我想通过使用库中的Geod类来计算两个lon / lat点之间的距离pyproj。
Geod
pyproj
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])
我收到以下错误:
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)
此错误消息从何而来?
这两点仅相距几厘米。看起来pyproj/Geod不能很好地解决彼此靠近的点。这有点奇怪,因为在这样的距离下简单的平面几何体已经足够了。同样,该错误消息也有点可疑,因为这表明这两点是对立的,即在直径上是相反的,显然不是这种情况!OTOH,也许它所提到的对立点是在计算中以某种方式出现的中间点…不过,我还是会犹豫使用这种行为的库。
鉴于此缺陷,我怀疑pyproj还有其他缺陷。特别是,它可能会使用旧的Vincenty公式进行椭球测地线计算,众所周知,当处理近对映点时该算法不稳定,并且在大距离上并非特别精确。我建议使用CFF Karney的现代算法。
Karney博士是一个主要贡献者上测地线的维基百科文章,特别是对一个椭球大地测量,他的geographiclib可PyPI上,这样你就可以很容易地使用安装它pip。请参阅他的SourceForge网站以获取更多信息以及其他语言的geolib绑定。
pip
FWIW,这是一个使用geoliblib计算问题距离的简短演示。
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'])
输出
0.07345528623159624
该数字以米为单位,因此这两点相距73毫米多一点。
如果您希望看到geoliblib用于解决复杂的测地线问题,请参阅我去年在gist上使用Python 2/3源代码编写的math.stackexchange答案。
希望这不再是问题,因为pyproj现在使用来自geoliblib的代码。