目录

上一个主题

5. 矢量数据的空间分析:使用Shapely

下一个主题

5.2. Shapely 中几何要素的操作

关注公众号


常见问题

  1. Windows下的安装说明
  2. Jupyter免费在线实验环境
  3. 勘误与补充


>>> from helper import info; info()
页面更新时间: 2019-10-28 20:02:03
操作系统/OS: Linux-4.19.0-6-amd64-x86_64-with-debian-10.1
Python: 3.7.3

5.1. 空间数据模型

5.1.1. 数据模型

由Shapely提供的几何对象的基本类型包括点, 曲线,和曲面。 每一种都与三套平面的点 (可能无限)相联系。 一个要素的内部, 边界, 外部是相互排斥的,并且它们的并集与整个平面重叠。

  • 点的内部有一个点,边界没有点,点的外部包括其它所有点。点的拓扑维数是0。

  • 曲线有沿其长度的无穷多个点组成的内部集 (想象一个点在空间拖动), 包括2个端点的边界集, 和包括其它所有点的外部集. 曲线的拓扑维数是1。

  • 曲面有一个包含无穷个点内部集 (想象一条曲线在空间上拖动以覆盖一个区域), 有一个由一个或多个曲面组成的边界, 和包含所有其它点的外部集。包括那些所有在表面可能存在的洞.。曲面的拓扑维数是2。

点是由点类实现的; 曲线是由直线和环线类实现的; 曲面是由多边形类实现的。 Shapely不执行平滑的 ( 即有连续的切线)曲线。 所有的曲线必须近似线性样条曲线。 所有的圆斑必须近似线性样条范围内的地区。

点的集合由一多样点类来实现, 曲线的集合由一多样线类来实现, 曲面的集合由一多样多边形来实现. 这些集合在计算上作用不是很大, 但对某些种类的功能建模非常有用。例如,一个Y形线 可以很好地由多样线建模为一个整体。

标准的数据模型对几何对象的特定类型 有特定的附加约束, 这将在接下来的章节里讨论。

5.1.2. 一些术语

  • 坐标, 在已经定义好的精确度模型中可以被精确地描绘的空间中的一点。

  • 精确计算, 在操作过程中支持所有数字的数值计算,通常要用到大量复杂的算法。

  • 节点, 相同或不同的几何图形内两条线相交的点。这点不必用坐标表示出来,因为大部分的相交点的输出计算比输入计算需要更加严格的精确度。

  • 节点的计算, 计算出在一个或多个几何体相交之处形成的节点的过程。

  • 非坐标, 不可被描绘成一个坐标的点。

  • 数值稳定性, 算法的稳定性取决于它的输出结果中误差的最大范围。如果算法的误差范围是很小的,那么可以认为它是稳定的。

  • 点, R3中任意一点,一般来说,可以无限地描绘出来。

  • 真交点, 两条线段相交所形成的唯一交点,并且是在这两条线段上的点。

  • 健壮计算, 一种对于所有的输入都必定会输出正确答案的数值计算方法,通常需要有特别处理舍入误差而设的算法。

  • SFS, 开放式地理信息系统协会OGC的简单要素实现规范

  • 分辨率单位, 在已经定义好的精确度模型中的最小的可以表示出来的距离。

  • 顶点, 几何体中的边角点,这些点的坐标都可以确定的话,就可以确定一个几何体的位置。

5.1.3. Shapely的基本技术特征

  1. 得到二维形状的联合分析,交叉分析,或者差异分析的结果;

  2. 决定二维形状是否相交,接触,或者一个中包含另一个。

以及更多。

关系

空间数据模型有一组几何对象之间关系的自然语言 – 包含,相交,叠置,相切等– 和一个理论框架用于理解它们 使用它们的组成点集的相互交叉的3x3矩阵。

运算符

承接JTS技术规格, 本手册将区别构造(缓冲,凸壳)和 理论集操作(相交,联合等)。

坐标系

虽然地球是不平的,但是有很多分析类的问题能够通过 试验的和真正的算法,将地球特征转化到直角平面, 接着将结果转为地理坐标来解决。这种做法与传统的精确纸质地图一样古老。

Shapely不支持坐标系统转换。 所有2个或更多特征的操作假定特征在同一个直角平面。

空间操作与计算

GEOS主要支持的操作和计算:

空间关系计算,主要支持的计算

OGC Filter定义了一些xml标签,从而实现了查询Feature集合的接口。 程序可以将它转换为查询语言,检索Feature集合后返回有效的Feature子集合。

  • 相等(Equals)几何形状拓扑上相等。

  • 脱节(Disjoint)几何形状没有共有的点。

  • 相交(Intersects)几何形状至少有一个共有点(区别于脱节)

  • 接触(Touches)几何形状有至少一个公共的边界点,但是没有内部点。

  • 交叉(Crosses)几何形状共享一些但不是所有的内部点。

  • 内含(Within)几何形状A的线都在几何形状B内部。

  • 包含(Contains)几何形状B的线都在几何形状A内部(区别于内含)

  • 重叠(Overlaps)几何形状共享一部分但不是所有的公共点,而且相交处有他们自己相同的区域。

一些情况下,确切定义有些微妙。 你可以参考JTS技术说明书来完全确定在任何情况下的返回情况。

以上的运算返回的都是 true 或者 false

另外一种是空间叠加分析操作。主要有下面几个操作:

  • 缓冲区分析(Buffer)包含所有的点在一个指定距离内的多边形和多多边形。

  • 凸壳分析(ConvexHull)包含几何形体的所有点的最小凸壳多边形,(就是外包多边形啦)

  • 交叉分析(Intersection)交叉操作就是多边形AB中所有共同点的集合。

  • 联合分析(Union)AB的联合操作就是AB所有点的集合。

  • 差异分析(Difference)AB形状的差异分析就是A里有B里没有的所有点的集合。

  • 对称差异分析(SymDifference)AB形状的对称差异分析就是位于A中或者B中但不同时在AB中的所有点的集合

另外还支持多边形化,连接有向线段,压出节点等等操作。

5.1.4. 快速上手

用法

Shapely版本,GEOS库版本和GEOS C API版本均可通过shapely访问。version,

shapely.geos.geos_version_string, and shapely.geos.geos_capi_version.

>>> from osgeo import ogr
>>> import shapely
>>>
>>> print(shapely.__version__)
>>> import shapely.geos
>>>
>>> print(shapely.geos.geos_version)
>>> print(shapely.geos.geos_version_string)
1.6.4
(3, 7, 1)
3.7.1-CAPI-1.11.1 27a5e771

缓冲一个点

>>> from shapely.geometry import Point
>>> point = Point(10, 10)
>>> pt_buf = point.buffer(5)
>>> print(pt_buf.wkt)
POLYGON ((15 10, 14.97592363336098 9.509914298352198, 14.90392640201615 9.02454838991936, 14.78470167866104 8.54857661372769, 14.61939766255643 8.086582838174554, 14.40960632174178 7.643016315870014, 14.15734806151273 7.222148834901992, 13.86505226681369 6.828033579181776, 13.53553390593274 6.464466094067266, 13.17196642081823 6.134947733186318, 12.77785116509802 5.842651938487276, 12.35698368412999 5.590393678258227, 11.91341716182545 5.380602337443569, 11.45142338627232 5.215298321338958, 10.97545161008065 5.096073597983849, 10.49008570164781 5.024076366639017, 10.00000000000001 5, 9.509914298352205 5.024076366639015, 9.024548389919367 5.096073597983846, 8.548576613727697 5.215298321338953, 8.086582838174561 5.380602337443563, 7.643016315870021 5.59039367825822, 7.222148834901997 5.842651938487268, 6.82803357918178 6.134947733186308, 6.464466094067269 6.464466094067256, 6.134947733186321 6.828033579181766, 5.842651938487278 7.222148834901982, 5.590393678258229 7.643016315870005, 5.380602337443569 8.086582838174545, 5.215298321338958 8.548576613727683, 5.096073597983849 9.024548389919353, 5.024076366639016 9.509914298352191, 5 9.999999999999995, 5.024076366639015 10.4900857016478, 5.096073597983847 10.97545161008064, 5.215298321338954 11.45142338627231, 5.380602337443565 11.91341716182545, 5.590393678258224 12.35698368412999, 5.842651938487273 12.77785116509801, 6.134947733186315 13.17196642081823, 6.464466094067261 13.53553390593274, 6.828033579181771 13.86505226681368, 7.222148834901985 14.15734806151272, 7.643016315870007 14.40960632174177, 8.086582838174545 14.61939766255643, 8.548576613727679 14.78470167866104, 9.024548389919348 14.90392640201615, 9.509914298352184 14.97592363336098, 9.999999999999986 15, 10.49008570164779 14.97592363336098, 10.97545161008062 14.90392640201616, 11.45142338627229 14.78470167866105, 11.91341716182543 14.61939766255644, 12.35698368412997 14.40960632174179, 12.77785116509799 14.15734806151274, 13.17196642081821 13.8650522668137, 13.53553390593272 13.53553390593276, 13.86505226681367 13.17196642081825, 14.15734806151271 12.77785116509804, 14.40960632174176 12.35698368413002, 14.61939766255642 11.91341716182548, 14.78470167866103 11.45142338627235, 14.90392640201615 10.97545161008068, 14.97592363336098 10.49008570164784, 15 10.00000000000004, 15 10))

有关Shapely提供的空间操作和谓词的更多示例, 请参阅tests.txt和Predicates.txt。 或参见Point.txt,LineString?.txt等几何API的例子。

与shapfile

>>> # Generated in chapter 1.
>>> datasource = ogr.Open('/gdata/GSHHS_c.shp')
>>> layer = datasource.GetLayer(0)
>>> feature = layer.GetFeature(0)
>>> # geom = feature.GetGeometryRef()
>>> polygon = feature.GetGeometryRef()
>>>
>>> # print(polygon.ExportToWkt())
>>>
>>> buf = polygon.Buffer(2)
>>> buf.ExportToWkt()[:40]
'POLYGON ((-10.9650568460992 37.422419294'

Numpy接口

所有的Shapely几何实例都提供了Numpy数组接口:

>>> from numpy import asarray
>>> a = asarray(point)
>>> a.size
>>> a.shape
(2,)

Numpy数组也可以被调整为Shapely的点和线串(linestrings)。

>>> from shapely.geometry import asLineString
>>> a = asarray([[1.0, 2.0], [3.0, 4.0]])
>>> line = asLineString(a)
>>> print(line.wkt)
LINESTRING (1 2, 3 4)

性能

Shapely 使用 GEOS 库来完成所有的操作。 GEOS 是使用 C++ 来完成的并且得到广泛的应用, 可以相信它已经得到了高度的优化。

shapely.speedups模块包含用C编写的性能增强功能。 当Python在安装过程中访问编译器时,它们会自动安装。

您可以检查加速功能是否与可用属性一起安装。 默认情况下禁用构造函数加速。 你可以启用加速调用enable()。 你也可以使用disable()恢复默认的实现。

>>> from shapely import speedups
>>> speedups.available
>>>
True
>>> speedups.enable()