按形状以平面单位(例如平方米)计算多边形面积

2022-04-11 00:00:00 python shapely polygon area

问题描述

我正在使用Python3.4和Shapely 1.3.2从一系列经纬度坐标对创建一个Polygon对象，我将这些坐标对转换为熟知的文本字符串，以便对其进行解析。这样的多边形可能如下所示：

POLYGON ((-116.904 43.371, -116.823 43.389, -116.895 43.407, -116.908 43.375, -116.904 43.371))

由于Shapely不处理任何投影并实现笛卡尔空间中的所有几何对象，因此在该多边形上调用Area方法，如下所示：

poly.area

以平方度为单位给出了那个多边形的面积。要以平方米这样的平面单位计算面积，我想我必须使用不同的投影来转换多边形的坐标(哪一个？)。

我多次阅读到pyproj库应该提供实现这一点的方法。使用pyproj，有没有办法将一个形状整齐的多边形物体转换成另一个投影，然后计算出面积？

我对我的多边形做了一些其他工作(不是你现在想的那样)，只有在某些情况下，我才需要计算面积。

到目前为止，我只找到了这个例子： http://all-geo.org/volcan01010/2012/11/change-coordinates-with-pyproj/

这意味着将每个Polygon对象拆分为其外环和内环(如果存在)，获取坐标，将每对坐标转换为另一个投影并重新构建Polygon对象，然后计算其面积(那么它到底是什么单位？)这看起来像是一个解决方案，但不太实用。

有更好的主意吗？

解决方案

好了，我终于用matplotlib库的底图工具包完成了。我会解释它是如何工作的，也许这会在某个时候对某人有所帮助。

1. 在您的系统上下载并安装matplotlib库。 http://matplotlib.org/downloads.html

对于Windows二进制文件，我建议使用此页面： http://www.lfd.uci.edu/~gohlke/pythonlibs/#matplotlib 注意提示：

需要NumPy、Dateutil、PYTZ、PYPARING、Six和可选Pillow， Pycairo、tornado、wxpython、pyside、pyqt、ghost script、miktex、ffmpeg、 Mencoder、avconv或Imagemagick。

因此，如果您的系统上还没有安装，您必须下载并安装NumPy、Dateutil、PYTZ、PYPARSING和Six，matplotlib才能正常工作(对于Windows用户：可以在页面上找到它们，对于Linux用户，可以使用python包管理器"pip"来完成此任务)。

2. 下载并安装matplotlib的"底图"工具包。不是来自 http://matplotlib.org/basemap/users/installing.html 或者对于Windows二进制文件，也可以从以下位置获得： http://www.lfd.uci.edu/~gohlke/pythonlibs/#basemap

3. 在您的Python代码中进行投影：

#Import necessary libraries
from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np

#Coordinates that are to be transformed
long = -112.367
lat = 41.013

#Create a basemap for your projection. Which one you use is up to you.
#Some examples can be found at http://matplotlib.org/basemap/users/mapsetup.html
m = Basemap(projection='sinu',lon_0=0,resolution='c')

#Project the coordinates:
projected_coordinates = m(long,lat)

输出：

投影坐标 (10587117.191355567,14567974.064658936)

就这么简单。现在，当使用投影坐标来构建一个带有Shapely的多边形，然后通过Shapely的面积方法计算面积时，您将获得以平方米为单位的面积(根据您使用的投影)。要得到平方公里，除以10^6。

编辑：我努力不仅要转换单个坐标，而且要转换整个几何体对象(如多边形)，因为这些对象已经作为形状对象提供，而不是通过它们的原始坐标。这意味着要向

编写大量代码

获取多边形外环的坐标
确定多边形是否有洞，如果有，请分别处理每个洞
变换外环和任何孔的每对坐标
将整个物体放回原处，并使用投影坐标创建一个多边形对象
这仅适用于多边形...为多面体和几何集合添加更多循环

然后，我偶然发现了这部分精美的文档，它使事情变得容易得多： http://toblerity.org/shapely/manual.html#shapely.ops.transform

当设置投影图时，例如如上所述：

m = Basemap(width=1,height=1, resolution='l',projection='laea',lat_ts=50,lat_0=50,lon_0=-107.)

然后，可以使用此投影通过以下方式变换任何形状良好的几何对象：

from shapely.ops import transform
projected_geometry = transform(m,geometry_object)

相关文章