如何使用python计算地球表面上多边形的面积？

假设您以GeoJSON格式表示了科罗拉多州

{"type": "Polygon",  "coordinates": [[   [-102.05, 41.0],    [-102.05, 37.0],    [-109.05, 37.0],    [-109.05, 41.0] ]]}

所有坐标均为经度，纬度。您可以使用pyproj投影坐标，并使用Shapely查找任何投影多边形的面积：

co = {"type": "Polygon", "coordinates": [    [(-102.05, 41.0),     (-102.05, 37.0),     (-109.05, 37.0),     (-109.05, 41.0)]]}lon, lat = zip(*co['coordinates'][0])from pyproj import Projpa = Proj("+proj=aea +lat_1=37.0 +lat_2=41.0 +lat_0=39.0 +lon_0=-106.55")

这是一个以关注区域为中心并包围其关注区域的等面积投影。现在制作新的投影GeoJSON表示形式，将其转换为Shapely几何对象，并采用以下区域：

x, y = pa(lon, lat)cop = {"type": "Polygon", "coordinates": [zip(x, y)]}from shapely.geometry import shapeshape(cop).area  # 268952044107.43506

这是与被调查区域非常接近的近似值。对于更复杂的特征，您需要沿顶点之间的边缘进行采样，以获取准确的值。上面有关日期变更线等的所有警告均适用。如果您仅对区域感兴趣，则可以在投影之前将特征从日期线移开。