Angular Diameter Distances

In gravitational lensing studies, it is necessary to calculate angular diameter distances in function of redshift. On this page, you can perform such calculations.

The program

z1	Omega_m	w
z2	Omega_r
h	Omega_v

Angular diameter distance is: ?? Mpc

h determines H₀: H₀ = h 100 km s^-1Mpc^-1.
w specifies the equation of state for the so-called dark energy. For a Λ value which is truly a cosmological constant, this value is -1.
Omega_v, Omega_m, Omega_r specify contributions to the energy density.
z1 and z2: redshifts of objects between which the angular diameter distance should be calculated.

The calculation

The following PDF file derives the required calculations:

angdiam.pdf (98,652 bytes)

Comments, questions and error reports can be sent to the following e-mail address: jori.liesenborgs@gmail.com

The code

The JavaScript code is based on the following Java code. The routine can easily be adapted to your preferred programming language. Note that you have to make sure that z1 < z2 and W_m + W_r + W_v + W_k = 1, as no checks are performed to verify this. The value returned is the calculated distance in units of 1 Mpc.

final static double TH = 1.0/100000.0;
final static double c = 299792458.0;

static double calculateAngularDiameterDistance(double z1, double z2, 
   double h, double W_m, double W_r, double W_v, 
   double W_k, double w, int num)
{
   double sum = 0;
   double R2 = 1.0/(1.0+z1);
   double R1 = 1.0/(1.0+z2);
   double dR = (R2-R1)/((double)num);
   double R = R1;
   for (int i = 0 ; i < num ; i++, R += dR)
   {
      double term1 = W_v*Math.pow(R,1.0-3.0*w);
      double term2 = W_m*R;
      double term3 = W_r;
      double term4 = W_k*R*R;

      double val1 = 1.0/Math.sqrt(term1+term2+term3+term4);

      term1 = W_v*Math.pow(R+dR,1.0-3.0*w);
      term2 = W_m*(R+dR);
      term3 = W_r;
      term4 = W_k*(R+dR)*(R+dR);

      double val2 = 1.0/Math.sqrt(term1+term2+term3+term4);

      sum += ((val1+val2)/2.0)*dR; // trapezium rule
   }

   double A = 0;

   if (W_k == 0)
      A = sum;
   else if (W_k > 0)
      A = (1.0/Math.sqrt(W_k))*Math.sinh(Math.sqrt(W_k)*sum);
   else // W_k < 0
      A = (1.0/Math.sqrt(-W_k))*Math.sin(Math.sqrt(-W_k)*sum);

   double result = c*A*((1.0/(1.0+z2))*TH)/h;
   return result;
}