# Distance between point and line segments

How would you determine the shortest distance between a point and one or more segments?

For example, what is the shortest distance between the point and the two segments below? Clearly the point is closer to the right line segment. That means the shortest distance between the point and the line segments is equal to the distance between the point and the right segment.

point = {{8, 15}}
lines = {{{20, 10}, {11, 27}}, {{11, 27}, {1, 27}}}
Show[Graphics[{Thick, Line@lines}],
Graphics[{PointSize[Large], Pink, Point@point}]]


Can you think of a way that also works with an arbitrary number of segments?

-
@YvesKlett Your link refers to the theory shortest distance. I would like to know how you can calculate distance in practice (with Mathematica of course). – sjdh Feb 24 '13 at 9:02
When you say "arbitrary number of lines", do you mean line, line segment, or ray? – Silvia Feb 24 '13 at 12:10
@Silvia I mean the union of the line segments – sjdh Feb 24 '13 at 12:35

In version 10 (now available publicly through the Programming Cloud), use RegionNearest

point = {8, 15};
line = Line[{{{20, 10}, {11, 27}}, {{11, 27}, {1, 27}}}];

np = RegionNearest[line, point]
(* {5663/370, 6981/370} *)

Graphics[{PointSize[Large], {line, Point[point]}, {Red, Point[np]}}]


-
OR since the damn Programming cloud has not worked (Server error) since it was released, just use the same RegionNearest in v9 as I did here – RunnyKine Jun 24 '14 at 15:11

Using the parametric equation of the lines.

point = {8, 15}
lines = {{{20, 10}, {11, 27}}, {{11, 27}, {1, 27}}}

distance[{start_, end_}, pt_]:= Module[{param = ((pt - start).(end - start))/Norm[end - start]^2},
EuclideanDistance[pt, start + Clip[param , {0, 1}] (end - start)]];

Min[distance[#, point] & /@ lines]

(*159/Sqrt[370]*)


Plotting isodistance lines

l  = Table[{Cos[2 Pi n/5], Sin[2 Pi n/5]}, {n, 1, 6}];
l1 = Partition[Riffle[l, l[[2 ;;]]], 2];
Quiet@ContourPlot[Min[distance[#, {x, y}] & /@ l1], {x, -2, 2}, {y, -2, 2}]


-

Not as elegant as belisarius' first answer but was a bit faster; now slower than the new one. A bit of explanation: If the point on the line segment $PQ$ closest to $X$ is on $PQ$, then the exterior (supplementary) angles of the triangle $PXQ$ at $P$ and $Q$ will be greater than a right angle and the cosines will be negative. [Edit: To clarify, cosAngles is actually a list of Dot products, which are just positive numbers times the cosines; and we only use the sign.]

lines = RandomReal[{0, 10}, {10, 2, 2}];
point = RandomReal[{0, 10}, {2}];

cosAngles = Dot @@@ Partition[Differences[{#1, #2, #3, #1}], 2, 1] &;
nearest[{P_, Q_}, X_] := Module[{points},
If[And @@ Negative[cosAngles[X, P, Q]],
P + Projection[X - P, Q - P], First@Nearest[{P, Q}, X]]
];
findNearest[lines_, point_] :=
Nearest[nearest[#, point] & /@ lines, point];
distance[lines_, point_] :=
Norm[First@findNearest[lines, point] - point];


Here's the output:

point
findNearest[lines, point]
distance[lines, point]

{7.48559, 3.5353}

{{8.20107, 4.23693}}

1.0021


Here's a picture of what's going on:

g := Module[{distPts, nearPt},
distPts = nearest[#, point] & /@ lines;
nearPt = First@Nearest[distPts, point];
Graphics[{Line@lines, Blue, Thin, Line[{point, #} & /@ distPts],
Point[nearest[#, point] & /@ lines], Thick, Line[{point, nearPt}],
Red, Point@point}]
]

g


Same lines, different point. The nearest point is an end point of a segment.

Block[{point = {8, 9}},
Print[distance[lines, point]];
g
]

1.46937


-
I think my new solution is faster :) – Dr. belisarius Feb 24 '13 at 23:06
@belisarius Touché! Nice work. About twice as fast, or a little better, on a 1000 lines. – Michael E2 Feb 24 '13 at 23:16
Actually 60% faster (1.3/0.8). I really liked the other approach (it seemed more "elegant"), but it was really slow – Dr. belisarius Feb 24 '13 at 23:18
@belisarius Well I can't find a different idea that's faster. As a math teacher, I admire your geometric distance function more than using MinValue. It's faster still if you skip Module and insert the param formula directly. – Michael E2 Feb 25 '13 at 11:28

The numbers shown are the normal distances. 8.266 is the smaller distance. ref above for the equation of the normals.

x0 = 8; y0 = 15;
point = {{x0, y0}};
(*x1=8;y1=20; x2=10;y2=22; x3=12;y3=28*)
x1 = 20; y1 = 10; x2 = 11; y2 = 27; x3 = 1; y3 = 27;
lines = {{{x1, y1}, {x2, y2}}, {{x2, y2}, {x3, y3}}};
v1 = {y2 - y1, -(x2 - x1)};
r1 = {x1 - x0, y1 - y0};
v2 = {y3 - y2, -(x3 - x2)};
r2 = {x3 - x0, y3 - y0};

d1 = Abs@Dot[v1, r1]/Norm[v1];
d2 = Abs@Dot[v2, r2]/Norm[v2];

Labeled[
Graphics[{
{Thick, Line@lines}, {PointSize[Large], Pink, Point@point}
}, Axes -> True, GridLines -> {Range[5, 20, 1], Range[10, 30, 1]},
GridLinesStyle -> LightGray],
Grid[{{"distance to first line ", N@d1},
{"distance to second line ", N@d2}}, Alignment -> Left]
]


-
This works for the line segments given by the OP, but doesn't seem to work in some other case, such as when the line segments are changed to x1 = 8; y1 = 20; x2 = 10; y2 = 22; x3 = 12; y3 = 28; – Mike Z. Feb 24 '13 at 14:32

This solution I am about to present is effectively a blend of Michael's and belisarius's approaches. To wit,

PointLineDistance[pt_, {s1_, s2_}] :=
With[{tp = s1 - pt}, EuclideanDistance[tp, Projection[tp, s2 - s1]]]

segs = {{{20, 10}, {11, 27}}, {{11, 27}, {1, 27}}};
nf = Nearest[segs -> Automatic, DistanceFunction -> PointLineDistance];

PointLineDistance[{8, 15}, Extract[segs, nf[{8, 15}]]]
159/Sqrt[370]


Trying to reproduce bel's isodistance plot on my box gives something rather different:

segs = Partition[Table[Through[{Cos, Sin}[2 Pi n/5]], {n, 5}] // N, 2, 1, 1];
nf = Nearest[segs -> Automatic, DistanceFunction -> PointLineDistance];

ContourPlot[PointLineDistance[{x, y}, Extract[segs, nf[{x, y}]]], {x, -2, 2}, {y, -2, 2},
AspectRatio -> Automatic, ColorFunction -> "ThermometerColors"]


Finally, here's a three dimensional test:

BlockRandom[SeedRandom[42, Method -> "MKL"]; (* for reproducibility *)
segs = RandomReal[{-1, 1}, {50, 2, 3}]];

nf = Nearest[segs -> Automatic, DistanceFunction -> PointLineDistance];

Graphics3D[{Line[segs],
{Directive[Red, AbsoluteThickness[4]], Line[Extract[segs, nf[{0, 0, 0}]]]},
{Directive[Red, AbsolutePointSize[6]], Point[{0, 0, 0}]}}]


-
I think that nf in the first block finds the segment whose extended line is closest to a point. But the problem is about distance to a segment only, not to its extended line, isn't it? – BoLe Apr 19 '13 at 9:01
How then do you define "distance to a segment"? – J. M. Apr 19 '13 at 9:02
Something like Min[EuclideanDistance[p1 + t (p2 - p1), p3]] with 0 < t < 1. – BoLe Apr 19 '13 at 9:32
p1 = {{8, 15}};
polyline = {{{20, 10}, {11, 27}}, {{11, 27}, {1, 27}}};


Discretizing the segments:

(* decrease dt for better result *)
data = With[{dt = .1},
Flatten[polyline /. r : {{_, _}, {_, _}} :>
Table[s[[1]] + t (s[[2]] - s[[1]]), {s, r}, {t, 0., 1, dt}], 1]];

p2 = Nearest[data, p][[1]]
(* {{15.5, 18.5}} *)

EuclideanDistance @@ Join[p1, p2]
(* 8.27647 *)

Graphics[{
Thick, Line@polyline,
PointSize[Large], Pink, Point[p1~Join~p2]}]


I think this should work in 3D too and with any number of joined or disjoint segments.

-

Another approach to this problem is to use the image processing functions and operate directly on the image of the lines. Here we import the two lines from the OPs question, and take the DistanceTransform. Finding the value of the distance transform at the orange point gives the distance (in pixels) of that point from that line.

img = Import["http://i.stack.imgur.com/Inf9t.png"];
point = Round[First@ImageValuePositions[img, RGBColor[1, 0.5, 0.5]]];
distImg = DistanceTransform[img];

91.418