# How to check if a given point is inside some closed path in MetaPost

If there is some arbitrary closed path (without self intersections) in MetaPost, is there a simple way to check if a given point falls inside or outside of it?

For example:

path p;
pair q, r;
p := (-3cm, 0) .. (-4cm, 1cm) .. (2cm, 3cm) .. (3cm, 1cm) .. (-1cm, -1cm) .. (1cm, -2cm) .. (-2cm, -1cm) .. cycle;
q := (0, 0);
r := (0, 1cm);
draw p;
fill (fullcircle scaled 1mm) shifted q withcolor red;
fill (fullcircle scaled 1mm) shifted r withcolor green;

In this case I want to get 'false' for 'q', the red one, and 'true' for 'r', the green one regarding their relation to 'p'.

• One ad hoc poposal: take one point of the path (or better a point inside the contour), connect it through a straight line (say) with the point in question, and then use this method to count the intersections of the straight line with the path. If it is even, the point is inside the contour, otherwise it is outside.
– user121799
Mar 13, 2018 at 0:21
• If the point is on path, this might not work. For example: imagine a point on the "north" of the "southern peninsula" on the picture above; the number of intersections will be even for 'r' and odd for 'q'. And to find a point inside is a problem on its own (shouldn't be too difficult though). And even if the same thing is done the other way around (as Marcel Krüger have suggested below), intersections in MP do not seem to be a very reliable tool for such a task. Mar 13, 2018 at 7:49

Based on Get all intersection points in MetaPost and implementing the even-odd rule:

vardef is_inside(expr p,q) =
begingroup;
save cut, inside;
path cut;boolean inside;
cut := q -- (ulcorner p + (-2,2));
inside := false;
forever:
cut := cut cutbefore p;
exitif length cuttings = 0;
cut := subpath(4epsilon, length cut) of cut;
inside := not inside;
endfor;
inside
endgroup
enddef;

beginfig(1) ;
path p;
p := (-0.5cm,0cm) .. (5cm,5cm) .. (0cm,3cm) .. (-7cm,9cm) .. cycle;
z0 = (2cm,3.5cm);
draw p;
draw z0 withpen pencircle scaled 2;
draw origin withpen pencircle scaled 1;

if is_inside(p,z0):
label.lrt("IN", origin);
else:
label.lrt("OUT", origin);
fi
endfig ;

An alternative rule is the non-zero winding rule. This is an implementation of windingnumber for MetaPost by Bogusław Jackowski. I extended the inside operator to allow testing if a point is inside a path:

vardef mock_arclength(expr p) = % |p| -- Bézier segment
% |mock_arclength(p)>=arclength(p)|
length((postcontrol 0 of p)-(point 0 of p)) +
length((precontrol 1 of p)-(postcontrol 0 of p)) +
length((point 1 of p)-(precontrol 1 of p))
enddef;
vardef windingangle(expr p,q) = % |p| -- point, |q| -- Bézier segment
save a,b,v;
a=length(p-point 0 of q); b=length(p-point 1 of q);
if min(a,b)<2eps: % MP is not the master of precision, we’d better stop now
errhelp "It is rather not advisable to continue. Will return 0.";
errmessage "windingangle: point unsafely near Bézier segment (dist="
& decimal(min(a,b)) & ")";
0
else:
v:=mock_arclength(q); % |v| denotes both length and angle
if (v>=a) and (v>=b): % possibly too long Bézier arc
windingangle(p, subpath (0,1/2) of q)+windingangle(p, subpath (1/2,1) of q)
else:
v:=angle((point 1 of q)-p)-angle((point 0 of q)-p);
if v>180: v:=v-360; fi
if v<-180: v:=v+360; fi
v
fi
fi
enddef;
vardef windingnumber (expr p,q) = % |p| -- point, |q| -- Bézier spline
save a; a:=0;
for t:=1 upto length(q):
a:=a+windingangle(p, subpath(t-1,t) of q);
endfor
a/360
enddef;
tertiarydef a inside b =
if path a: % |and path b|; |a| and |b| must not touch each other
begingroup
save a_,b_; (a_,b_)=
(windingnumber(point 0 of a,b), windingnumber(point 0 of b,a));
(abs(a_-1)<eps) and (abs(b_)<eps)
endgroup
elseif pair a: % |and path b|; |a| must not lie on |b|
begingroup
abs(windingnumber(a,b))>eps
endgroup
else: % |numeric a and pair b|
begingroup
(a>=xpart b) and (a<=ypart b)
endgroup
fi
enddef;

This is based on a topological invariant, so it should be more reliable. At least it should fail if it can't compute the right result.

• This solution as it is doesn't seem to be very reliable. Take, for example, these path and point: p := (0, 0) -- (0, 10cm) -- (10cm, 10cm) -- (10cm, 0cm) -- cycle; z0 = (5cm, 5cm); 'z0' is clearly inside, but the output is 'OUT' (it seems to have something to do with passing through angles). or here: p := (0, 0) -- (0, 10cm) -- (1cm, 10cm) -- (4cm, 6cm) -- (10cm, 6cm) -- (10cm, 0cm) -- cycle; z0 = (5cm, 7cm); 'z0' is outside, but the output is 'IN'. Plus, it feels that there could be some problems if the cutting line is tangent to something. Mar 13, 2018 at 7:29
• @SergeySlyusarev I added another solution which should be more reliable Mar 13, 2018 at 10:37