Distance between a point and a line
... easy for the eye but a challenge in math & code
Published: 2016-06-05

Some problems are easy to solve using pen and paper, but suddenly become hard when you have to express things in mathematics. In this blog we will discuss one such case:
What is the shortest distance between a point and a line (or line segment, or ray)?
And what are the coordinates of the point on the line that is closest?
We will go from idea to code in four steps:
- We sketch the solution for a line, a line segment and a ray.
- We use linear algebra to define the vector equations that solve the problem in an N-dimensional space.
- We implement the equations in M-code
- We turn the code snippets into a robust MATLAB function
The resulting MATLAB function is available on GitHub.
Drawing the solution
We solve the problem graphically, by drawing it in 2D. Two points, A and B, define the line, line segment or ray. P is the independent point. The figure below shows how we construct the distance (red, dashed) and closest point (C). We consider three cases:
- A line through A and B (the top row)
- A line segment AB (the middle row)
- A ray which starts at A and passes through B (the bottom row)
The difference between these cases becomes clear when we consider the intersection point of the line through A and B and the perpendicular through P. The closest point can be the intersection point, or A, or B:
- For a line, the intersection point is the closest point.
- For a line segment, if the intersection point is on line segment AB, this is the closest point (middle column). If the intersection point is outside segment AB, the closest point is the end of the segment (left and right column).
- For a ray, if the intersection point is on the ray, this is the closest point (middle and right column). If the intersection point is not on the ray (behind source A), the closest point is source A (left column).
This approach provides insight, but is limited: we need a ruler to measure the distance, an it only works in 2D.
Putting it into math
The next step is to define a set of vector equations that solve the problem in N dimensions.
Let \(\vec{A}\) and \(\vec{B}\) be two points on line \(\vec{L}\). \(\vec{P}\) is the independent point, which is not on the line.
The vector equation of the line is:
\(\vec{L}(t) = \vec{A} + t\vec{M}\) (Eq. 1)
with \(t\) the running parameter and \(\vec{M}\) the direction vector:
\(\vec{M} = \vec{B} – \vec{A}\)
The running parameter can take different values for a line, line segment and ray:
- A line is infinite, so \(t \in (-\infty, \infty)\)
- A line segment is finite, so \(t \in \left[0, 1\right]\)
- A ray is semi-infinite, so \(t \in [0, \infty)\)
We need to find the orthogonal projection of our independent point \(\vec{P}\) on line \(\vec{L}\). This can be done by using dot products, as follows:
\(\vec{L}(t_0) = \vec{A} + \left[\frac{(\vec{P}-\vec{A}) \cdot \vec{M}}{\vec{M} \cdot \vec{M}}\right] \vec{M}\) (Eq. 2)
From Eq. 1 and Eq. 2 follows that \(t_0\), the value of the running parameter at the intersection of \(\vec{L}\) and its perpendicular through \(\vec{P}\) is:
\(t_0 = \frac{(\vec{P}-\vec{A}) \cdot \vec{M}}{\vec{M} \cdot \vec{M}}\)
Now that we have all the parameters, we can calculate the shortest distance \(d\).
In case of a line:
\(d = \left|\vec{P} – (\vec{A} + t_0\vec{M})\right|\)
In case of a line segment:
\(d = \left\{ \begin{array}{ll} \left| \vec{P} – \vec{A} \right| & t_0 \leq 0 \\ \left| \vec{P} – (\vec{A} + t_0\vec{M})\right| & 0 < t_0 < 1 \\ \left| \vec{P} – \vec{B} \right| & t_0 \geq 1 \end{array}\right.\)
And in case of a ray:
\(d = \left\{ \begin{array}{ll} \left| \vec{P} – \vec{A} \right| & t_0 \leq 0 \\ \left| \vec{P} – (\vec{A} + t_0\vec{M})\right| & t_0 > 0 \end{array}\right.\)
Drafting the implementation
We only need a few lines of code to implement the linear algebra equations MATLAB.
- Direction vector \(\vec{M}\):
M = B - A;
- Running parameter at the orthogonal intersection \(t_0\):
t0 = dot(M, P - A) / dot(M, M);
- Intersection point:
intersectPnt = A + t0 * M;
- Distance \(d\) between independent point and closest point:
d = norm(P-C);
Implementing a function
The last step involves coding a robust, documented, and readable MATLAB function. At least three input arguments are required: the points A and B that define a line and test point P. The optional fourth input argument specifies the line type: 'line' (the default), 'segment' or 'ray'. The function returns up to three outputs: distance d, closest point C, and running parameter at the orthogonal intersection t0. The function definition line reflects this:
function [d, C, t0] = distancePoint2Line(A, B, P, varargin)
% distancePoint2Line Calculate the distance between a point and a line
% D = distancePoint2Line(A, B, P) returns the distance from point P to the
% line through A and B in an N-dimensional space.
%
% D = distancePoint2Line(A, B, P, LineType) returns the distance from
% point P to
% - a line through A and B,
% - a ray which starts at A and passes through B, or
% - a line segment AB.
% The definition of distance D depends on LineType, see below.
%
% [D, C, t0] = distancePoint2Line(A, B, P, ..) returns in addition the
% closest point C and the running parameter t0 that defines the intersection
% point of the line through A,B and the perpendicular through P. The
% definition of closest point depends on LineType, see below.
Get the distancePoint2Line() function on GitHub to inspect the full implementation and to use it for your work.