How to get equally spaced points on a line in Shapely

Question

I'm trying to (roughly) equally space the points of a line to a predefined distance.

It's ok to have some tolerance between the distances but as close as possible would be desirable.

I know I could manually iterate through each point in my line and check the p1 distance vs p2 and add more points if needed.

But I wondered if anyone knows if there is a way to achieve this with shapely as I already have the coords in a LineString.

Sleuthhound

Answer 1

One way to do that is to use interpolate method that returns points at specified distances along the line. You just have to generate a list of the distances somehow first. Taking the input line example from Roy2012's answer:

import numpy as np
from shapely.geometry import LineString
from shapely.ops import unary_union

line = LineString(([0, 0], [2, 1], [3, 2], [3.5, 1], [5, 2]))

Splitting at a specified distance:

distance_delta = 0.9
distances = np.arange(0, line.length, distance_delta)
# or alternatively without NumPy:
# points_count = int(line.length // distance_delta) + 1
# distances = (distance_delta * i for i in range(points_count))
points = [line.interpolate(distance) for distance in distances] + [line.boundary[1]]
multipoint = unary_union(points)  # or new_line = LineString(points)

Note that since the distance is fixed you can have problems at the end of the line as shown in the image. Depending on what you want you can include/exclude the [line.boundary[1]] part which adds the line's endpoint or use distances = np.arange(0, line.length, distance_delta)[:-1] to exclude the penultimate point.

Also, note that the unary_union I'm using should be more efficient than calling object.union(other) inside a loop, as shown in another answer.

Splitting to a fixed number of points:

n = 7
# or to get the distances closest to the desired one:
# n = round(line.length / desired_distance_delta)
distances = np.linspace(0, line.length, n)
# or alternatively without NumPy:
# distances = (line.length * i / (n - 1) for i in range(n))
points = [line.interpolate(distance) for distance in distances]
multipoint = unary_union(points)  # or new_line = LineString(points)

Stymie

Wow, that's perfect. – Sleuthhound

I've just realised that points = [line.interpolate(distance) for distance in distances] + [line.boundary[1]] flips the axis. I think you need boundry[0] – Sleuthhound

@LewisMorris What do you mean by "flips the axis"? I cannot reproduce the issue with the given example. I'm also not sure how changing boundary[1] to boundary[0] would solve the problem since boundary[0] should already be included as the first value of distances is zero. – Stymie

line.boundary isn't subscriptable anymore. You need line.boundary.geoms[-1]. – Bullnecked

Oh, and if Linestring is a loop, line.boundary won't have any geometry, so you should check line.boundary is truthy before getting the last point. – Bullnecked

Answer 2

You can use the shapely substring operation:

from shapely.geometry import LineString
from shapely.ops import substring

line = LineString(([0, 0], [2, 1], [3,2], [3.5, 1], [5, 2]))

mp = shapely.geometry.MultiPoint()
for i in np.arange(0, line.length, 0.2):
    s = substring(line, i, i+0.2)
    mp = mp.union(s.boundary)

The result for this data is given below. Each circle is a point.

Navarino

Answer 3

You can use shapely.segmentize:

line = gpd.read_file("line.geojson")
line = shapely.segmentize(line,max_segment_length=100)

The points on the line won't strictly be equally spaced, but if max_segment_length is much higher than the number of points in the original line, the distances will be very close to each other.

Betancourt