Calculating the shortest distance between two lines (line segments) in 3D

Question

I have two line segments: X1,Y1,Z1 - X2,Y2,Z2 And X3,Y3,Z3 - X4,Y4,Z4

I am trying to find the shortest distance between the two segments.

I have been looking for a solution for hours, but all of them seem to work with lines rather than line segments.

Any ideas how to go about this, or any sources of furmulae?

Unknown Unknown · Accepted Answer · 2009-03-31T17:36:02

I'll answer this in terms of matlab, but other programming environments can be used. I'll add that this solution is valid to solve the problem in any number of dimensions (>= 3).

Assume that we have two line segments in space, PQ and RS. Here are a few random sets of points.

> P = randn(1,3)
P =
     -0.43256      -1.6656      0.12533

> Q = randn(1,3)
Q =
      0.28768      -1.1465       1.1909

> R = randn(1,3)
R =
       1.1892    -0.037633      0.32729

> S = randn(1,3)
S =
      0.17464     -0.18671      0.72579

The infinite line PQ(t) is easily defined as

PQ(u) = P + u*(Q-P)

Likewise, we have

RS(v) = R + v*(S-R)

See that for each line, when the parameter is at 0 or 1, we get one of the original endpoints on the line returned. Thus, we know that PQ(0) == P, PQ(1) == Q, RS(0) == R, and RS(1) == S.

This way of defining a line parametrically is very useful in many contexts.

Next, imagine we were looking down along line PQ. Can we find the point of smallest distance from the line segment RS to the infinite line PQ? This is most easily done by a projection into the null space of line PQ.

> N = null(P-Q)
N =
     -0.37428     -0.76828
       0.9078     -0.18927
     -0.18927      0.61149

Thus, null(P-Q) is a pair of basis vectors that span the two dimensional subspace orthogonal to the line PQ.

> r = (R-P)*N
r =
      0.83265      -1.4306

> s = (S-P)*N
s =
       1.0016     -0.37923

Essentially what we have done is to project the vector RS into the 2 dimensional subspace (plane) orthogonal to the line PQ. By subtracting off P (a point on line PQ) to get r and s, we ensure that the infinite line passes through the origin in this projection plane.

So really, we have reduced this to finding the minimum distance from the line rs(v) to the origin (0,0) in the projection plane. Recall that the line rs(v) is defined by the parameter v as:

rs(v) = r + v*(s-r)

The normal vector to the line rs(v) will give us what we need. Since we have reduced this to 2 dimensions because the original space was 3-d, we can do it simply. Otherwise, I'd just have used null again. This little trick works in 2-d:

> n = (s - r)*[0 -1;1 0];
> n = n/norm(n);

n is now a vector with unit length. The distance from the infinite line rs(v) to the origin is simple.

> d = dot(n,r)
d =
       1.0491

See that I could also have used s, to get the same distance. The actual distance is abs(d), but as it turns out, d was positive here anyway.

> d = dot(n,s)
d =
       1.0491

Can we determine v from this? Yes. Recall that the origin is a distance of d units from the line that connects points r and s. Therefore we can write dn = r + v(s-r), for some value of the scalar v. Form the dot product of each side of this equation with the vector (s-r), and solve for v.

> v = dot(s-r,d*n-r)/dot(s-r,s-r)
v =
       1.2024

This tells us that the closest approach of the line segment rs to the origin happened outside the end points of the line segment. So really the closest point on rs to the origin was the point rs(1) = s.

Backing out from the projection, this tells us that the closest point on line segment RS to the infinite line PQ was the point S.

There is one more step in the analysis to take. What is the closest point on the line segment PQ? Does this point fall inside the line segment, or does it too fall outside the endpoints?

We project the point S onto the line PQ. (This expression for u is easily enough derived from similar logic as I did before. Note here that I've used \ to do the work.)

> u = (Q-P)'\((S - (S*N)*N') - P)'
u =
      0.95903

See that u lies in the interval [0,1]. We have solved the problem. The point on line PQ is

> P + u*(Q-P)
ans =
      0.25817      -1.1677       1.1473

And, the distance between closest points on the two line segments was

> norm(P + u*(Q-P) - S)
ans =
        1.071

Of course, all of this can be compressed into just a few short lines of code. But it helps to expand it all out to gain understanding of how it works.

Calculating the shortest distance between two lines (line segments) in 3D

6 Answers