File Exchange

## Distance from points to polyline or polygon

version 2.0.0.0 (6 KB) by
Find minimum distances from points to a polyline or to a closed polygon

Updated 04 Nov 2015

Editor's Note: This file was selected as MATLAB Central Pick of the Week

The following files are included:
p_poly_dist.m - Compute the distances from each one of a set of np points p(1), p(2),... p(np) on a 2D plane to a polyline or a closed polygon. Polyline is defined as a set of nv-1 segments connecting nv ordered vertices v(1), v(2), ..., v(nv). The polyline can optionally be treated as a closed polygon.
Distance from point j to a segment k is defined as a distance from this point to a straight line passing through vertices v(k) and v(k+1), when the projection of point j on this line falls INSIDE of segment k; and to the closest of v(k) or v(k+1) vertices, when the projection falls OUTSIDE of segment k. Distance from point j to a polyline is defined as a minimum of this point's distances to all segments.
In a case when the projected point fall OUTSIDE of ALL polyline segments, the distance to a closest vertex of the polyline is returned
test_p_poly_dist.m - a simple unit test for p_poly_dist. Plots the results of a call to p_poly_dist function (see help for usage example)
p_poly_dist1.m – previous version (ver. 1.0) of p_poly_dist.m
Notes on version 2.0:
Oct 2, 2015 - version 2.0 (Michael Yoshpe). The original ver. 1.0 function was completely redesigned. The main changes introduced in Ver. 2.0 are:
1. Distances to polyline (instead of a closed polygon in Ver. 1.0) are returned. The polyline can optionally be treated as a closed polygon.
2. Distances from multiple points to the same polyline can be found.
3. The algorithm for finding the closest points is now based on coordinate system transformation. The new algorithm avoids numerical problems that Ver. 1.0 algorithm could experience in "ill-conditioned" cases.
4. Many optional output variables were added. In particular, the closest points on polyline can be returned.

### Cite As

Michael Yoshpe (2021). Distance from points to polyline or polygon (https://www.mathworks.com/matlabcentral/fileexchange/12744-distance-from-points-to-polyline-or-polygon), MATLAB Central File Exchange. Retrieved .

Michael Yoshpe

Hello, Rebecca, thank you very much for your rating!
I haven't tested the function with self-intersected polygons, but I think it should work, because the algorithm computes the distance to each segment independently, and then finds the minimum. It won't work for multiply connected polygon though, simply because there is currently no code inside the function that would separate the polygons, like in Matlab inpolygon function.
Regards
Michael

Rebecca Jiang

This function is great -- it is efficient, accurate, and very useful.

One question: does the polygon have to be simply connected for this function to compute the distance accurately?

Thank you!

Gerald Koudijs

Great function, it save me a lot of time. As I needed did for a lot of data points, I needed to smartly select the points and crop the polyline to make to improve the computation time. I think the function itself could be optimized a bit as well to improve the computation time. Functionally however perfect.

Michael Yoshpe

Hello, Wai Kei
I am glad you've found my function useful, but I don't think I fully understand what do you mean by saying "why need input the same dimensions of points and polygon?". As stated in the description, the function was specifically designed to find the distances of MULTIPLE point to a SINGLE polyline or polygon. Please see the description of the input variables:
% xp - vector of points X coordinates (1 X np)
% yp - vector of points Y coordinates (1 X np)
% xv - vector of polygon vertices' X coordinates (1 X nv)
% yv - vector of polygon vertices' Y coordinates (1 X nv)

The np and nv DO NOT have to be of the same value.
I suggest you have a look at test_p_poly_dist function. Please uncomment the lines after usage example and run the code, i think it might help you.

Best regards

Michael

wai kei

Important code, and help me work a lot.
but one question, why need input the same dimensions of points and polygon? If I want to lots of points, but only one polygon, do I need to split the points to the same dimensions of polygon,then use the command?

thanks a lot!

Juan Marrugo

This code was very useful to me. I needed the distance from a fault line to a set of points. At the beginning I thought of measuring the distance to the vertices in the fault line, but then I realized I needed to account for the closest distance to be in the segments as well. This code saved me a lot of time, thanks a lot. JC

Songqiu

Dear Micheal,
The code is very useful. But I want your comment on another code written for the same purpose.
The code name is "dpoly". It also find the distance between point(s) and a polygon. You can find it here
https://people.sc.fsu.edu/~jburkardt/m_src/dist_plot/dpoly.m

Can you please comment on this?
I would greatly appreciate it if you kindly give me some feedback on this.

Regards

Samuel Londner

Very useful function, fast and exact.

foboo

Dear Michael,
Thanks so much for this great function!

I have one question/feature request.
Based on the current output of the function, I can derive what I call the "closest line segment function", i.e., the line segment on which the closest point has been located.
Let me denote this function as
f(t) = a + t*(b-a)
where "t" is between 0 and 1, and where the points "a" and "b" can already be deduced from the current output of the function.

However, I am also looking for the value of "t" which corresponds to the smallest distance (as computed).

Currently, I am trying to compute value of "t" via
t = ([x*, y*] - a )/ (b-a),
where x*,y* are the closest points on the polyline (as returned by the function).
This is cumbersome, as it gives me two possible solutions for "t".

Would be possible that the function itself also returns the value of "t" corresponding to the smallest distance?

foboo

I did not mean
t = (distance - a )/ (b-a).

t = ([x*, y*] - a )/ (b-a),
where x*,y* are the closest points on the polyline (as returned by the function).

snoop58

Hi Michael,
thanks to your reply, I started digging and found that my polyline generator does not always make consecutive connected points. Moreover, the function trisurf() that I have used to display the polyline, plots a confusing edge-only graph, even when the points are not consecutive!
Anyhow, after passing a convhull() filter over my polylines, I was able to apply your function like a charm. Thank you!

Michael Yoshpe

Hello, snoop58
Thank you very much for your comment.
Unfortunately, I could not reproduce the problem you mentioned entirely, since I could only load the polylines without the test code (the code itself is in pdf file instead of .m file).
Still, I have a few remarks.
First, setting the optional find_in_out flag to 'true' treats any OPEN polyline as A CLOSED one. I quote from the function's help section:
% find_in_out - logical flag. When true, the polyline is treated as a
% closed polygon, and each input point is checked wether it is inside or
% outside of this polygon. In such case, an open polyline is automatically
% closed by adding a last point identical to a first one.

So, using find_in_out = 'true', effectively adds a segment to the polyline which connects last point to the first one, and the polylines in your figures are plotted as OPEN without this closing segment.

Second, I could not reproduce the polylines themselves as they are plotted in your document. For example, I used the following command to plot the A polyline for the wrong case:

S2 = polylines.S2;
plot(S2{1}.vertices(:,1), S2{1}.vertices(:,2), 'm-')

and the line I got is very different from what you plotted (try this command yourself). Note that xv,yv input vectors are treated as the coordinates of consecutive connected points.

I am sorry I could not be of more help to you. I suggest you read the help section again and check the input data. If you still have problems with the results, please give both the data AND the code as .m function, so I could try to debug it further.

Best regards

Michael

snoop58

Hello Michael. Great work, thank you.
The function is very nice and does exactly what I want; however, it fails in some cases. I attach a PDF and mat file, which contain a couple of screenshots of correct and incorrect distance calculation, a mat file with polyline triangulation, and a code to run the data.
https://www.dropbox.com/s/dg77yepl0btsykg/distance_test.pdf?dl=0
https://www.dropbox.com/s/nt17douftlyhj8u/polylines.mat?dl=0

In short: there are 2 sets of 4 open polylines and a single point in fixed location.
In the first case, the point is outside all 4 polylines, and the shortest distances to each polyline are calculated correctly.
In the second case, the point is outside top & bottom polylines (distance is correct) but inside the left & right polylines: the min distance is wrong, as well as projection coordinates.
I have tried various polyline combinations, and I always get wrong distance calculation when the left/right polylines cover the point of interest. I have also tried to keep the polylines open or closed (find_in_out = true/false) with identical results.

Michael Yoshpe

Hello, Philippe
I ran the unit test test_p_poly_dist (included with my submission) on Matlab R2015a, and it was OK.

Please let me know which version of Matlab caused the problems you observed, and what were the inputs, specifically the coordinates of vertices (xv, yv), and points (xp, yp).

Best regards

Michael

Philippe Lebel

also, on line 291 (or so...)

d_min(is_vertex) = dpv_min(is_vertex);

Philippe Lebel

Hi, I found out about some bugs that might have been introduced by new matlab releases.

When trying to compute multiple points, the function did not work at all.

here are the fixes i had to implement in order to make it work:

at line 219:
if length(xp)<= 1
Ppr1 = Ppr1 - a1';
Ppr2 = Ppr2 - a2';
else
Ppr1 = Ppr1 - a1;
Ppr2 = Ppr2 - a2;
end

then at line 316:
I_cr_min_temp = I_cr_min(~is_vertex);
if ~isempty(I_cr_min_temp)
I_cr_min_temp = I_cr_min_temp + [0:length(I_cr_min_temp)-1]'*(nv - 1);
xc_temp = xc';
yc_temp = yc';
x_d_min = xc_temp(I_cr_min_temp);
y_d_min = yc_temp(I_cr_min_temp);
end

Turgut S.

David S

Hi Michael, thanks for this very useful function. Here is an example of a group of polygons that form "donut" holes. Interestingly, the code correctly parses the right donut hole, but at the left one, it adds a line connecting the outer polygon to the inner one. Is there a way to avoid this? https://www.dropbox.com/s/9yxilcfx635bvvi/p_poly_dist_example.m?dl=0

Michael Yoshpe

Charles, you might want to have a look on a p_spoly_dist function I just added (http://www.mathworks.com/matlabcentral/fileexchange/52734-p-spoly-dist).
I think it's pretty close to what you want, EXCEPT it DOES NOT for now return an in/out indication for a point in a closed spherical polygon. I simply could not find a reliable algorithm for this task.

Pablo Nilo

Very very useful, perfect for geographical data filtering and screening. Thanks a lot!!

Just a question from an Earth Observartion investigator: Could it be faster? We need to deal with satellite tracks (120.000 records each orbit).

charles

very helpful, thank you. however, is 2D. is there a comparable function for finding the distance on surface of a sphere, from point on the sphere's surface to polygon on the sphere's surface? point and vertices given in lat-long?

Ilya

To Eric Schmitz's comment:
why we need "clear id"? By the way, it consumes quite a lot of time.

Trey

Sorry, my previous comment is in regards to Bastian's update, not the original by Michael Yoshpe

Trey

Hi, thanks for the code. It's fast, but the 23 July 2013 code has a bug.

There's an 'X' in

uY = [Y(2:end, :); X(1, :)] - Y;

that should be 'Y':

uY = [Y(2:end, :); Y(1, :)] - Y;

Bastian

I've changed the code a little bit. So it can be used with the same syntax as the MATLAB-Funktion inpolygon.

function d = distancepolygon(x, y, X, Y)
m = size(X, 1) * size(X, 2);
n = size(x, 1) * size(x, 2);

X = X(:) * ones(1, n);
Y = Y(:) * ones(1, n);

x = ones(m, 1) * x(:)';
y = ones(m, 1) * y(:)';

% projection
uX = [X(2:end, :); X(1, :)] - X;
uY = [Y(2:end, :); X(1, :)] - Y;

LAMBDA = ((x - X) .* uX + (y - Y) .* uY) ./ (uX.^2 + uY.^2);

PX = X + LAMBDA .* uX;
PY = Y + LAMBDA .* uY;

% distance to projections points which is inside a segment
PD = NaN(size(x));
b = ((LAMBDA > 0) & (LAMBDA < 1));
PD(b) = sqrt((PX(b) - x(b)).^2 + (PY(b) - y(b)).^2);

% distance to polygon point
D = sqrt((X - x).^2 + (Y - y).^2);

% distance
d = min([D; PD]);

% change sign if the point is inside the polygon
b = inpolygon(x(1, :), y(1, :), X(:, 1), Y(:, 1));
d(b) = -d(b);

Xavier Nieves

Hi, I have been trying to create a function similar to this one... if I were to have a random polygon but always have as an input the x,y coordinates of the vertices's and and the x,y coordinate of a point. can I use this function to determine the minimum distance of the point to the polygon, and if I can should express the inputs in what way... ?? nice function

Cristian D

Thanks for this function! it was very usefull!

Alejandro Weinstein

I added the suggested test from Eric Schmitz, and also return the point of the polygon closest to the point. It is available here:

Eric Schmitz

For the case where a polygon rib is either horizontal or vertical, the test for checking if the projected point is inside a segment may not work all the time because of numerical accuracy. I added a test to force either the x or y coordinate of projected point to be set equal to the corresponding rib x or y coordinate (before the section computing idx_x and idx_y):

id = find(diff(xv)==0);
xp(id)=xv(id);
clear id
id = find(diff(yv)==0);
yp(id)=yv(id);

Jesse Collins

this function would benefit greatly from also returning the point on the given polygon closest to the point from which the distance is to be found.

Kardi Teknomo

##### MATLAB Release Compatibility
Created with R12
Compatible with any release
##### Platform Compatibility
Windows macOS Linux