version 3.0.0.0 (138 KB) by
Darren Engwirda

A fast test to determine point inclusion for general polygonal geometries.

INPOLY computes the intersection between a set of points and a general polygonal region in the plane, returning the 'inside', 'outside' and 'boundary' status for each vertex. General non-convex and multiply-connected polygonal regions can be handled. INPOLY is intended as a fast replacement for MATLAB's default INPOLYGON routine.

See POLYDEMO to get started with a set of example problems:

polydemo(1); % a simple example

polydemo(2); % multiply-connected domains

polydemo(3); % speed comparison

INPOLY implements a sorted 'crossing-number' test designed to achieve fast performance for complex inputs. Given a configuration with N points and M polygon edges, INPOLY runs in approximately O((N+M)*LOG(N)) time on average, improving on the O(N*M) scaling of naive implementations.

Darren Engwirda (2021). INPOLY: A fast points-in-polygon test (https://github.com/dengwirda/inpoly), GitHub. Retrieved .

Created with
R2018b

Compatible with any release

**Inspired:**
FINDPOLY: A fast points-in-polygons test, The Barycentric Fixed-Mass method for estimating fractal dimensions, Maximum Inscribed Circle using Voronoi Diagram, Flow Cytometry GUI for Matlab, Maximum Inscribed Circle using Distance Transform, Fast Inpolygon in MEX

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!Create scripts with code, output, and formatted text in a single executable document.

Seras VictoriaErfanul HuqdkRueben MendelsbergNice work, much faster than built in function.

If using polygons with multiple regions separated by [nan nan], make sure to the last entry of the input point list is also [nan nan] or it will join the first and last point in the list (which are different regions and should not be closed)

Gerald KoudijsDarren Engwirda@Yue: If you have a large number of polygons (and points) the findpoly library might be of use, it is often much faster than calling inpoly in a loop:

https://www.mathworks.com/matlabcentral/fileexchange/74806-findpoly-a-fast-points-in-polygons-test

If you have any questions please feel free to discuss via GitHub, email, etc.

Yue HuangThanks Darren. I'll spend some time doing a bit more testing. I have a few polygons but they're all separated. I've got it working by running it for each of the polygons instead of all polygons together (not sure if this is some limitation I've omitted).

Thanks again for this fantastic work!

Darren Engwirda@Yue: it's hard to comment on the behaviour you're describing without an actual example --- if you believe you've found an error please create an "issue" via GitHub. Please also note that inpoly() does not support self-intersecting inputs like MATLAB's inpolygon.

Yue HuangI've tried almost every toolbox that finds points inside polygon and I'd say this is amongst the fastest. My only problem is the results are not fully correct, some parts of the shape fail with the points' distribution all inverted. This hasn't happened when I tried the same data with MATLAB's built-in inpolygon() or another File Exchange function inpolygons(), both of which are much slower. I'm not sure if it's the fast algorithm that causes loopholes

Arthur de JongWorks great! Had to detect >300K points in a polygon. With isinterior, Matlab consumed all available RAM and took forever. With inpoly, it runs in <1 sec. Thanks!

Eugenio GrabovicWondering if a gpu version would be achiveable for this kind of alghorithm.

chung-cheng leeExcellent function, thank you.

ha haPlease show me: How to use you INPOLY funnction? I try to follow your guide in https://github.com/dengwirda/inpoly. But i still can't know how to use it. Thanks

Steven KoxChad GreeneGreat function! I find it performs the test in about 10% of the time compared to inpolygon, and in all of my tests, inpoly2 and inpolygon give exactly the same results.

I typically use inpolygon to get logical masks of 2D data, so for usability I wrote this wrapper function, which lets me use the same input syntax as inpolygon:

function in = inpoly3(x,y,xv,yv,varargin)

in = reshape(inpoly2([x(:) y(:)],[xv yv],varargin{:}),size(x));

Kotya KarapetyanFaez AlkadiHi Darren,

What if I want to find points detect pints with in tolerance of (0.25) from the Edge?

How can I do that?

Thank you.

Lars Corbijn van WillenswaardDoes not correctly work with points on the boundary. For example

[in, on] = inpoly([1,1], [0,0; 0,1; 1,1; 1,0]) => in = 0, on = 0

[in, on] = inpoly([0,0], [0,0; 0,1; 1,1; 1,0]) => in = 1, on = 0

[in, on] = inpoly([0,0.5], [0,0; 0,1; 1,1; 1,0]) => in = 1, on = 0

[in, on] = inpoly([1,0.5], [0,0; 0,1; 1,1; 1,0]) => in = 0, on = 0

Behaviour of points on the boundary in the mulitpoint version depends on the the point before it:

[in, on] = inpoly([0.5, 0.5; 1,1], [0,0; 0,1; 1,1; 1,0]) -> in = [1,1], on = [0,1]

[in, on] = inpoly([1.5, 1; 1,1], [0,0; 0,1; 1,1; 1,0]) -> in = [0,0], on = [0,0]

[in, on] = inpoly([1.5, 0.5; 1,1], [0,0; 0,1; 1,1; 1,0]) -> in = [0,1], on = [0,1]

Lars Corbijn van WillenswaardNibodh BoddupallilowrieFabientic

[in,on] = inpolygon(xq, yq, M0(:,1), M0(:,2));

toc

"Elapsed time is 143.513056 seconds."

tic

[cn,on] = inpoly(horzcat(xq, yq),M0);

toc

"Elapsed time is 0.577500 seconds."

Many thanks \o/ \o/

John MorganHad some issues with detecting points on the polygon. Fixed by changing line 184 from

on(j) = on(j) || (abs((y2-Y)*(x1-X)-(y1-Y)*(x2-X))<tol);

to

on(j) = on(j) || (abs((y2-Y)*(x1-X)-(y1-Y)*(x2-X))<TOL);

As far as I can tell, in some cases tol was being set to zero rather than a sufficiently small number.

Daeyong KimAjayVery dependable code, thanks.

Any recommendations for a 3-D point-in-polygon test, including non-convex cases?

LucaHi,

running the following code, with N = 1000000, never gives inpoly faster than inpolygon.For "small" N (N=100-1000), they alternate.

However, when running the examples provided in polydemo.m, inpoly is consistently faster than inpolygon.

Am I doing anything wrong? A graphical check (commented in the code) gives what seems to be right, so I think I'm using the function as I should. Anyway, I don't see any improvement in speed, the contrary (for large N).

%%% sample code %%%%%%

for k=1:10

small_box = [ 1/3 1/3; 1/3 2/3; 2/3 2/3; 2/3 1/3];

N = 1000000;

A = rand(N,1);

B = rand(N,1);

% figure; hold on;

% plot(A,B,'.')

% patch([small_box(:,1); nan],[small_box(:,2); nan],ones(5,1))

tic; in = inpolygon(A,B,small_box(:,1),small_box(:,2)); t1=toc;

tic; in2 = inpoly([A B],small_box); t2=toc;

time_in_sec = struct('inpolygon',t1,'inpoly',t2);

if t2<t1

fprintf('inpoly is %f seconds faster\n',t1-t2)

else

fprintf('inpolygon is %f seconds faster\n',t2-t1)

end

% plot(A(in),B(in),'.r')

% plot(A(in2),B(in2),'.g')

end

%%%%%%%%%

Mario CASTRO GAMANice, used it in a problem with an unstructured grid and worked perfect.

JosephNot sure when matlab changed this but in matlab 2012b the inpolygon function works just slightly faster than this one. Still great work!

Edison LeeHi was wondering could you tell me from which conference/journal paper this function algorithm is based from?

Thanks.

Michelle TadmorThanks!!

Woody WongLuke WinslowGreat and fast little tool. Much faster than matlab inpolygon. My only issue is that I wish it natively understood the typical "GIS" format for polygons which includes a NaN separated list of polygons (NaNs separate the major outline from the 'islands').

Of course you can just handle this with the edges field though, so for future reference, here's my simple create edges code for NaN separated GIS objects. 'shp' variable is n Nx2 matrix of latitudes and longitudes.

shpEnd = find(isnan(shp(:,1)));

shpEnd = vertcat(0,shpEnd);

edges = nan(length(shp(:,1))-length(shpEnd),2);

count = 1;

for j=1:length(shpEnd)-1

endCount = count+length((shpEnd(j)+1:shpEnd(j+1)-2));

edges(count:endCount,:) = [(shpEnd(j)+1:shpEnd(j+1)-2)' ...

(shpEnd(j)+2:shpEnd(j+1)-1)';shpEnd(j+1)-1 shpEnd(j)+1];

count = endCount+1;

end

Svenindeed the on edge test is not robust, points on the edge will not be detected as on edge or as inside polygon. So watch out with that!

I am using the build in inpolygon routine now, works for me.

Bruno LuongThis is a quality code

Luigi GiaccariThank You,

it was very helpfull, impressive code.

If I can suggest an improvement:

a mex version will be the top of performance.

Dag LindboVery nice routine!

Is it reasonable to look for more efficient method for the particular case when the points to query at can be assumed to lie on a cartesian grid? E.g. imagine

in = inpoly(p,node);

with p = {X, Y}, where [X Y] = meshgrid(x,y)

Thanks!

Lili WanIt is a good job, but I still find some problems when detecting a point of a polygon lies on the polygon edge. My test run in Matlab R2006a. Suppose node be the points of a polygon,

"[in,on] = inpoly(node,node);" can get right results, but after running

"[in,on] = inpoly(node(1,:),node);", on is false.

Armin MüllerVery fast, very accurate. Good jobb, Darren!

Joseph MarksMy review below is too harsh.

After examining several runs on many different "ultra-concave" polygons, I have found inpoly to be very good.

If the error problem I reported earlier is real, it is very rare. It could have been due to some other factor in my software including my own bug.

inpoly is very much faster than inpolygon for large test vectors.

Joseph Marksinpoly provides a very large speed increase for large polygon problems!

Unfortunately, I too noticed a bug.

I am checking a very large rectangle grid's points to see if they are in a "massively concave polygon" -- think the outline of a "robot with arms, legs, etc."

The inpoly algorithm "incorrectly" *ADDS* a "shock of hair" to the "robot" -- obviously I am speaking metaphorically here -- and hence I am stuck using the *MUCH SLOWER* "inpolygon".

Has anyone reported an error like this to you?

Is there anything I can do?

Is there some middle ground between the two -- for example, by and large, I am dealing with concave polygons that just have an outside boundary -- no hole or anything.

Darren EngwirdaSmall bug (as noted below) fixed.

Users don't beware, inpoly should work in all cases.

Further bugs can be emailed if necessary...

Alex StorerVery fast, but fails on some cases. Users beware! Perhaps older versions are more robust?

Michael MIt's an incredible speed up compared with the poor inpolygon function delivered with MATLAB! Must be O(M*log(N)) ;-)

Matt K.Perfectly suited for my needs of finding points on the boundary of a polygon

Urs (us) SchwarzQUOTE

The reason that cnect is defined on its own is so that the domain can be multiply connected (polygon with "islands").

QUOTE

please, note that i said at run-time, i still feel that you should make the third arg optional

us

Darren Engwirda (The author)The reason that cnect is defined on its own is so that the domain can be multiply connected (polygon with "islands").

If you assume that cnect can be built by taking consecutive nodes this is no longer possible.

Thanks though, I will update it to flag boundary points as us mentioned.

John D'ErricoWhile I prefer the edge list implementation this code uses as opposed to Matlab's polygon, it would be easy enough to generate the connectivity assuming consecutive points on the polygon as us points out.

Regardless, this code is indeed fast and nice.

Urs (us) Schwarzvery nice (and indeed fast) snippet with a good help section and economic implementation of the crossing number test

minor comments:

- the third arg CNECT should be computed automatically (if not defined at run-time) on the assumption that the user-defined points are connected consecutively; this behavior could/should be mentioned in the help section

- unfortunately, unlike INPOLYGON it does not (yet) distinguish between IN and ON the polygon

- the help section should give a pointer to INPOLYGON

us