如何将自相交LineString分割为非自相交子线串？

问题背景 #

我有一组定义LineString的坐标，该LineString存在自相交情况：

coordinates = [
    [0, 3],
    [0, 5],
    [4, 5],
    [4, 0],
    [0, 0],
    [0, 5],
    [2, 5]
]

需要将其分割为更小的子线串，满足以下要求：

• 子线串数量最少
• 各子线串的坐标数量尽可能均衡

本案例的预期输出为：

line0 = [
    [0, 3],
    [0, 5],
    [4, 5],
    [4, 0]
]
line1 = [
    [4, 0],
    [0, 0],
    [0, 5],
    [2, 5]
]

我的尝试 #

目前我通过Shapely的LineString构建了相交矩阵来检测交点：

from shapely.geometry import LineString
from itertools import combinations
import numpy as np

def get_intersection_matrix(coordinates):
    linestrings = [
        (ix, LineString([c0, c1]))
        for ix, (c0, c1) in enumerate(zip(coordinates[:-1], coordinates[1:]))
    ]
    M = np.zeros((len(linestrings), len(linestrings)))
    for (ix0, ls0), (ix1, ls1) in combinations(linestrings, 2):
        if abs(ix0 - ix1) == 1: # 忽略相邻连接的线段
            continue
        if ls0.intersects(ls1):
            M[ix0, ix1], M[ix1, ix0] = 1, 1
    return M

该函数输出的“相交矩阵”如下：

>>> get_intersection_matrix(coordinates)
array([[0, 0, 0, 0, 1, 1],
       [0, 0, 0, 0, 1, 1],
       [0, 0, 0, 0, 0, 0],
       [0, 0, 0, 0, 0, 0],
       [1, 1, 0, 0, 0, 0],
       [1, 1, 0, 0, 0, 0]])

该矩阵可解读为：

• 第0、1段线段与第4、5段线段相交
• 第4、5段线段与第0、1段线段相交

另外，我认为相交簇的数量可指示子线串的数量：相交簇数量 + 1。

解决方案 #

思路分析 #

从相交矩阵能看出，第0、1段和第4、5段属于交叉冲突组，不能放在同一个子线串中；而第2、3段是无冲突的过渡线段，它的端点可以作为分割点。选择该过渡线段的终点[4,0]（坐标索引3）作为分割点，既能让两个子线串内部无自相交，又满足子线串数量最少、坐标数量均衡的要求。

代码实现 #

from shapely.geometry import LineString
from itertools import combinations
import numpy as np

def get_intersection_matrix(coordinates):
    linestrings = [
        (ix, LineString([c0, c1]))
        for ix, (c0, c1) in enumerate(zip(coordinates[:-1], coordinates[1:]))
    ]
    M = np.zeros((len(linestrings), len(linestrings)), dtype=int)
    for (ix0, ls0), (ix1, ls1) in combinations(linestrings, 2):
        if abs(ix0 - ix1) == 1:
            continue
        if ls0.intersects(ls1):
            M[ix0, ix1], M[ix1, ix0] = 1, 1
    return M

def find_split_points(coordinates):
    M = get_intersection_matrix(coordinates)
    # 找出无相交的过渡线段索引
    no_intersections = np.where(M.sum(axis=1) == 0)[0]
    
    split_points = []
    for seg_ix in no_intersections:
        split_ix = seg_ix + 1
        # 检查该分割点是否能分隔冲突组
        left_segs = set(range(seg_ix + 1))
        right_segs = set(range(seg_ix + 1, len(M)))
        cross_intersect = False
        for l in left_segs:
            for r in right_segs:
                if M[l, r] == 1:
                    cross_intersect = True
                    break
            if cross_intersect:
                break
        if not cross_intersect:
            split_points.append(split_ix)
    
    return split_points

def split_linestring(coordinates):
    split_points = find_split_points(coordinates)
    # 取第一个有效分割点，保证子线串数量最少
    if split_points:
        split_ix = split_points[0]
        return [coordinates[:split_ix + 1], coordinates[split_ix:]]
    # 无自相交时直接返回原线串
    return [coordinates]

# 测试执行
coordinates = [
    [0, 3],
    [0, 5],
    [4, 5],
    [4, 0],
    [0, 0],
    [0, 5],
    [2, 5]
]

result = split_linestring(coordinates)
for i, line in enumerate(result):
    print(f"line{i} = {line}")

运行输出：

line0 = [[0, 3], [0, 5], [4, 5], [4, 0]]
line1 = [[4, 0], [0, 0], [0, 5], [2, 5]]

内容的提问来源于stack exchange，提问作者Joost Döbken