Split lines at points using Shapely

Question

I'm trying to split a shapefile of geometry type 'line' at a shapefile of geometry type 'point. These points are snapped to the lines and are the start points/end points of other line shapefiles.

Something like this : The points are snapped to the black lines.

splitAtPoints

As of now I'm working with Geopandas/Shapely/Fiona. I've searched for a native function in Shapely that performs the split; however, I've had no luck in finding a solution. This is also true when searching through Shapely's official docs. This is odd since the split is among the basic geoprocessing tools.

My goal is to give a function of the two shapefiles (lines & points) and split the lines at points (like ArcGIS' tool). The split lines should then be stored in a different shapefile. This seems easy to do but the lack of documentation is a bump in the road.

Any ideas on how I can do this in python using Geopandas/Shapely?

Answer 1

You have many solutions and I use here a simple example

1) the easiest way

from shapely.geometry import Point, LineString
line = LineString([(1,2),(2,4),(4,5)])
point = Point(2,4)

First, you must determine if the point is within the line (Determine if shapely point is within a linestring/multilinestring)

line.distance(point) < 1e-8
True
print LineString([line.coords[0],point.coords[:][0]])
LINESTRING (1 2, 2 4)
print LineString([point.coords[:][0], line.coords[-1]])
LINESTRING (2 4, 4 5)

2) with the shapely function split

from shapely.ops import split
result = split(line, point)
result.wkt
'GEOMETRYCOLLECTION (LINESTRING (1 2, 2 4), LINESTRING (2 4, 4 5))'

3) from Get the vertices on a LineString either side of a Point

from shapely.geometry import Point,LineString 
def split(line_string, point):
 coords = line_string.coords
 j = None 
 for i in range(len(coords) - 1):
 if LineString(coords[i:i + 2]).intersects(point):
 j = i
 break 
 assert j is not None 
 # Make sure to always include the point in the first group
 if Point(coords[j + 1:j + 2]).equals(point):
 return coords[:j + 2], coords[j + 1:]
 else:
 return coords[:j + 1], coords[j:]
line1,line2 = split(line,point)
line1 = LineString(line1)
line2 = LineString(line2)
print line1, line2
LINESTRING (1 2, 2 4) LINESTRING (2 4, 4 5)

4) from Shapely Split LineStrings at Intersections with other LineStrings

# First coords of line (start + end)
coords = [line.coords[0], line.coords[-1]] 
# Add the coords from the points
coords += point.coords
# Calculate the distance along the line for each point
dists = [line.project(Point(p)) for p in coords]
# sort the coordinates
coords = [p for (d, p) in sorted(zip(dists, coords))]
lines = [LineString([coords[i], coords[i+1]]) for i in range(len(coords)-1)]
for lin in lines:
 print lin
LINESTRING (1 2, 2 4)
LINESTRING (2 4, 4 5)

5) you can also examine and adapt

• cut.py by Sean Gillies
• Efficient algorithms to split and join lines?
• Shapely: Cut a piece from a linestring at two cutting points
• ....

• GeoSal: Thank you @gene for your detailed answer! At first I was looking for a native function in Shapely that does that.
• gene: The function split
• GeoSal: The split function in Shapely is not recognized apparently. I tried importing : from shapely.ops import split But I get an unresolved error. I install the split package and nothing happens still. Very odd. I already have the points and lines in shapefiles so all I need to do is the split. The point are already projected on the lines.
• gene: look at the github version of shapely ops.py line 329
• Chau: The split function could be nice, but no matter how I try to snap the point onto the line, it won't split the line at that point.