shapely

name	shapely
description	Use when performing computational geometry in Python — intersection, union, buffer, convex hull, simplification. Shapely: Python bindings for GEOS, the C++ geometry engine that powers PostGIS.
tags	["python","geometry","geos","wkt","wkb","geojson","spatial","numpy"]

项目地址： 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

Type	Class	Description
Point	shapely.Point	Single coordinate (x, y [, z [, m]])
LineString	shapely.LineString	Ordered sequence of 2+ points
LinearRing	shapely.LinearRing	Closed, simple LineString
Polygon	shapely.Polygon	Exterior ring + optional holes
MultiPoint	shapely.MultiPoint	Collection of Points
MultiLineString	shapely.MultiLineString	Collection of LineStrings
MultiPolygon	shapely.MultiPolygon	Collection of Polygons
GeometryCollection	shapely.GeometryCollection	Heterogeneous collection

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)}

Constructive Operations

All available as shapely.<function>(geometry, ...).

Function	Description
buffer(geom, distance, quad_segs=8, cap_style="round", join_style="round", mitre_limit=5.0, single_sided=False)	Minkowski sum/difference with circle
offset_curve(geom, distance, ...)	Parallel curve at given distance
centroid(geom)	Geometric center
point_on_surface(geom)	Point guaranteed inside geometry
boundary(geom)	Topological boundary
convex_hull(geom)	Minimum convex enclosure
concave_hull(geom, ratio)	Concave enclosure (GEOS >= 3.11)
envelope(geom)	Axis-aligned bounding box
oriented_envelope(geom)	Minimum-area rotated rectangle
minimum_bounding_circle(geom)	Smallest enclosing circle
simplify(geom, tolerance, preserve_topology=True)	Douglas-Peucker simplification
snap(geom, reference, tolerance)	Snap vertices to reference
segmentize(geom, max_segment_length)	Densify by adding vertices
reverse(geom)	Reverse coordinate order
normalize(geom)	Canonical coordinate ordering
make_valid(geom, method="linework")	Repair invalid geometry
build_area(geom)	Polygonize from linework
clip_by_rect(geom, xmin, ymin, xmax, ymax)	Fast rectangle clip (may produce invalid output)
delaunay_triangles(geom, tolerance=0.0, only_edges=False)	Delaunay triangulation
constrained_delaunay_triangles(geom)	Constrained Delaunay (GEOS >= 3.11)
voronoi_polygons(geom, tolerance=0.0, extend_to=None)

Spatial Predicates

All return bool (or bool array). Available as shapely.<function>(a, b) or a.<method>(b).

Function	Description
contains(a, b)	B entirely inside A (boundary excluded)
contains_properly(a, b)	B inside A with no shared boundary
contains_xy(geom, x, y)	Fast point-in-geometry test
covered_by(a, b)	No point of A outside B
covers(a, b)	No point of B outside A
crosses(a, b)	A and B spatially cross
disjoint(a, b)	No shared space
dwithin(a, b, distance)	Within given distance
equals(a, b)	Topologically equal
equals_exact(a, b, tolerance)	Structurally equal within tolerance
intersects(a, b)	Share any space
intersects_xy(geom, x, y)	Fast point-intersects test
overlaps(a, b)	Partial overlap, same dimension
touches(a, b)	Only boundaries touch
within(a, b)	A entirely inside B
relate(a, b)	DE-9IM intersection matrix string
relate_pattern(a, b, pattern)	Test DE-9IM pattern match
is_valid(geom)	Well-formed geometry
is_valid_reason(geom)	Description if invalid
is_empty(geom)	Empty geometry test
is_simple(geom)	No self-intersections
is_ring(geom)	Closed and simple
is_closed(geom)	First == last point
is_ccw(geom)

Set-Theoretic Operations

Function	Description
intersection(a, b, grid_size=None)	Shared geometry
intersection_all(geometries)	N-way intersection
union(a, b, grid_size=None)	Merge two geometries
union_all(geometries)	N-way union
difference(a, b, grid_size=None)	A minus B
symmetric_difference(a, b, grid_size=None)	XOR of A and B
coverage_union(a, b)	Optimized union for non-overlapping polygons
coverage_union_all(geometries)	N-way coverage union

All binary operations support optional grid_size for precision snapping.

Measurements

Function	Description
area(geom)	Area of polygon(s)
length(geom)	Length of line / polygon perimeter
distance(a, b)	Minimum Cartesian distance
hausdorff_distance(a, b, densify=None)	Hausdorff distance
frechet_distance(a, b, densify=None)	Frechet distance
bounds(geom)	(minx, miny, maxx, maxy) per geometry
total_bounds(geometries)	Combined bounds of all geometries
minimum_clearance(geom)	Smallest move to invalidate
minimum_bounding_radius(geom)	Radius of minimum bounding circle

STRtree Spatial Index

from shapely import STRtree

tree = STRtree(geometries, node_capacity=10)

# Bounding-box query — returns indices
idx = tree.query(geometry)                            # shape (n,)
idx = tree.query(geom_array)                          # shape (2, n)

# With predicate filter
idx = tree.query(geom, predicate="intersects")
# predicates: intersects, within, contains, overlaps, crosses,
#             touches, covers, covered_by, contains_properly, dwithin

# Distance-based query
idx = tree.query(geom, predicate="dwithin", distance=5.0)

# Nearest neighbor
idx = tree.nearest(geometry)                          # single nearest

# All nearest with options
result = tree.query_nearest(
    geometry,
    max_distance=10,
    return_distance=True,   # returns (indices, distances)
    exclusive=False,
    all_matches=True
)

Geometry Accessors

shapely.get_type_id(geom)                  # GeometryType enum int
shapely.get_dimensions(geom)               # 0=point, 1=line, 2=polygon
shapely.get_coordinate_dimension(geom)     # 2, 3, or 4
shapely.get_x(point)                       # x coordinate
shapely.get_y(point)                       # y coordinate
shapely.get_z(point)                       # z coordinate
shapely.get_num_points(linestring)
shapely.get_point(linestring, index)       # extract nth point
shapely.get_exterior_ring(polygon)
shapely.get_interior_ring(polygon, index)
shapely.get_num_interior_rings(polygon)
shapely.get_geometry(collection, index)
shapely.get_num_geometries(collection)
shapely.get_parts(geom, return_index=False)
shapely.get_rings(polygon, return_index=False)
shapely.get_coordinates(geom, include_z=False, return_index=False)
shapely.set_coordinates(geom, new_coords)  # modifies in-place!
shapely.count_coordinates(geom)
shapely.get_srid(geom) / shapely.set_srid(geom, srid)
shapely.get_precision(geom) / shapely.set_precision(geom, grid_size)
shapely.force_2d(geom) / shapely.force_3d(geom, z=0)

Affine Transforms

from shapely.affinity import affine_transform, rotate, scale, skew, translate

affine_transform(geom, [a, b, d, e, xoff, yoff])   # 2D 6-element matrix
affine_transform(geom, matrix_12)                    # 3D 12-element matrix
rotate(geom, angle, origin="center")                 # degrees CCW
scale(geom, xfact=1.0, yfact=1.0, zfact=1.0, origin="center")
skew(geom, xs=0, ys=0, origin="center")              # degrees
translate(geom, xoff=0.0, yoff=0.0, zoff=0.0)

Coordinate Transforms

# Vectorized transform (fixed coordinate count)
shapely.transform(geom, lambda coords: coords * 2)
shapely.transform(geom, func, include_z=True)

# Per-ring transform (may change coordinate count)
shapely.transform_coordseq(geom, func)

Linear Referencing

shapely.line_interpolate_point(line, distance, normalized=False)
shapely.line_locate_point(line, point, normalized=False)
shapely.line_merge(multilinestring, directed=False)
shapely.shared_paths(a, b)
shapely.shortest_line(a, b)

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

---

常见问题

问题	解决
面积/距离为 0 或极小	CRS 为地理坐标系（度），需要投影到投影坐标系再计算
contains 对边界上的点返回 False	使用 covers() 或 intersects() 代替
几何操作产生 TopologyException	几何无效，用 make_valid() 修复后再操作
OOP 和函数式 API 结果不一致	buffer 默认 quad_segs 不同（OOP=16，函数式=8），显式指定即可
None 参与运算报错	使用 shapely.is_missing() 先检查，或过滤掉 None 再运算
STRtree 查询返回空数组	检查几何是否 None/Empty；确认查询几何与索引几何的 CRS 一致
WKB 序列化丢失 LinearRing	WKB 规范不支持 LinearRing，会自动转为 LineString
to_geojson 丢失 Z 坐标	GeoJSON 标准不支持 Z，改用 WKB 或自定义序列化

---

相关技能

• 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/