This is a litle benchmark of three Python libraries to compute geodesic distances:
Geodesic distance computation with Pyproj (python interface to PROJ4 C library).
https://github.com/jswhit/pyproj
PROJ4 C library (>= v.4.9.0) routines used to compute geodesic distances are a simple transcription from C. Karney Geographiclib C++ Library.
Geodesic distance calculations using Charles Karney geodesic algorithms: C. F. F. Karney, Algorithms for Geodesics, J. Geodesy 87(1), 43–55 (Jan. 2013)
https://trac.osgeo.org/proj/browser/trunk/proj/src/geodesic.h
PROJ4 C library (< v.4.9.0) routines used to compute geodesic distances are: Paul D. Thomas, Spheroidal Geodesics, Reference Systems, and Local Geometry. U.S. Naval Oceanographic Office, p. 162 Engineering Library 526.3 T36s (1970)
Default values used (WGS84): https://github.com/jswhit/pyproj/blob/master/lib/pyproj/__init__.py
Therefore:
from pyproj import Geod
def grcrcl1(lon_1, lat_1, lon_2, lat_2):
g = Geod(ellps='WGS84')
az, az_inv, dist = g.inv(lon_1, lat_1, lon_2, lat_2)
return dist
Geodesic distance computation with Pygc.
https://github.com/axiom-data-science/pygc
Geodesic distance calculations using Vincenty's formulae: http://en.wikipedia.org/wiki/Vincenty%27s_formulae
Default values used (WGS84): https://github.com/axiom-data-science/pygc/blob/master/pygc/gc.py
Therefore:
from pygc import great_distance
def grcrcl2(startlong, startlat, endlong, endlat):
res = great_distance(start_latitude=startlat, start_longitude=startlong,
end_latitude=endlat, end_longitude=endlong)
return float(res['distance'])
Geodesic distance computation with GeographicLib (python interface to GeographicLib C++ library).
http://geographiclib.sourceforge.net/html/other.html#python
Geodesic distance calculations using Charles Karney geodesic algorithms: C. F. F. Karney, Algorithms for Geodesics, Journal of Geodesy 87(1), 43–55 (Jan. 2013)
Default values used (WGS84): http://geographiclib.sourceforge.net/html/Constants_8hpp_source.html
Therefore:
from geographiclib.geodesic import Geodesic
def grcrcl3(lon_1, lat_1, lon_2, lat_2):
res = Geodesic.WGS84.Inverse(lat_1,lon_1,lat_2,lon_2)
return res['s12']
Function to print results from geodesic distance functions:
def printResults(input_coords):
lon_1, lat_1, lon_2, lat_2 = input_coords
dist1 = grcrcl1(lon_1, lat_1, lon_2, lat_2)
dist2 = grcrcl2(lon_1, lat_1, lon_2, lat_2)
dist3 = grcrcl3(lon_1, lat_1, lon_2, lat_2)
print "Distances for: lon_1 {0}, lat_1 {1}, lon_2 {2}, lat_2 {3}\n".format(lon_1, lat_1, lon_2, lat_2)
print "Pyproj___________________________{:,.10f} m\n".format(dist1)
print "Pygc_____________________________{:,.10f} m\n".format(dist2)
print "GeographicLib____________________{:,.10f} m\n".format(dist3)
print "Difference Pyproj,Pygc___________{:.15f} m\n".format(abs(dist1 - dist2))
print "Difference Pyproj,GeographicLib__{:.15f} m\n".format(abs(dist1 - dist3))
print "Difference GeographicLib,Pygc____{:.15f} m\n".format(abs(dist3 - dist2))
print "------------------------------------------\n"
Creating test data. Test data is a list of lists. Each inner list contains:
data = [
[-3.6,40.5,-118.4,33.9],
[-6.,37.,-145.,11.],
[-150.,37.,140.,11.],
[-50.,7.,40.,11.],
[-100.,80.,-140.,30.]
]
Printing three functions results with created test data:
map(printResults, data)
Distances for: lon_1 -3.6, lat_1 40.5, lon_2 -118.4, lat_2 33.9 Pyproj___________________________9,406,468.4536480606 m Pygc_____________________________9,406,468.4522715993 m GeographicLib____________________9,406,468.4536480606 m Difference Pyproj,Pygc___________0.001376461237669 m Difference Pyproj,GeographicLib__0.000000000000000 m Difference GeographicLib,Pygc____0.001376461237669 m ------------------------------------------ Distances for: lon_1 -6.0, lat_1 37.0, lon_2 -145.0, lat_2 11.0 Pyproj___________________________13,189,756.0717643406 m Pygc_____________________________13,189,756.0719679464 m GeographicLib____________________13,189,756.0717643406 m Difference Pyproj,Pygc___________0.000203605741262 m Difference Pyproj,GeographicLib__0.000000000000000 m Difference GeographicLib,Pygc____0.000203605741262 m ------------------------------------------ Distances for: lon_1 -150.0, lat_1 37.0, lon_2 140.0, lat_2 11.0 Pyproj___________________________7,510,240.7866536863 m Pygc_____________________________7,510,240.7865020484 m GeographicLib____________________7,510,240.7866536863 m Difference Pyproj,Pygc___________0.000151637941599 m Difference Pyproj,GeographicLib__0.000000000000000 m Difference GeographicLib,Pygc____0.000151637941599 m ------------------------------------------ Distances for: lon_1 -50.0, lat_1 7.0, lon_2 40.0, lat_2 11.0 Pyproj___________________________9,871,047.1972853243 m Pygc_____________________________9,871,047.1963992771 m GeographicLib____________________9,871,047.1972853225 m Difference Pyproj,Pygc___________0.000886047258973 m Difference Pyproj,GeographicLib__0.000000001862645 m Difference GeographicLib,Pygc____0.000886045396328 m ------------------------------------------ Distances for: lon_1 -100.0, lat_1 80.0, lon_2 -140.0, lat_2 30.0 Pyproj___________________________5,853,723.4893777100 m Pygc_____________________________5,853,723.4893133771 m GeographicLib____________________5,853,723.4893777100 m Difference Pyproj,Pygc___________0.000064332969487 m Difference Pyproj,GeographicLib__0.000000000000000 m Difference GeographicLib,Pygc____0.000064332969487 m ------------------------------------------
[None, None, None, None, None]
As it was expected there are no difference betwen Pyproj and GeographicLib outcomes. However there are a minimal difference between Pygc and Pyproj/Geographiclib.