next up previous
Next: Longitude difference Up: Geodesics on a spheroid Previous: Geodesics on a spheroid

Geodesic arc length

From Fig. 3 and Eqn. (17) we see that:

\begin{displaymath}
\cos\alpha=a(1-e^2\cos^2\theta)^{1/2}\frac{\mathrm{d}\theta}{\mathrm{d}s}=\pm(\cos^2\theta-\sin^2\alpha_0)^{1/2}/\cos\theta
\end{displaymath} (18)

and thus
\begin{displaymath}
\frac{\mathrm{d}s}{\mathrm{d}\theta}=\pm
a\frac{\cos\theta(1-e^2\cos^2\theta)^{1/2}}{(\cos^2\alpha_0-\sin^2\theta)^{1/2}}
\end{displaymath} (19)

with the sign being that of $\cos\alpha$. We now substitute $\sin\theta=\sin\sigma\cos\alpha_0$, where $a\sigma$ is the arc-length along the great circle on the auxiliary sphere, measured from where it crosses the equator in a Northerly direction, obtaining:
\begin{displaymath}
\frac{\mathrm{d}s}{\mathrm{d}\sigma}=a(1-e^2\cos^2\theta)^{1/2}=b(1+u^2\sin^2\sigma)^{1/2}
\end{displaymath} (20)

where $u^2\equiv e^2\cos^2\alpha_0/(1-e^2)$. By expanding in a power series in $u^2$ and integrating term by term[3,5], we obtain $s$ in the form:
$\displaystyle s/b =$   $\displaystyle \sigma(1+u^2/4-3u^4/64+5u^6/256-175u^8/16384+\ldots)$  
    $\displaystyle -\sin 2\sigma (u^2/8)(1-u^2/4+15u^4/128-35u^6/512+\ldots)$  
    $\displaystyle -\sin 4\sigma (u^4/256)(1-3u^2/4+35u^4/64-\ldots)$  
    $\displaystyle -\sin 6\sigma (u^6/3072)(1-5u^2/4+\ldots)$  
    $\displaystyle -\sin 8\sigma(5u^8/131072)(1-\ldots)$  
    $\displaystyle - \ldots$ (21)

The distance between two points $s(\sigma_2, \alpha_0)-s(\sigma_1,
\alpha_0)$, on a geodesic arc is best obtained, after differencing Eqn. (21), by using the identity $\sin( 2n\sigma_2) -\sin(
2n\sigma_1)=2\cos (2n\sigma_m) \sin n(\sigma_2-\sigma_1)$, where $\sigma_m=(\sigma_1+\sigma_2)/2$. This avoids excessive loss of significant digits when the two points are close together.

Vincenty[5] has rearranged a subset of the resulting equations into nested forms more suitable for computation:


$\displaystyle \tan\sigma_1$ $\textstyle =$ $\displaystyle \tan\phi_1/\cos\alpha_1$  
$\displaystyle \sin\alpha_0$ $\textstyle =$ $\displaystyle \cos\phi_1\sin\alpha_1$  
$\displaystyle u^2$ $\textstyle =$ $\displaystyle e^2 \cos^2\alpha_0/(1-e^2)$  
$\displaystyle A$ $\textstyle =$ $\displaystyle 1 +\frac{u^2}{16384}(4096+u^2(-768+u^2(320-175u^2)))$  
$\displaystyle B$ $\textstyle =$ $\displaystyle \frac{u^2}{1024}(256+u^2(-128+u^2(74-47u^2)))$  
$\displaystyle \sigma_m$ $\textstyle =$ $\displaystyle \sigma_1 +\sigma/2$  
$\displaystyle \Delta\sigma$ $\textstyle =$ $\displaystyle B\sin\sigma (\cos
2\sigma_m+\frac{B}{4}(\cos\sigma (-1+2\cos^2 2\sigma_m)$  
  $\textstyle \quad$ $\displaystyle -\frac{B}{6}\cos
2\sigma_m (-3 +4\sin^2\sigma)(-3+4\cos^2 2\sigma_m))))$  
$\displaystyle s$ $\textstyle =$ $\displaystyle bA(\sigma- \Delta\sigma)$ (22)

In equations (22) the origins of $\sigma$ and $s$ have been shifted from the equator to the initial point ($1$), where the reduced latitude is $\phi_1$ and the azimuth of the geodesic is $\alpha_1$.


next up previous
Next: Longitude difference Up: Geodesics on a spheroid Previous: Geodesics on a spheroid
Ed Williams
2002-03-21