[Proj] Differential behavior of geodesics

Charles Karney charles.karney at sri.com
Mon Jul 4 15:11:02 EST 2011

I've worked out some relationships between reduced lengths and geodesic
scales which I give here...

The story so far: Consider a geodesic from point 1 (azimuth azi1) to
point 2 with length s12.  Consider a second geodesic of length s12
starting from point 1 with azimuth azi1 + dazi1.  The end of this
geodesic lies a distance m12 * dazi1 from point 2.  It is easy to show
that m12 + m21 = 0.  m12 is called the "reduced length".

Consider a third geodesic of length s12 parallel to the first at point 1
but separated a distance dt1 from it.  At point 2 the geodesics are
separated by M12 * dt1.  M21 is defined likewise.  M12 and M21 are
called the "geodesic scales".

Together m12, M12, and M21 completely define the behavior of all
geodesics near the original one.

Derivatives holding point 1 and azi1 fixed and extending the length of
the geodesic:

    dm12/ds2 = M21
    dM12/ds2 = -(1 - M12*M21)/m12

Addition rules, points 1, 2, and 3 all lie on the same geodesic:

    m23 = m13*M12 - m12*M13
    M23 = M21*M13 + (1 - M12*M21) * m13/m12
    M32 = M31*M12 + (1 - M13*M31) * m12/m13

Hence we get Christoffel's addition rule (Theory of Geodesic Triangles,
Art 16, Eq 5) (points 1, 2, 3, and 4 all lie on the same geodesic):

    m12*m34 + m13*m42 + m14*m23 = 0

In the GIS domain, there an interest in defining "buffer regions" around
shapes.  If point 1 is on the original "reference" boundary, then a
point (point 2) on the buffer boundary is obtained by traveling a
distance s12 along a geodesic normal to the reference curve.  The buffer
boundary is therefore a geodesic parallel of the reference curve.

If the geodesic curvature of the reference curve at point 1 is K1, then
the geodesic curvature at point 2 on the buffer boundary is

    K2 = (M21*K1 - (1 - M12*M21)/m12) / (m12*K1 + M12)

If the reference curve has a corner in it, then K1 -> inf, and

    K2 = M21 / m12

is the curvature of a geodesic circle.  If K1 -> 0 then

    K2 = - (1 - M12*M21) / (M12*m12)

is the curvature of a curve parallel to a geodesic.  The test K2 -> inf
can to used to detect when the buffer boundary self-intersects and can
therefore be simplified.

These results are all general; they apply to all reasonable two
dimensional surfaces (not just ellipsoids of revolution).

Charles Karney <charles.karney at sri.com>
SRI International, Princeton, NJ 08543-5300
Tel: +1 609 734 2312

More information about the Proj mailing list