python地理处理包——Shapely介绍及用户手册

本文主要是基于shapely官方文档 ，根据个人理解进行整理

shapely主要是在笛卡尔平面对几何对象进行操作和分析，它是一个BSD授权的Python包。Shapely不关心数据格式或坐标系，但可以很容易地与这些文件包集成。

性能

Shapely中所有的操作都是使用GEOS库。GEOS是用C++写的，也被用在许多应用程序中，你可以期待所有操作都是高度优化的。创建新的几何图形与许多坐标,然而,涉及一些开销,可能会减缓你的代码。

用法

下面是一个典型的例子，通过点缓冲来构建一个缓冲圆。

>>> from shapely.geometry import Point
>>> patch = Point(0.0, 0.0).buffer(10.0)
>>> patch
<shapely.geometry.polygon.Polygon object at 0x...>
>>> patch.area
313.65484905459385

 用户手册

• 空间数据模型

Shapely实现的几何对象的基本类型是点、曲线和曲面。每一个都与平面上的三组（可能是无限的）点相关联。要素的内部集、边界集和外部集是互斥的，它们的并集与整个平面重合，其数据模型如下：

1、点有一个正好由一个点组成的内部集，一个完全没有点的边界集，以及一个所有其他点的外部集。点的拓扑维数为0；

2、曲线有一个由沿其长度无限多个点组成的内部集（想象一个在空间中拖动的点），一个由其两个端点组成的边界集，以及一个由所有其他点组成的外部集。曲线的拓扑维数为1；

3、曲面有一个内部集，它由内部无限多个点组成（想象一条曲线在空间中拖动以覆盖一个区域），一个边界集由一条或多条曲线组成，以及所有其他点的外部集，包括可能存在于曲面中的孔内的点。曲面的拓扑维数为2；

点类型由point类实现；curve由LineString和LinearRing类实现；surface由Polygon类实现。Shapely实现不平滑（即具有连续切线）曲线。所有曲线必须用线性样条曲线近似。所有的圆形面片必须由线性样条线包围的区域来近似。

点集合由MultiPoint类实现，曲线集合由MultiLineString类实现，曲面集合由MultiPolygon类实现。这些集合在计算上并不重要，但对某些类型的特征建模非常有用。例如，Y形线要素可以很好地通过多线串作为一个整体进行建模。

• 关系

空间数据模型伴随着一组几何对象之间的自然语言关系（包含、相交、重叠、接触等），以及一个使用其组成点集3的相互交集的3x3矩阵来理解它们的理论框架：DE-9IM。这是一种拓扑模型，用于描述两个几何图形空间关系的一种标准，学过地理课程的一般都了解

• 操作

根据JTS技术规范，本手册将区分构造（缓冲区、凸壳）和集合论操作（交集、并集等）。

• 坐标系

Shapely不支持坐标系转换。对两个或多个特征的所有操作都假定这些特征存在于同一个笛卡尔平面上。

• 几何对象

几何对象是以典型的Python方式创建的，使用类本身作为实例工厂。

Point、LineString和LinearRing的实例最重要的属性是确定其内部、边界和外部点集的有限坐标序列。一个线串可以由2个点确定，但包含无限多个点。坐标序列是不变的。构造实例时可以使用第三个z坐标值，但对几何分析没有影响。所有操作都在x-y平面上进行。

在所有构造函数中，数值被转换为float类型。换句话说，点（0，0）和点（0.0，0.0）生成几何等效的实例。在构造实例时，Shapely不会检查拓扑的简单性或有效性。

注意：Shapely是一个平面几何库，在几何分析中忽略了平面上方或下方的高度z。对于用户来说，这里有一个潜在的陷阱：只有z轴不同的坐标元组之间没有区别，它们的应用可能会导致令人惊讶的无效几何对象。例如，LineString（[（0，0，0），（0，0，1）]）不返回单位长度的垂直线，而是平面中长度为零的无效直线。类似地，多边形（[（0，0，0），（0，0，1），（1，1，1）]）不是由闭合环包围的，因此是无效的。

• 一般对象和方法

object.area

返回对象面积(float) 

object.bounds

返回对象边界(minx, miny, maxx, maxy) 的元组(float) 

object.length

返回对象长度(float)

object.minimum_clearance

返回可以移动节点以生成无效几何体的最小距离。

这可以被认为是一个几何体鲁棒性的度量，其中最小间隙值越大，表示几何体越坚固。如果几何图形（例如点）不存在最小间隙，则将返回 math.infinity。该方法需要GEOS 3.6 以上版本。

>>> from shapely.geometry import Polygon
>>> Polygon([[0, 0], [1, 0], [1, 1], [0, 1], [0, 0]]).minimum_clearance
1.0

object.geom_type

　　返回对象的几何类型

>>> Point(0, 0).geom_type
'Point'

object.distance(other)

　　返回到另一个几何对象的最小距离

>>> Point(0,0).distance(Point(1,1))
1.4142135623730951

object.hausdorff_distance(other) （1.6.0版本以上）

　　返回到另一个几何对象的Hausdorff距离（float）。两个几何体之间的Hausdorff距离是任意一个几何体上的一个点与另一个几何体上的最近点之间的最远距离

>>> point = Point(1, 1)
>>> line = LineString([(2, 0), (2, 4), (3, 4)])
>>> point.hausdorff_distance(line)
3.605551275463989
>>> point.distance(Point(3, 4))
3.605551275463989

object.representative_point()

　　返回一个可以保证在几何对象内的廉价计算点

注意：这通常与质心不同。

>>> donut = Point(0, 0).buffer(2.0).difference(Point(0, 0).buffer(1.0))
>>> donut.centroid.wkt
'POINT (-0.0000000000000001 -0.0000000000000000)'
>>> donut.representative_point().wkt
'POINT (-1.5000000000000000 0.0000000000000000)'

• Points

 class Point(coordinates)

　　点构造函数采用位置坐标值或点元组参数。

>>> from shapely.geometry import Point
>>> point = Point(0.0, 0.0)
>>> q = Point((0.0, 0.0))

　　1、点的面积和长度属性为0

>>> point.area
0.0
>>> point.length
0.0

　　2、点的边界是均为0的元组

>>> point.bounds
(0.0, 0.0, 0.0, 0.0)

　　3、点的坐标值可以使用coords, x, y, 和z 获取

>>> list(point.coords)
[(0.0, 0.0)]
>>> point.x
0.0
>>> point.y
0.0

　　4、坐标可以被切片（1.2.14版本新特性）

>>> point.coords[:]
[(0.0, 0.0)]

　　5、Point构造函数还接受另一个Point实例，从而生成一个副本

>>> Point(point)
<shapely.geometry.point.Point object at 0x...>

• LineString

class LineString(coordinates)

LineString构造函数采用2个或更多（x，y[，z]）点元组的有序序列

构造的LineString对象表示点之间的一个或多个连接的线性样条曲线。允许按顺序重复点，但可能会导致性能损失，应避免。线串可以交叉（即复杂而不简单）

图1 左边是一个简单的LineString，右边是一个复杂的LineString。每个点的（MultiPoint）边界显示为黑色，描述这些线的其他点显示为灰色　

　　1、LineString的面积为0，长度不为0

>>> from shapely.geometry import LineString
>>> line = LineString([(0, 0), (1, 1)])
>>> line.area
0.0
>>> line.length
1.4142135623730951

　　2、LineString的边界是 (minx, miny, maxx, maxy) 的元组

>>> line.bounds
(0.0, 0.0, 1.0, 1.0)

　　3、LineString的坐标值可以使用coords获取

>>> len(line.coords)
2
>>> list(line.coords)
[(0.0, 0.0), (1.0, 1.0)]

　　4、坐标可以被切片（1.2.14版本新特性）

>>> point.coords[:]
[(0.0, 0.0), (1.0, 1.0)]
>>> point.coords[1:]
[(1.0, 1.0)]

　　5、LineString构造函数还接受另一个LineString实例，从而生成一个副本

>>> LineString(line)
<shapely.geometry.linestring.LineString object at 0x...>

 　　6、也可以使用一系列混合点实例或坐标元组来构造LineString。各个坐标将复制到新对象中

>>> LineString([Point(0.0, 1.0), (2.0, 3.0), Point(4.0, 5.0)])
<shapely.geometry.linestring.LineString object at 0x...>

• LineRings

class LinearRing(coordinates)

LinearRing构造函数采用（x，y[，z]）点元组的有序序列

通过在第一个和最后一个索引中传递相同的值，可以显式地闭合序列。否则，将通过将第一个元组复制到最后一个索引来隐式闭合序列。与LineString一样，允许有序序列中的重复点，但可能会导致性能损失，因此应该避免。线迹不能交叉，也不能单点接触

图2 左边是有效的LinearRing，右边是无效的自相交LinearRing。描述环的点以灰色显示。环的边界是空的

注意：Shapely不会阻止这种环的产生，但在对其进行操作时会引发异常。

　　1、LinearRing的面积为0，长度不为0

>>> from shapely.geometry.polygon import LinearRing
>>> ring = LinearRing([(0, 0), (1, 1), (1, 0)])
>>> ring.area
0.0
>>> ring.length
3.4142135623730949

　　2、LinearRing的边界是 (minx, miny, maxx, maxy) 的元组

>>> ring.bounds
(0.0, 0.0, 1.0, 1.0)

　　3、LinearRing的坐标值可以使用coords获取

>>> len(ring.coords)
4
>>> list(ring.coords)
[(0.0, 0.0), (1.0, 1.0), (1.0, 0.0), (0.0, 0.0)]

　　4、坐标可以被切片（1.2.14版本新特性）

>>> point.coords[:]
[(0.0, 0.0), (1.0, 1.0)]
>>> point.coords[1:]
[(1.0, 1.0)]

　　5、LinearRing构造函数还接受另一个LinearRing实例，从而生成一个副本

>>> LinearRing(ring)
<shapely.geometry.polygon.LinearRing object at 0x...>

• Polygons

class Polygon(shell[, holes=None])

Polygon构造函数接受两个位置参数。第一个是（x，y[，z]）点元组的有序序列，其处理方式与LinearRing完全相同。第二个是一个可选的无序的环状序列，指定了特征的内部边界或“洞”

有效多边形的环不能相互交叉，且只能接触一个点。同样，Shapely不会阻止无效特性的创建，但是在操作它们时会引发异常

 a)有效，其中一个内环在一个点上与外环接触；b)无效，因为它的内环在多个点上与外环接触；c)无效，因为它的外环和内环沿一条线接触；d)无效，因为它的内环沿着一条线接触

　1、Polygon的面积和长度不为0

>>> from shapely.geometry import Polygon
>>> polygon = Polygon([(0, 0), (1, 1), (1, 0)])
>>> polygon.area
0.5
>>> polygon.length
3.4142135623730949

　2、Polygon的边界是(minx, miny, maxx, maxy)数组

>>> from shapely.geometry import Polygon
>>> polygon = Polygon([(0, 0), (1, 1), (1, 0)])
>>> polygon.area
0.5
>>> polygon.length
3.4142135623730949

　3、Polygon的构成环通过exterior 和 interiors 属性获取

>>> list(polygon.exterior.coords)
[(0.0, 0.0), (1.0, 1.0), (1.0, 0.0), (0.0, 0.0)]
>>> list(polygon.interiors)
[]

　4、Polygon的构造接受LineString 和LinearRing实例

>>> coords = [(0, 0), (1, 1), (1, 0)]
>>> r = LinearRing(coords)
>>> s = Polygon(r)
>>> s.area
0.5
>>> t = Polygon(s.buffer(1.0).exterior, [r])
>>> t.area
6.5507620529190334

  5、矩形的构造可以使用shapely.geometry.box()函数（1.2.9新版本特性）

shapely.geometry.box(minx, miny, maxx, maxy, ccw=True) 根据提供的边界框值生成矩形，默认情况下按逆时针顺序

>>> from shapely.geometry import box
>>> b = box(0.0, 0.0, 1.0, 1.0)
>>> b
<shapely.geometry.polygon.Polygon object at 0x...>
>>> list(b.exterior.coords)
[(1.0, 0.0), (1.0, 1.0), (0.0, 1.0), (0.0, 0.0), (1.0, 0.0)]

 5、要获得具有已知方向的多边形，使用shapely.geometry.polygon.orient()函数（1.2.10新版本特性）

shapely.geometry.polygon.orient(polygon,sign=1.0)

返回给定多边形的正确定向副本。返回结果有符号区域将具有给定的符号。符号1.0表示其外环的坐标方向将为逆时针方向

• Collections

几何对象的异构集合可能是由一些Shapely操作造成的。例如，两个LineStrings 可以沿一条线在一个点相交。为了表示这些结果，Shapely提供了类似冻结集的、不可变的几何对象集合。集合可以是同质（MultiPoint 等）或异构的

>>> a = LineString([(0, 0), (1, 1), (1,2), (2,2)])
>>> b = LineString([(0, 0), (1, 1), (2,1), (2,2)])
>>> x = a.intersection(b)
>>> x
<shapely.geometry.collection.GeometryCollection object at 0x...>
>>> from pprint import pprint
>>> pprint(list(x))
[<shapely.geometry.point.Point object at 0x...>,
 <shapely.geometry.linestring.LineString object at 0x...>]

 a)一条绿色和黄色的线，沿着一条线以及一个点相交；b)交点（蓝色）是一个包含一个LineString 和一个Point的集合

  1、GeometryCollection的成员可以通过geoms属性或迭代器协议in或list()进行访问

>>> pprint(list(x.geoms))
[<shapely.geometry.point.Point object at 0x...>,
 <shapely.geometry.linestring.LineString object at 0x...>]
>>> pprint(list(x))
[<shapely.geometry.point.Point object at 0x...>,
 <shapely.geometry.linestring.LineString object at 0x...>]

  2、集合也可以切片（1.2.14新版本特性）

>>> from shapely.geometry import MultiPoint
>>> m = MultiPoint([(0, 0), (1, 1), (1,2), (2,2)])
>>> m[:1].wkt
'MULTIPOINT (0.0000000000000000 0.0000000000000000)'
>>> m[3:].wkt
'MULTIPOINT (2.0000000000000000 2.0000000000000000)'
>>> m[4:].wkt
'GEOMETRYCOLLECTION EMPTY'

注意：如果可能的话，最好使用下面描述的同类集合类型之一。

• Collections of Points

class MultiPoint(points)

MultiPoint 构造函数接受一系列（x，y[，z]）点元组

  1、MultiPoint 面积和长度为0

>>> from shapely.geometry import MultiPoint
>>> points = MultiPoint([(0.0, 0.0), (1.0, 1.0)])
>>> points.area
0.0
>>> points.length
0.0

  2、其 x-y 边界范围是(minx, miny, maxx, maxy) 元组

>>> points.bounds
(0.0, 0.0, 1.0, 1.0)

  3、GeometryCollection的成员可以通过geoms属性或迭代器协议in或list()进行访问

>>> import pprint
>>> pprint.pprint(list(points.geoms))
[<shapely.geometry.point.Point object at 0x...>,
 <shapely.geometry.point.Point object at 0x...>]
>>> pprint.pprint(list(points))
[<shapely.geometry.point.Point object at 0x...>,
 <shapely.geometry.point.Point object at 0x...>]

  4、构造函数还接受另一个MultiPoint 实例或无序的点实例序列，从而生成副本

>>> MultiPoint([Point(0, 0), Point(1, 1)])
<shapely.geometry.multipoint.MultiPoint object at 0x...>

• Collections of Points

class MultiLineString(lines)

MultiLineString构造函数接受一系列类似于直线的序列或对象

 a)不连接的MultiLineString；b)复杂的 MultiLineString  

  1、MultiLineString 的面积和长度为0

>>> from shapely.geometry import MultiLineString
>>> coords = [((0, 0), (1, 1)), ((-1, 0), (1, 0))]
>>> lines = MultiLineString(coords)
>>> lines.area
0.0
>>> lines.length
3.4142135623730949

  2、其 x-y 边界范围是(minx, miny, maxx, maxy) 元组

>>> lines.bounds
(-1.0, 0.0, 1.0, 1.0)

  3、它的成员是LineString的实例，可以通过geoms属性或迭代器协议in或list()进行访问

>>> len(lines.geoms)
2
>>> pprint.pprint(list(lines.geoms))
[<shapely.geometry.linestring.LineString object at 0x...>,
 <shapely.geometry.linestring.LineString object at 0x...>]
>>> pprint.pprint(list(lines))
[<shapely.geometry.linestring.LineString object at 0x...>,
 <shapely.geometry.linestring.LineString object at 0x...>]

  4、构造函数还接受多行字符串的另一个实例或无序的LineString实例序列，从而生成副本

>>> MultiLineString(lines)
<shapely.geometry.multilinestring.MultiLineString object at 0x...>
>>> MultiLineString(lines.geoms)
<shapely.geometry.multilinestring.MultiLineString object at 0x...>

• Collections of Polygons

class MultiPolygon(polygons)

MultiPolygon构造函数采用一系列外部环和孔列表元组：[（（a1，…，aM），[（b1，…，bN），…]），…]

 

 a)有效的有两个成员的MultiPolygon；b)无效MultiPolygon，因为它的成员接触到无限多个点（沿着一条线）

  1、其 x-y 边界范围是(minx, miny, maxx, maxy) 元组

>>> polygons.bounds
(-1.0, -1.0, 2.0, 2.0)

  2、它的成员是Polygon的实例，可以通过geoms属性或迭代器协议in或list()进行访问

>>> len(polygons.geoms)
3
>>> len(polygons)
3

• Empty features

  1、它是一个点集与空集重合；不是没有，而是类似于集（[]）。可以通过不带参数地调用各种构造函数来创建空要素。空要素几乎不支持任何操作。

>>> line = LineString()
>>> line.is_empty
True
>>> line.length
0.0
>>> line.bounds
()
>>> line.coords
[]

  2、可以设置空要素的坐标，之后几何图形不再为空

>>> line.coords = [(0, 0), (1, 1)]
>>> line.is_empty
False
>>> line.length
1.4142135623730951
>>> line.bounds
(0.0, 0.0, 1.0, 1.0)

• Coordinate sequences

  1、描述geometry的坐标序列表示为CoordinateSequence对象，这些序列不应该直接序列化，但可以从现有的geometry的Geometry.coords属性访问

>>> line = LineString([(0, 1), (2, 3), (4, 5)])
>>> line.coords
<shapely.coords.CoordinateSequence object at 0x00000276EED1C7F0>

  2、坐标序列可以被索引、切片和迭代，就像它们是一个坐标元组的列表一样。

>>> line.coords[0]
(0.0, 1.0)
>>> line.coords[1:]
[(2.0, 3.0), (4.0, 5.0)]
>>> for x, y in line.coords:
...     print("x={}, y={}".format(x, y))
...
x=0.0, y=1.0
x=2.0, y=3.0
x=4.0, y=5.0

  3、多边形的外部和每个内环都有一个坐标序列

>>> poly = Polygon([(0, 0), (0, 1), (1, 1), (0, 0)])
>>> poly.exterior.coords
<shapely.coords.CoordinateSequence object at 0x00000276EED1C048>

  4、多部分几何图形没有坐标序列。其坐标序列存储在组成的几何图形上

>>> p = MultiPoint([(0, 0), (1, 1), (2, 2)])
>>> p[2].coords
<shapely.coords.CoordinateSequence object at 0x00000276EFB9B320>

• 线性参考方法

使用一维参照系统指定沿线性要素（例如LineStrings 和MultiLineStrings ）的位置非常有用。Shapely支持基于长度或距离的线性参考，计算沿几何对象到给定点的投影的距离，或沿对象给定距离的点的距离（需GEOS3.2.0及以上版本支持）

object.interpolate(distance[, normalized=False])

       返回沿线性几何对象指定距离的点

如果normalized  参数为True，距离将被解释为几何对象长度的比例部分

>>> ip = LineString([(0, 0), (0, 1), (1, 1)]).interpolate(1.5)
>>> ip
<shapely.geometry.point.Point object at 0x740570>
>>> ip.wkt
'POINT (0.5000000000000000 1.0000000000000000)'
>>> LineString([(0, 0), (0, 1), (1, 1)]).interpolate(0.75, normalized=True).wkt
'POINT (0.5000000000000000 1.0000000000000000)'

object.project(other[, normalized=False])

       返回沿此几何对象到另一个对象最近的点的距离。

 如果normalized  参数为True，则返回对象长度的标准化距离，project()方法是interpolate()的逆方法。

 >>> LineString([(0, 0), (0, 1), (1, 1)]).project(ip)

1.5
>>> LineString([(0, 0), (0, 1), (1, 1)]).project(ip, normalized=True)
0.75

例如，可以使用线性参照方法在指定距离处剪切直线

def cut(line, distance):
    # Cuts a line in two at a distance from its starting point
    if distance <= 0.0 or distance >= line.length:
        return [LineString(line)]
    coords = list(line.coords)
    for i, p in enumerate(coords):
        pd = line.project(Point(p))
        if pd == distance:
            return [
                LineString(coords[:i+1]),
                LineString(coords[i:])]
        if pd > distance:
            cp = line.interpolate(distance)
            return [
                LineString(coords[:i] + [(cp.x, cp.y)]),
                LineString([(cp.x, cp.y)] + coords[i:])]

>>> line = LineString([(0, 0), (1, 0), (2, 0), (3, 0), (4, 0), (5, 0)])
>>> pprint([list(x.coords) for x in cut(line, 1.0)])
[[(0.0, 0.0), (1.0, 0.0)],
 [(1.0, 0.0), (2.0, 0.0), (3.0, 0.0), (4.0, 0.0), (5.0, 0.0)]]
>>> pprint([list(x.coords) for x in cut(line, 2.5)])
[[(0.0, 0.0), (1.0, 0.0), (2.0, 0.0), (2.5, 0.0)],
 [(2.5, 0.0), (3.0, 0.0), (4.0, 0.0), (5.0, 0.0)]]

• 谓词和关系

几何对象中解释的类型的对象提供标准谓词作为属性（对于一元谓词）和方法（对于二元谓词）。无论是一元还是二元谓词，都返回True或False

object.interpolate(distan

• 一元谓词

  标准的一元谓词被实现为只读属性属性

object.has_z

      如果特征不仅具有的x和y坐标，还具有三维（或所谓的2.5D）几何图形的z坐标，则返回True

>>> Point(0, 0).has_z
False
>>> Point(0, 0, 0).has_z
True

 object.is_ccw （1.2.10版本）

      如果坐标按逆时针顺序（以正向符号包围区域）返回True。此方法仅适用于线性化对象

>>> LinearRing([(1,0), (1,1), (0,0)]).is_ccw
True

环可以通过以下方式实现反转：

>>> ring = LinearRing([(0,0), (1,1), (1,0)])
>>> ring.is_ccw
False
>>> ring.coords = list(ring.coords)[::-1]
>>> ring.is_ccw
True

object.is_empty

      如果要素的内部和边界（在点集术语中）与空集重合，则返回True

>>> Point().is_empty
True
>>> Point(0, 0).is_empty
False

注意：在操作符模块的attrgetter()函数的帮助下，诸如is_empty之类的一元谓词可以很容易地用作内置filter()或itertools.ifilter()

>>> from operator import attrgetter
>>> empties = filter(attrgetter('is_empty'), [Point(), Point(0, 0)])
>>> len(empties)
1

 object.is_ring

      如果要素是闭合的简单LineString，则返回True。一个封闭要素的边界与空集合重合。

>>> LineString([(0, 0), (1, 1), (1, -1)]).is_ring
False
>>> LinearRing([(0, 0), (1, 1), (1, -1)]).is_ring
True

此属性适用于LineString和LinearRing实例，但对其他实例没有意义

 object.is_simple

      如果要素不自相交，则返回True （只对LineStrings 和LinearRings有意义）

>>> LineString([(0, 0), (1, 1), (1, -1), (0, 1)]).is_simple
False

object.is_valid

      如果要素是“有效的”，则返回True

一个有效的LinearRing 不能在一个点上与自己交叉或相接。一个有效的Polygon 不能有任何重叠的外环或内环。有效的MultiPolygon 不能集合任何重叠的多边形。对无效要素的操作可能会失败

>>> MultiPolygon([Point(0, 0).buffer(2.0), Point(1, 1).buffer(2.0)]).is_valid
False

上面的两个点非常接近，缓冲区操作产生的多边形（在下一节中解释）会重叠

注意：is_valid谓词可用于编写验证装饰器，该装饰器可确保从构造函数函数只返回有效对象

from functools import wraps
def validate(func):
    @wraps(func)
    def wrapper(*args, **kwargs):
        ob = func(*args, **kwargs)
        if not ob.is_valid:
            raise TopologicalError(
                "Given arguments do not determine a valid geometric object")
        return ob
    return wrapper

>>> @validate
... def ring(coordinates):
...     return LinearRing(coordinates)
...
>>> coords = [(0, 0), (1, 1), (1, -1), (0, 1)]
>>> ring(coords)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "<stdin>", line 7, in wrapper
shapely.geos.TopologicalError: Given arguments do not determine a valid geometric object