Relay-Version: version B 2.10 5/3/83; site utzoo.UUCP
Path: utzoo!watmath!clyde!caip!think!mit-eddie!genrad!decvax!ucbvax!jade!eris!mwm
From: mwm@eris.berkeley.edu (Mike (Don't have strength to leave) Meyer)
Newsgroups: sci.math
Subject: Re: point inside a polygon (triangle)
Message-ID: <1482@jade.BERKELEY.EDU>
Date: Mon, 20-Oct-86 03:57:46 EDT
Article-I.D.: jade.1482
Posted: Mon Oct 20 03:57:46 1986
Date-Received: Tue, 21-Oct-86 22:51:05 EDT
References: <8900031@osiris>
Sender: usenet@jade.BERKELEY.EDU
Reply-To: mwm@eris.UUCP (Mike (Don't have strength to leave) Meyer)
Organization: Missionaria Phonibalonica
Lines: 88

In article <8900031@osiris> kaden@osiris.CSO.UIUC.EDU writes:
>
>		Can we generalise to determine whether (p1,p2) is
>		inside a polygon? Use the VECTOR approach, ONLY, as
>		it is clearly the most superior, elegant and non-
>		paleozoic approach to the problem. I will not even
>		scan any other precambrian attempts.
>
>
>						Thankyou.


I'm not sure what you mean by the "VECTOR" approach, but my solution
to the triangle problem (which I MAILED to the asker) generalizes
quite nicely to an arbitrary polygon. It's also generally faster than
any of the solutions I've seen posted.

The heart of the solution is a theorom which states that any ray from
a point crosses a closed curve an odd number of times iff the point is
inside the curve. To use this, we chose a computationally "nice" ray,
the one in the positive X direction from the point. The heart of the
algorithm is the routine crosses, below.

/*
 * Just to make things easy
 */
typedef struct	{
	float x, y ;
	} point ;

/*
 * Crosses determines if a ray in the positve x direction from p crosses
 * the line from v1 to v2. Start by running lots of tests to avoid
 * interpolation.
 */
crosses(p, v1, v2) point p, v1, v2; {

	if (p . y <= v1 . y && p . y < v2 . y)	/* to far in -y direction */
		return 0 ;
	if (p . y >= v1 . y && p . y > v2 . y)	/* to far in +y direction */
		return 0 ;
	if (p . x >= v1 . x && p . x > v2 . x)	/* to far in +x direction */
		return 0 ;
	if (p . x <= v1 . x && p . x < v2 . x)	/* Hit */
		return 1 ;
	/*
	 * We now have p in the box defined by [v1, v2). This means
	 * that 1) the box isn't a line, so v1 . y != v2 . y; and
	 * 2) we have to interpolate. To bad.
	 */
	if ((p . y - v1 . y) * (v1 . x - v2 . x) / (v1 . y - v2 . y) > 0)
		return 1 ;
	return 0 ;
	}

Now, let's define an N-gon by an array of points v[1] through v[N]
which are the vertices of the N-gon as you follow the boundary
(clockwise, counterclockwise, who cares?). Define v[0] equal to v[N]
to make my life easier. The the following code will test if p is
inside of v.

inside(p, v, n) point p, *v; unsigned n; {
	register unsigned i, crossings = 0;

	for (i = 0; i < n; i++)
		crossings += crosses(p, v[i], v[i + 1]) ;
	return crossings & 1 ;
	}

which returns true if crossings is odd, indicating that p is inside of
v.

Just to let you exercise a little, I'll point out that there's a bug
in the above code. There are lots of ways to fix it, and I still
haven't decided on the best. I'd like to see other peoples fixes
before posting anything.

Oh, yeah - the triangle version of inside is particularly elegant, being:

inside(p, v1, v2, v3) point p, v1, v2, v3; {

	return crosses(p, v1, v2) ^ crosses(p, v2, v3) ^ crosses(p, v3, v1) ;
	}

Of course, by adding some assumptions about the n-gon (triangle), you
can elimate the bug, and allow some nice speed improvements.

	<mike