Angular Diameter Distances

In gravitational lensing studies, it is necessary to calculate angular diameter distances in function of redshift. On this page, you'll find a Java applet which enables you to perform such calculations.

The applet

The applet uses the Open Source Physics package. If your Java settings allow it, pressing the button below will cause a window to appear in which several parameters kan be specified:

  • h determines H0: H0 = 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 and Omega_k specify contributions to the energy density. Their sum has to be equal to 1.
  • z1 and z2: redshifts of objects between which the angular diameter distance should be calculated. Note that z2 should be larger than z1.
  • Calculating the angular diameter distance requires an integral to be calculated and steps specifies the number of integration steps used to calculate this integral numerically. The more steps, the more accurate the result, but the longer it takes to calculate.

The calculation

The following PDF file derives the required calculations:

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

The code

The following Java code is used in the calculation. 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;
}