Welcome to OGeek Q&A Community for programmer and developer-Open, Learning and Share
Welcome To Ask or Share your Answers For Others

Categories

0 votes
483 views
in Technique[技术] by (71.8m points)

python - Geod ValueError : undefined inverse geodesic

I want to compute the distance between two lon / lat points by using Geod class from pyproj library.

from pyproj import Geod

g = Geod(ellps='WGS84')
lonlat1 = 10.65583081724002, -7.313341167341917
lonlat2 = 10.655830383300781, -7.313340663909912

_, _, dist = g.inv(lonlat1[0], lonlat1[1], lonlat2[0], lonlat2[1])

I get the following error :

ValueError                                Traceback (most recent call last)
<ipython-input-5-8ba490aa5fcc> in <module>()
----> 1 _, _, dist = g.inv(lonlat1[0], lonlat1[1], lonlat2[0], lonlat2[1])

/usr/lib/python2.7/dist-packages/pyproj/__init__.pyc in inv(self, lons1, lats1, lons2, lats2, radians)
    558         ind, disfloat, dislist, distuple = _copytobuffer(lats2)
    559         # call geod_inv function. inputs modified in place.
--> 560         _Geod._inv(self, inx, iny, inz, ind, radians=radians)
    561         # if inputs were lists, tuples or floats, convert back.
    562         outx = _convertback(xisfloat,xislist,xistuple,inx)

_geod.pyx in _geod.Geod._inv (_geod.c:1883)()

ValueError: undefined inverse geodesic (may be an antipodal point)

Where does this error message come from ?

See Question&Answers more detail:os

与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
Welcome To Ask or Share your Answers For Others

1 Reply

0 votes
by (71.8m points)

Those two points are only a few centimetres apart. It looks like pyproj / Geod doesn't cope well with points which are that close together. That's a bit strange, since simple plane geometry is more than adequate at such distances. Also, that error message is a bit suspicious, since it's suggesting that the two points are antipodal, i.e., diametrically opposite, which is clearly not the case! OTOH, maybe the antipodal point it mentions is some intermediate point that arises somehow in the calculation... Still, I'd be rather hesitant in using a library that behaves like this.

Given this defect, I suspect that pyproj has other flaws. In particular, it probably uses the old Vincenty's formulae for its ellipsoid geodesic calculations, which is known to be unstable when dealing with near-antipodal points, and not particularly accurate over large distances. I recommend using the modern algorithms of C. F. F. Karney.

Dr Karney is a major contributor to the Wikipedia articles on geodesics, in particular Geodesics on an ellipsoid, and his geographiclib is available on PyPi, so you can easily install it using pip. See his SourceForge site for further information, and geographiclib binding in other languages.

FWIW, here's a short demo of using geographiclib to compute the distance in your question.

from geographiclib.geodesic import Geodesic

Geo = Geodesic.WGS84

lat1, lon1 = -7.313341167341917, 10.65583081724002
lat2, lon2 = -7.313340663909912, 10.655830383300781

d = Geo.Inverse(lat1, lon1,  lat2, lon2)
print(d['s12'])

output

0.07345528623159624

That figure is in metres, so those two points are a little over 73mm apart.


If you'd like to see geographiclib being used to solve a complex geodesic problem, please see this math.stackexchange answer I wrote last year, with Python 2 / 3 source code on gist.


Hopefully, this is no longer an issue, since pyproj now uses code from geographiclib.


与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
OGeek|极客中国-欢迎来到极客的世界,一个免费开放的程序员编程交流平台!开放,进步,分享!让技术改变生活,让极客改变未来! Welcome to OGeek Q&A Community for programmer and developer-Open, Learning and Share
Click Here to Ask a Question

...