在与其他LineString的交点处匀称地拆分LineString

AJG*_*519 5 python shapely geopandas

我有一组LineString,这些LineString与其他LineString相交,并且我想在这些交点将LineString分成单独的段。我有一个解决方案,但我认为这不是最好的方法。

假设我们正在处理一个LineString:

>>> import shapely
>>> from shapely.geometry import *
>>> import geopandas as gpd
>>> 
>>> MyLine=LineString([(0,0),(5,0),(10,3)])
>>> MyLine
<shapely.geometry.linestring.LineString object at 0x1277EEB0>
>>> 
Run Code Online (Sandbox Code Playgroud)

与该LineString相交的2行:

>>> IntersectionLines=gpd.GeoSeries([LineString([(2,1.5),(3,-.5)]), LineString([(5,.5),(7,.5)])])
>>> IntersectionLines
0    LINESTRING (2 1.5, 3 -0.5)
1     LINESTRING (5 0.5, 7 0.5)
dtype: object
>>> 
Run Code Online (Sandbox Code Playgroud)

看起来像这样: 在此处输入图片说明

我可以得到交点如下:

>>> IntPoints=MyLine.intersection(IntersectionLines.unary_union)
>>> IntPointCoords=[x.coords[:][0] for x in IntPoints]
>>> IntPointCoords
[(2.75, 0.0), (5.833333333333333, 0.5)]
>>> 
Run Code Online (Sandbox Code Playgroud)

然后,我提取起点和终点,并用这些点和相交点创建对,这些对点将用于形成线段:

>>> StartPoint=MyLine.coords[0]
>>> EndPoint=MyLine.coords[-1]
>>> SplitCoords=[StartPoint]+IntPointCoords+[EndPoint]
>>> 
>>> def pair(list):
...     for i in range(1, len(list)):
...         yield list[i-1], list[i]
... 
>>> 
>>> SplitSegments=[LineString(p) for p in pair(SplitCoords)]     
>>> SplitSegments
[<shapely.geometry.linestring.LineString object at 0x127F7A90>, <shapely.geometry.linestring.LineString object at 0x127F7570>, <shapely.geometry.linestring.LineString object at 0x127F7930>]
>>> 

>>> gpd.GeoSeries(SplitSegments)
0                      LINESTRING (0 0, 2.75 0)
1    LINESTRING (2.75 0, 5.833333333333333 0.5)
2      LINESTRING (5.833333333333333 0.5, 10 3)
dtype: object
>>> 
Run Code Online (Sandbox Code Playgroud)

但是,我的方法存在一个问题,就是我知道哪个交点应该与LineString起点连接,以及哪个交点应该与LineString终点配对。如果沿线以与起点和终点相同的顺序列出交点,则此程序将起作用。我想可能会出现并非总是这样的情况吗?有没有一种方法可以概括这种方法,或者完全有更好的方法?

jor*_*ris 5

这是一种更通用的方法:计算所有点(线的起点和终点+您要分割的点)沿线的距离,对这些点进行排序,然后以正确的顺序生成线段。一起发挥作用:

def cut_line_at_points(line, points):

    # First coords of line (start + end)
    coords = [line.coords[0], line.coords[-1]]

    # Add the coords from the points
    coords += [list(p.coords)[0] for p in points]

    # Calculate the distance along the line for each point
    dists = [line.project(Point(p)) for p in coords]

    # sort the coords based on the distances
    # see http://stackoverflow.com/questions/6618515/sorting-list-based-on-values-from-another-list
    coords = [p for (d, p) in sorted(zip(dists, coords))]

    # generate the Lines
    lines = [LineString([coords[i], coords[i+1]]) for i in range(len(coords)-1)]

    return lines
Run Code Online (Sandbox Code Playgroud)

在您的示例上应用此功能:

In [13]: SplitSegments = cut_line_at_points(MyLine, IntPoints)

In [14]: gpd.GeoSeries(SplitSegments)
Out[14]:
0                      LINESTRING (0 0, 2.75 0)
1    LINESTRING (2.75 0, 5.833333333333333 0.5)
2      LINESTRING (5.833333333333333 0.5, 10 3)
dtype: object
Run Code Online (Sandbox Code Playgroud)

唯一的事情是,这不会保留原始线的拐角(但是您在问题中的示例也没有这样做,因此我不知道这是否是必需条件。可能,但会使它变得更复杂)


更新一个版本,该版本保持原始行中的角不变(我的方法是也保留一个0/1列表,该列表指示是否要拆分坐标):

def cut_line_at_points(line, points):

    # First coords of line
    coords = list(line.coords)

    # Keep list coords where to cut (cuts = 1)
    cuts = [0] * len(coords)
    cuts[0] = 1
    cuts[-1] = 1

    # Add the coords from the points
    coords += [list(p.coords)[0] for p in points]    
    cuts += [1] * len(points)        

    # Calculate the distance along the line for each point    
    dists = [line.project(Point(p)) for p in coords]    
?
    # sort the coords/cuts based on the distances    
    # see http://stackoverflow.com/questions/6618515/sorting-list-based-on-values-from-another-list    
    coords = [p for (d, p) in sorted(zip(dists, coords))]    
    cuts = [p for (d, p) in sorted(zip(dists, cuts))]          

    # generate the Lines    
    #lines = [LineString([coords[i], coords[i+1]]) for i in range(len(coords)-1)]    
    lines = []        

    for i in range(len(coords)-1):    
        if cuts[i] == 1:    
            # find next element in cuts == 1 starting from index i + 1   
            j = cuts.index(1, i + 1)    
            lines.append(LineString(coords[i:j+1]))            

    return lines
Run Code Online (Sandbox Code Playgroud)

应用于示例:

In [3]: SplitSegments = cut_line_at_points(MyLine, IntPoints)

In [4]: gpd.GeoSeries(SplitSegments)
Out[4]:
0                           LINESTRING (0 0, 2.75 0)
1    LINESTRING (2.75 0, 5 0, 5.833333333333333 0.5)
2           LINESTRING (5.833333333333333 0.5, 10 3)
dtype: object
Run Code Online (Sandbox Code Playgroud)