Geodesics
A geodesic is the shortest path between two points along a manifold. On the ellipsoid there are two classical problems, and each has its own function. Both are solved with the series of Karney (2013), on the ellipsoid attached to the datum of the coordinates.
Meshes.geodesicfwd — Function
geodesicfwd(p, ϕ, l)Solve the direct geodesic problem: return the point reached from the point p after walking the length l along the geodesic that leaves p with the azimuth ϕ, measured clockwise from the north.
The geodesic is the one of the ellipsoid attached to the datum of p, and is computed with the series of Karney (2013).
See also geodesicbwd.
Examples
p = Point(LatLon(0, 0))
geodesicfwd(p, 90u"°", 1000u"km")
geodesicfwd(p, 90, 1000000)References
- Karney, C. F. F. 2013. Algorithms for geodesics
Meshes.geodesicbwd — Function
geodesicbwd(p₁, p₂)Solve the inverse geodesic problem: return the azimuth at the point p₁ of the shortest geodesic connecting p₁ to the point p₂, measured clockwise from the north.
The length of that geodesic is GeodesicDistance, and the azimuth at the other end is geodesicbwd(p₂, p₁).
See also geodesicfwd.
Examples
p₁ = Point(LatLon(0, 0))
p₂ = Point(LatLon(0, 1))
geodesicbwd(p₁, p₂)References
- Karney, C. F. F. 2013. Algorithms for geodesics
Meshes.geodesictangent — Function
geodesictangent(p, ϕ)Unit vector tangent to the ellipsoid at the point p, pointing along the azimuth ϕ, measured clockwise from the north.
The vector is expressed in the geocentric Cartesian coordinates of the datum of p, and is the direction in which geodesicfwd walks.
See also geodesicazimuth.
Examples
p = Point(LatLon(0, 0))
geodesictangent(p, 90u"°")
geodesictangent(p, 90)Meshes.geodesicazimuth — Function
geodesicazimuth(p, v)Azimuth of the vector v at the point p of the ellipsoid, measured clockwise from the north. Any component of v along the normal of the ellipsoid is ignored, and the azimuth is undefined if nothing is left.
See also geodesictangent.
Examples
p = Point(LatLon(0, 0))
geodesicazimuth(p, Vec(0, 1, 0))The direct problem walks a given length along a given azimuth. Plain numbers are taken to be degrees and meters:
p = Point(LatLon(-33.8688, 151.2093))
geodesicfwd(p, 45, 1000000)Point with GeodeticLatLon{WGS84Latest} coordinates
├─ lat: -27.28175431149266°
└─ lon: 158.34048229360366°The inverse problem recovers the azimuth that connects two points. Its length is GeodesicDistance, so the two together describe the geodesic completely:
sydney = Point(LatLon(-33.8688, 151.2093))
london = Point(LatLon(51.5074, -0.1278))
geodesicdist = GeodesicDistance()
ϕ = geodesicbwd(sydney, london)
d = geodesicdist(sydney, london)
ϕ, d(-40.70736698782122°, 1.698929577054045e7 m)Walking that azimuth for that length arrives at the second point:
geodesicfwd(sydney, geodesicbwd(sydney, london), geodesicdist(sydney, london))Point with GeodeticLatLon{WGS84Latest} coordinates
├─ lat: 51.5074°
└─ lon: -0.12780000000000769°The azimuth at the far end points back along the same geodesic, and is obtained by swapping the arguments:
geodesicbwd(london, sydney)60.361685700155675°Tangent vectors
An azimuth at a point can also be expressed as a unit vector tangent to the ellipsoid there, in geocentric Cartesian coordinates:
geodesictangent(sydney, geodesicbwd(sydney, london))Vec(-0.05612928822352242 m, 0.7750343623232548 m, 0.6294213534840605 m)The conversion goes both ways, and any component of the vector along the normal of the ellipsoid is ignored:
geodesicazimuth(sydney, geodesictangent(sydney, 45))45.00000000000001°