Examples
This implementation of the geodesic parts of GeographicLib allows you to do the following:
- Compute azimuths and distances between two know points: use
inverse
. - Compute the end point from some starting point, a direction and distance: use
forward
orforward_deg
. - Calculate waypoints along a known great circle path: use
GeodesicLine
andwaypoints
. - Compute the area or perimeter of a polyon on the ellipsoid: use
Polygon
andproperties
.
These examples show how to do each of these in practice.
(Note that these examples are taken from the original documentation.)
Ellipsoids
By default, great circle calculations use the WGS84 ellipsoid. You therefore do not need to specify the first ellipse
argument to the functions unless you want to specify your own ellipsoid. However, to specify your own ellipsoid on which to compute great circles, this is easily done:
julia> using GeographicLib
julia> semimajor_radius_m = 6378_388;
julia> flattening = 1/297.0;
julia> ellipse = Geodesic(semimajor_radius_m, flattening)
Geodesic: a: 6.378388e6 f: 0.003367003367003367
Distance between two points
If you know the location of two points in longitude and latitude, use the inverse
function to find the distance between them, the azimuth from the first point and backazimuth from the second point to the first:
julia> lon1, lat1 = 174.81, -41.32; # Wellington, New Zealand
julia> lon2, lat2 = -5.50, 40.96; # Salamanca, Spain
julia> inverse(ellipse, lon1, lat1, lon2, lat2)
(azi = 161.1046750858827, baz = -161.21158298105982, dist = 1.9960331457114425e7, angle = 179.6197913140605)
This returns a named tuple; the distance in m is given by the .dist
field, .angle
gives the angular distance on the equivalent sphere, .azi
gives the forward azimuth and .baz
the backazimuth.
Moving a set distance from one point
If you know the starting point and want to know the end point from moving a certain distance along a certain azimuth, use forward
.
julia> forward(ellipse, 115.74, -32.06, 225, 2000e3) # 2000 km southwest of Perth, Australia
(lon = 98.22644722385655, lat = -43.646489349577195, baz = 55.854896821905356, dist = 2.0e6, angle = 18.003186633674353)
If you want to know the end point a certain angular distance away, use forward_deg
:
julia> forward_deg(ellipse, 115.74, -32.06, 225, 20) # 20° southwest of Perth
(lon = 95.90742651866789, lat = -44.744151411150156, baz = 57.47168951862605, dist = 2.221907360203232e6, angle = 20.0)
Computing waypoints
Consider the great circle path from Beijing Airport (116.6°E, 40.1°N) to San Fransisco Airport (122.4°W, 37.6°N). Compute waypoints and azimuths at intervals of 1000 km by creating a GeodesicLine
and then using waypoints
:
julia> line = GeodesicLine(ellipse, 116.6, 40.1, lon2=-122.4, lat2=37.6)
GeodesicLine: a: 6.378388e6 f: 0.003367003367003367 caps: 65439 lat1: 40.1 lon1: 116.6 azi1: 42.916274498536474 salp1: 0.6809289156056745 calp1: 0.7323495148438894 s13: 9.514421681743788e6 a13: 85.57941157162365
julia> waypoints(line, dist=1000e3)
11-element Vector{@NamedTuple{lon::Float64, lat::Float64, baz::Float64, dist::Float64, angle::Float64}}: (lon = 116.6, lat = 40.1, baz = -137.08372550146353, dist = 0.0, angle = 0.0) (lon = 125.4485667325373, lat = 46.37304779242446, baz = -131.00682481965458, dist = 1.0e6, angle = 8.998886872620938) (lon = 136.40637727605085, lat = 51.78759841880426, baz = -122.70667626241067, dist = 2.0e6, angle = 17.994676114113027) (lon = 149.93620398067395, lat = 55.924144706930434, baz = -111.75604881313774, dist = 3.0e6, angle = 26.987984341289916) (lon = 165.90470996945186, lat = 58.27453572396102, baz = -98.32022506798336, dist = 4.0e6, angle = 35.97966822059779) (lon = -176.97213147749545, lat = 58.43545906442418, baz = -83.71311081540205, dist = 5.0e6, angle = 44.970741063118155) (lon = -160.7345701291212, lat = 56.375367574145756, baz = -70.00410754285849, dist = 6.0e6, angle = 53.96227499942639) (lon = -146.83061753761507, lat = 52.459361612570746, baz = -58.670924994668155, dist = 7.0e6, angle = 62.9552977797524) (lon = -135.531489855743, lat = 47.19657788053033, baz = -50.01634694411308, dist = 8.0e6, angle = 71.95069362742086) (lon = -126.42034205680477, lat = 41.0241211786531, baz = -43.65845928486635, dist = 9.0e6, angle = 80.94911728584405) (lon = -122.4, lat = 37.59999999999999, baz = -41.10962796089933, dist = 9.514421681743788e6, angle = 85.57941157162367)
If you just want some number of points along the way, use the n
keyword argument:
julia> waypoints(line, n=10)
10-element Vector{@NamedTuple{lon::Float64, lat::Float64, baz::Float64, dist::Float64, angle::Float64}}: (lon = 116.6, lat = 40.1, baz = -137.08372550146353, dist = 0.0, angle = 0.0) (lon = 126.0126111856647, lat = 46.709034579241305, baz = -130.59740008727215, dist = 1.0571579646381987e6, angle = 9.51314678540973) (lon = 137.8175199510943, lat = 52.33432908010009, baz = -121.59372287667333, dist = 2.1143159292763975e6, angle = 19.022865841326787) (lon = 152.52050213344623, lat = 56.46815877414928, baz = -109.608496596686, dist = 3.1714738939145957e6, angle = 28.52991214184784) (lon = 169.78803691430426, lat = 58.512936752962524, baz = -95.01255785530321, dist = 4.228631858552795e6, angle = 38.03532928611469) (lon = -172.15587006256067, lat = 58.06322940031036, baz = -79.61693059230413, dist = 5.285789823190993e6, angle = 47.54033622959181) (lon = -155.6674221274284, lat = 55.21594792543654, baz = -65.81230407755629, dist = 6.342947787829191e6, angle = 57.04619601194877) (lon = -142.01246922996296, lat = 50.48799877187304, baz = -54.90015839311627, dist = 7.40010575246739e6, angle = 66.55407998760973) (lon = -131.12820808268523, lat = 44.465527833955235, baz = -46.856268007749804, dist = 8.45726371710559e6, angle = 76.06494142269099) (lon = -122.4, lat = 37.59999999999999, baz = -41.10962796089933, dist = 9.514421681743788e6, angle = 85.57941157162367)
Or set the distance between points in terms of angle:
julia> waypoints(line, angle=10)
10-element Vector{@NamedTuple{lon::Float64, lat::Float64, baz::Float64, dist::Float64, angle::Float64}}: (lon = 116.6, lat = 40.1, baz = -137.08372550146353, dist = 0.0, angle = 0.0) (lon = 126.55306724811369, lat = 47.02451056726837, baz = -130.20299081102675, dist = 1.1112708185226605e6, angle = 10.0) (lon = 139.1907543571846, lat = 52.837886783177815, baz = -120.50297777939181, dist = 2.222958662728229e6, angle = 20.0) (lon = 155.05072606130204, lat = 56.93534454913972, baz = -107.49355294926323, dist = 3.334963135588773e6, angle = 30.0) (lon = 173.53698469358338, lat = 58.62977452059097, baz = -91.81330485164187, dist = 4.447145746945927e6, angle = 40.0) (lon = -167.6599068447789, lat = 57.540243916824664, baz = -75.81180965562544, dist = 5.559346588636242e6, angle = 50.0) (lon = -151.1045662211484, lat = 53.91941776870715, baz = -62.09358473950738, dist = 6.671403561172048e6, angle = 60.0) (lon = -137.77647639652707, lat = 48.428003710944736, baz = -51.67995925549309, dist = 7.783171861221558e6, angle = 70.0) (lon = -127.29504443962117, lat = 41.70777455951728, baz = -44.236536959892305, dist = 8.894541409303235e6, angle = 80.0) (lon = -122.4, lat = 37.6, baz = -41.10962796089939, dist = 9.514421681743788e6, angle = 85.57941157162365)
Measuring areas
Measure the area of Antarctica by using a Polygon
and properties
:
julia> antarctica = [ ( -58,-63.1), (-74,-72.9), (-102,-71.9), (-102,-74.9), (-131,-74.3), (-163,-77.5), (163,-77.4), ( 172,-71.7), ( 140,-65.9), ( 113,-65.7), ( 88,-66.6), ( 59,-66.9), ( 25,-69.8), ( -4,-70.0), ( -14,-71.0), ( -33,-77.3), (-46,-77.9), ( -61,-74.7) ];
julia> lons = first.(antarctica);
julia> lats = last.(antarctica);
julia> p = Polygon(ellipse, lons, lats);
julia> properties(p)
(n = 18, perimeter = 1.683106789279071e7, area = 1.3662703680020125e13)
Equivalently, call add_point!
multiple times:
julia> p = Polygon(ellipse);
julia> for (lon, lat) in antarctica add_point!(p, lon, lat) end
julia> properties(p)
(n = 18, perimeter = 1.6831928215597793e7, area = 1.3664118444349844e13)