线条优美。长度单位

Question

我试图在Shapely中找到线串的长度（以米为单位），但似乎无法达到预期的结果。几乎可以肯定我在坐标系方面犯了一些错误，但我一直没能弄明白

下面是单行的一些简化代码：

from shapely.geometry import LineString
line = LineString([(12875996.563923, -3940011.116702), (12872802.929335, -3937989.118438)])
line.crs = 'EPSG:3857'
line.length

输出为3779.92米。然而，根据我的测量，它应该是约3159米

这是基于检查以下站点点的位置（两个环形交叉口的中心），然后在Google Earth中测量它们之间的距离

第1点： http://epsg.io/map#srs=3857&x=12875996.563923&y=-3940011.116702&z=14&layer=streets
第2点： http://epsg.io/map#srs=3857&x=12872802.929468&y=-3937989.118538&z=17&layer=streets

Answer 1

为了证明，shapely几何中的line.length只是一个欧几里德距离，可以这样计算：

import math
x1, y1 = (12875996.563923, -3940011.116702)
x2, y2 = (12872802.929335, -3937989.118438)
math.hypot(x2-x1,  y2-y1)  # 3779.91783790157

要使用pyproj模块获得正确的距离，请执行以下步骤：

import pyproj  #v2.4.1
x1, y1 = (12875996.563923, -3940011.116702)
x2, y2 = (12872802.929335, -3937989.118438)
lon1, lat1 = pyproj.Proj("epsg:3857")(x1, y1, inverse=True)
lon2, lat2 = pyproj.Proj("epsg:3857")(x2, y2, inverse=True)
_,_,dist_km = pyproj.Geod(ellps='WGS84').inv(lon1, lat1, lon2, lat2)
dist_km   # 3157.214113925091