使用 python 查找线和形状文件之间的交集

Question

我正在尝试根据 python 中的形状文件剪裁线条。我有一个有效的脚本，但速度很慢。有没有更快的方法来做到这一点？先决条件是我必须使用 python.

from shapely.geometry import LineString, shape, Point
from shapely.affinity import rotate
from itertools import chain
from shapely.ops import cascaded_union
from fiona import open
import numpy as np

# open the shp-file with land geometry
shoreline = open('land_boundary.shp')
shapes = [shape(f['geometry']) for f in shoreline]
mergedshorelines = cascaded_union(shapes)

# create an arbitrary line
x,y = 696346,6601295
x_end, y_end = 746345,6601295
line1 = LineString([(x, y), (x_end, y_end)])

# Creates lines from arbitrary line with 1 degree step
radii = [rotate(line1, deg, (x,y)) for deg in chain(np.mod(np.arange(0,360,1),360))]
mergedradii = cascaded_union(radii)

# the intersection between the two multilines is computed and the intersection point with the smallest distance is choosen
point =  Point(x,y)
points_ok = []

#----THIS IS THE SLOW PART-------
# clip lines with land boundary 
for line in mergedradii:
    if line.intersects(mergedshorelines):
        if line.intersection(mergedshorelines).type == "MultiPoint":
          # multiple points: select the point with the minimum distance
          a = {}
          for pt in line.intersection(mergedshorelines):
              a[point.distance(pt)] = pt
          points_ok.append(a[min(a.keys())])
        if line.intersection(mergedshorelines).type == "Point":
          # ok, only one intersection
          points_ok.append(line.intersection(mergedshorelines))
Shoreline_points = cascaded_union(points_ok) # coordinates of all intersection.

感谢任何意见！干杯! /比约恩

Answer 1

只需调用 intersection 一次并保存结果，而不是多次调用，可以节省大量时间。我在循环的顶部调用一次 intersection 方法并将其保存为变量，然后在逻辑中引用该变量，而不是再次运行它。

您还可以通过完全跳过 intersects 检查来节省时间，因为如果一条线不与合并的海岸线相交，它就不会是多点或点对象，而是一个 LineString。

您还可以优化脚本中寻找 MultiPoint 中最近点的部分。当试图找到离直线原点最近的点时，您只需要检查 MultiPoint 中的第一个和最后一个点，因为这些点是经过排序的。所以第一个或最后一个点将最接近原点。

for line in mergedradii:
    intersection = line.intersection(mergedshorelines)
    if intersection.type == "MultiPoint":
        # multiple points: select the point with the minimum distance
        first_pt = intersection[0]
        last_pt = intersection[-1]
        if point.distance(first_pt) < point.distance(last_pt):
            points_ok.append(first_pt)
        else:
            points_ok.append(last_pt)
    elif intersection.type == "Point":
        # ok, only one intersection
        points_ok.append(intersection)
Shoreline_points = cascaded_union(points_ok) # coordinates of all intersection.

在我的电脑上，这个脚本运行大约需要 40 秒，而原始脚本大约需要 190 秒。

使用 python 查找线和形状文件之间的交集

Find intersection between line and shape file with python

python

shapely