Python线要素与Shapefile裁剪求交性能优化方案咨询

优化方案 #

核心性能瓶颈原因 #

你的代码慢主要来自三个问题：

• 同一条线与海岸线的intersection重复计算了3次，每次计算都需要遍历复杂海岸线的所有顶点
• 没有用到空间索引，每条线的相交判断都是全量计算
• 纯Python for循环遍历360条线，没有用上向量化操作加速

具体优化措施 #

1. 缓存重复计算结果：单次计算相交结果后存为变量，重复调用即可
2. 引入空间索引加速相交判断：用shapely的STRtree封装海岸线要素，相交判断时间复杂度从O(n)降到O(logn)
3. 简化最近点计算逻辑：无需构建字典存储距离和点，直接用min函数指定距离为排序key，也可以用平方距离代替实际距离省略开根号运算，排序结果完全一致
4. 替换过时接口：cascaded_union已被弃用，改用性能更高的unary_union
5. 移除冗余代码：生成射线的chain和np.mod属于多余操作，np.arange(0,360,1)本身就生成0-359的角度序列
6. 可选：用geopandas向量化操作替代for循环：如果数据量更大，向量化操作可以进一步提升数倍性能

优化后代码示例 #

from shapely.geometry import LineString, shape, Point
from shapely.affinity import rotate
from shapely.ops import unary_union
from shapely.strtree import STRtree
from fiona import open
import numpy as np

# 读取岸线数据并构建空间索引，用with上下文自动关闭文件句柄
with open('land_boundary.shp') as shoreline:
    shapes = [shape(f['geometry']) for f in shoreline]
mergedshorelines = unary_union(shapes)
# 构建空间索引加速相交判断
shore_tree = STRtree([mergedshorelines])

# 生成基础射线
x,y = 696346,6601295
x_end, y_end = 746345,6601295
line1 = LineString([(x, y), (x_end, y_end)])
radii = [rotate(line1, deg, (x,y)) for deg in np.arange(0,360,1)]
mergedradii = unary_union(radii)

point = Point(x,y)
points_ok = []

for line in mergedradii.geoms:
    # 先用空间索引快速筛掉不相交的线
    if not shore_tree.query(line):
        continue
    # 仅计算一次相交结果
    intersect = line.intersection(mergedshorelines)
    if intersect.is_empty:
        continue
    if intersect.geom_type == "MultiPoint":
        # 用平方距离计算最近点，省略开根号运算
        nearest_pt = min(intersect.geoms, key=lambda p: (p.x - x)**2 + (p.y - y)**2)
        points_ok.append(nearest_pt)
    elif intersect.geom_type == "Point":
        points_ok.append(intersect)

Shoreline_points = unary_union(points_ok)

额外性能提升建议 #

如果岸线数据精度极高、顶点数极多，可以先对岸线做简化处理（用shapely的simplify方法），在允许的误差范围内减少顶点数，相交计算速度会有非常明显的提升。

内容的提问来源于stack exchange，提问作者björn almström