Shapely

项目地址： https://github.com/shapely/shapely

官方文档： https://shapely.readthedocs.io/

PyPI： https://pypi.org/project/shapely/

许可证： BSD-3-Clause

概述

Shapely 是 Python 计算几何的核心库，基于 GEOS（Geometry Engine - Open Source）的 C++ 几何引擎。它为 PostGIS、GeoPandas、PyQGIS 等众多开源 GIS 项目提供几何计算底层支持。核心能力：

• 几何对象模型：Point、LineString、Polygon、Multi*、GeometryCollection
• 空间关系判断：contains、intersects、within、touches、covers、disjoint、relate（DE-9IM）
• 集合运算：intersection、union、difference、symmetric_difference
• 构造操作：buffer、simplify、convex_hull、delaunay_triangles、voronoi_polygons、make_valid
• 几何度量：area、length、distance、hausdorff_distance、frechet_distance
• 空间索引：STRtree，加速批量空间查询
• 格式读写：WKT、WKB、GeoJSON
• 双重 API：函数式（向量化，支持 NumPy 广播）+ OOP（标量便捷）

环境要求： Python 3.9+，Shapely 2.0+

Geometry Types

Z coordinates are ignored for all spatial analysis — operations are performed in the x-y plane only.

Dual API Pattern

Shapely provides two APIs. Prefer the function-based API for arrays and performance; use the OOP API for scalar convenience.

Function-based (vectorized, NumPy ufunc)

All functions release the GIL during GEOS execution, support NumPy broadcasting, and handle arrays natively.

import shapely
import numpy as np

geoms = np.array([shapely.Point(0, 0), shapely.Point(1, 1)])
shapely.area(geoms)                    # array of floats
shapely.contains(polygon, geoms)       # array of bools
shapely.buffer(geoms, 1.0)            # array of polygons
shapely.distance(geoms, other)        # array of distances

OOP (method-based, scalar)

from shapely import Point, Polygon

p = Point(0, 0)
poly = Polygon([(0, 0), (1, 0), (1, 1), (0, 1)])
p.buffer(10)
p.distance(Point(1, 1))
poly.area                  # property
poly.contains(p)
poly.intersection(other)

Operator overloading

geom1 & geom2    # intersection
geom1 | geom2    # union
geom1 - geom2    # difference
geom1 ^ geom2    # symmetric_difference

Geometry Creation

From coordinates (OOP)

from shapely import Point, LineString, LinearRing, Polygon
from shapely import MultiPoint, MultiLineString, MultiPolygon, GeometryCollection

Point(0, 0)
Point(0, 0, 5)                                    # 3D
LineString([(0, 0), (1, 1), (2, 0)])
LinearRing([(0, 0), (1, 0), (1, 1), (0, 0)])      # must be closed
Polygon([(0, 0), (1, 0), (1, 1), (0, 0)])         # auto-closes shell
Polygon(shell_coords, [hole1_coords, hole2_coords])
MultiPoint([(0, 0), (1, 1)])
MultiLineString([[(0, 0), (1, 1)], [(2, 2), (3, 3)]])
MultiPolygon([poly1, poly2])
GeometryCollection([point, line, poly])

From coordinates (vectorized)

import shapely

shapely.points([[0, 1], [2, 3]])                  # array of Points
shapely.points(0, 1)                               # scalar Point
shapely.linestrings([[[0, 0], [1, 1]], [[2, 2], [3, 3]]])
shapely.linearrings([[0, 0], [1, 0], [1, 1]])     # auto-closes
shapely.polygons(shell, holes=[hole1, hole2])
shapely.box(xmin, ymin, xmax, ymax)               # rectangle(s)
shapely.multipoints([pt1, pt2])
shapely.multilinestrings(...)
shapely.multipolygons(...)
shapely.geometrycollections(...)

From serialization formats

# WKT
shapely.from_wkt("POINT (0 0)")
shapely.to_wkt(geom, rounding_precision=6, trim=True)

# WKB
shapely.from_wkb(b'\x01\x01...')
shapely.to_wkb(geom, hex=False, flavor="extended")

# GeoJSON
shapely.from_geojson('{"type":"Point","coordinates":[0,0]}')
shapely.to_geojson(geom, indent=None)

# GeoArrow ragged arrays
shapely.from_ragged_array(geometry_type, coords, offsets)
shapely.to_ragged_array(geometries)

From/to GeoJSON-like dict (__geo_interface__)

from shapely.geometry import shape, mapping

geom = shape({"type": "Point", "coordinates": [0, 0]})
d = mapping(geom)  # {"type": "Point", "coordinates": (0.0, 0.0)}

构造操作、空间谓词、集合运算、测量、空间索引、仿射与坐标变换及线性参考的完整 API 见 reference/operations-reference.md

Validation

shapely.is_valid(geom)
shapely.is_valid_reason(geom)
shapely.make_valid(geom, method="linework")   # or "structure"

Prepared Geometry

Caches spatial index for repeated predicate tests against the same geometry.

shapely.prepare(geom)            # in-place, returns True if newly prepared
shapely.destroy_prepared(geom)   # free cache
shapely.is_prepared(geom)

Missing Values

None represents missing geometries. Consistent behavior:

• Predicates return False
• Measurements return nan
• Constructive operations return None

Use shapely.is_missing(obj) to check, distinct from shapely.is_empty(geom).

Important Caveats

1. Z ignored in analysis: All spatial operations work in x-y plane only. Geometries with different z values may still intersect or be equal.

2. contains excludes boundary: contains(line, endpoint) is False. Use covers() or intersects() if you need boundary inclusion.

3. set_coordinates modifies in-place: Copy geometry first with .copy() if originals must be preserved.

4. clip_by_rect may produce invalid output: The fast clip is not guaranteed to yield valid topology.

5. WKB drops LinearRing: LinearRings become LineStrings during WKB serialization.

6. to_geojson drops Z: Third dimension is silently discarded; LinearRing outputs as null.

7. OOP vs function buffer defaults differ: Point(0,0).buffer(1) uses quad_segs=16; shapely.buffer(Point(0,0), 1) uses quad_segs=8.

8. set_precision changes vertex order: Returned geometry is in "mild canonical form" — vertex order should not be relied upon.

9. NaN coordinates: Creating geometries with NaN/Inf is allowed by default. Use handle_nan='error' or handle_nan='skip' in creation functions to control this.

10. Prepared state not preserved: After any operation producing a new geometry, prepare() must be called again.

11. GEOS version gates: concave_hull, constrained_delaunay_triangles, remove_repeated_points require GEOS >= 3.11. coverage_simplify, has_m, get_m require GEOS >= 3.12.

Common Patterns

Point-in-polygon test for many points

from shapely import STRtree, points

pts = points(np.random.rand(10000, 2))
tree = STRtree(pts)
idx = tree.query(polygon, predicate="intersects")
inside_points = pts[idx]

Batch union of many polygons

result = shapely.union_all(polygon_array)

Coordinate projection (e.g., with pyproj)

from pyproj import Transformer
transformer = Transformer.from_crs("EPSG:4326", "EPSG:3857", always_xy=True)

def project(coords):
    x, y = transformer.transform(coords[:, 0], coords[:, 1])
    return np.column_stack([x, y])

projected = shapely.transform(geom, project)

Extract all coordinates from mixed geometries

coords = shapely.get_coordinates(geom_array, include_z=False)
# coords shape: (N, 2)

Spatial join pattern

tree = STRtree(right_geoms)
left_idx, right_idx = tree.query(left_geoms, predicate="intersects")
# left_idx[i] and right_idx[i] are paired matches

AI 使用建议

推荐工作流

1. 创建几何：使用 shapely.points() / shapely.linestrings() / shapely.polygons() 向量化创建，或 Point() / LineString() / Polygon() OOP 方式
2. 空间运算：使用向量化函数（shapely.intersection() / shapely.buffer()）批量处理，比逐个 OOP 调用快 10-100 倍
3. 空间谓词：使用 shapely.contains() / shapely.intersects() 等函数式 API，返回 bool 数组
4. 空间索引：STRtree 查询 + predicate="intersects" 加速批量空间关系判断
5. 序列化：shapely.from_wkt() / shapely.to_geojson() 等进行格式互转

关键注意事项

• Z 坐标在分析中被忽略：所有空间运算仅在 x-y 平面进行
• contains 不包含边界：如需边界判断使用 covers() 或 intersects()
• 向量化函数 vs OOP 方法：函数式 API（shapely.area(geoms)）释放 GIL，支持 NumPy 广播，性能远优于 OOP（geom.area）
• set_precision 改变顶点顺序：返回的几何为"温和规范形式"，不应依赖顶点顺序
• Prepared 状态不持久：任何操作产生新几何后需重新 prepare()
• WKB 会丢弃 LinearRing：WKB 序列化时 LinearRing 变为 LineString

常见问题

相关技能

• geopandas — 基于 Shapely 的矢量数据分析库：../geopandas/SKILL.md
• pyqgis — QGIS Python 绑定（也使用 GEOS/QgsGeometry）：../pyqgis/SKILL.md
• jts — Java 几何计算（Shapely 的 GEOS 底层是 JTS 的 C++ 移植）：../jts/SKILL.md
• gdal — 命令行数据处理：../gdal/SKILL.md

Package Structure

shapely/
├── __init__.py              # Top-level re-exports
├── geometry/                # OOP geometry classes
│   ├── base.py              # BaseGeometry (all methods/properties)
│   ├── point.py, linestring.py, polygon.py
│   ├── multipoint.py, multilinestring.py, multipolygon.py
│   ├── collection.py        # GeometryCollection
│   └── geo.py               # shape(), mapping(), box()
├── constructive.py          # buffer, simplify, hull, etc.
├── predicates.py            # contains, intersects, etc.
├── set_operations.py        # union, intersection, difference, etc.
├── measurement.py           # area, length, distance, bounds
├── creation.py              # points(), linestrings(), polygons(), box()
├── coordinates.py           # transform, get/set_coordinates
├── io.py                    # from/to WKT, WKB, GeoJSON
├── strtree.py               # STRtree spatial index
├── linear.py                # line_interpolate_point, line_merge, etc.
├── affinity.py              # rotate, scale, skew, translate
├── validation.py            # make_valid
├── prepared.py              # PreparedGeometry
├── ops.py                   # Legacy: split, nearest_points, etc.
├── algorithms/              # polylabel, oriented_envelope fallback
└── plotting.py              # Matplotlib helpers (experimental)

参考资源

• GitHub 仓库： https://github.com/shapely/shapely
• 官方文档： https://shapely.readthedocs.io/
• API 参考： https://shapely.readthedocs.io/en/stable/reference.html
• PyPI： https://pypi.org/project/shapely/
• GEOS（底层 C++ 引擎）： https://libgeos.org/
• 上游中文教程： https://znlgis.github.io/gis/tutorial/shapely/