Projection of a point to a line segment Python Shapely

Question

I have a LineString defined by two points, so essentially a straight line segment, and I wanted to project a point on to it. I am aware of .project and .interpolate. However when the point is "outside" the segment, I don't want the closest point on the segment, but I want to extend the segment and draw a line going through the point and is orthogonal to the (extended) line segment. I want the coordinate of the projection.

For example if the point is "within" the segment

from shapely.geometry import Point
from shapely.geometry import LineString
point = Point(0.2, 0.5)
dist = LineString([(0, 1), (1, 1)]).project(point)
list(LineString([(0, 1), (1, 1)]).interpolate(dist).coords)

Anyone knows what to do when the point is outside of the segment?

Answer 1

It will be probably easiest to do this manually. If you denote the angle x - u - v as alpha, then

cos(alpha) = (v - u).(x - u) / (|x - u|*|v - u|)

where . denotes the dot product, and | | represents the Euclidean norm. The distance d of P(x) from u is therefore:

d = cos(alpha)*|x - u| = (v - u).(x - u) / |v - u|

Having calculated d, the projected point P(x) is then easily obtained as:

P(x) = u + d*(v - u)/|v - u|

The implementation:

import numpy as np
from shapely.geometry import Point
from shapely.geometry import LineString

point = Point(0.2, 0.5)
line = LineString([(0, 1), (1, 1)])

x = np.array(point.coords[0])

u = np.array(line.coords[0])
v = np.array(line.coords[len(line.coords)-1])

n = v - u
n /= np.linalg.norm(n, 2)

P = u + n*np.dot(x - u, n)
print(P) #0.2 1.

Answer 2

You may also consider using scikit-spatial as follows:

from skspatial.objects import Line

line = Line.from_points(point_a=[0,1], point_b=[1,1])
point = (1.2, 0.5) # example of a point not "within" the line segment

projected_point = line.project_point(point)


print(projected_point)

Output:

[1.2 1. ]