|
|
|
| Polygon Extreme Point Algorithms |
|
|
|
|
|
|
| Distance of a Polygon to a Line |
| Implementations |
|
|
|
|
| References |
It is often useful to find an extreme point of a 2D polygon. For example, the vertices with maximal or minimal x-coordinates or y-coordinates define a polygon's bounding box (see: Algorithm 8). More generally, one might want to find an extreme point in an arbitrary direction. For any set of n points, such as an n-vertex polygon, this is easily done in O(n) time by testing each point against the previously found extreme. However, for the special case of a convex polygon, an adaptation of binary searching can find the extreme point in only O(log n) time, as described by [O'Rourke, 1998]. We describe this algorithm in detail for an arbitrary query direction vector u. Also, we show how this algorithm can be applied to compute the distance from a convex polygon to a line.
We further note that the convex hull of a point set or polygon (see: Algorithm 10) is precisely the collection of extreme points in all possible directions. Thus, if several extreme point queries are expected for an arbitrary polygon, it may make sense to first compute its convex hull, and then do queries on this hull in O(log h) time, where h £ n is the number of hull vertices.
Let a 2D polygon W be given by n vertices V0,V1,...,Vn-1,Vn=V0 going counterclockwise (ccw) around the polygon. Also let ei be the i-th edge segment from vertex Vi to Vi+1 for i=0,n-1; and evi = (Vi+1-Vi) be the edge vector.
Next, let a direction be given by a vector u. We want to find the extreme, maximum and minimum, vertices of W in the direction u. That is, if the vertices Vi are orthogonally projected onto a line L in the direction u, then the extreme vertices are the ones whose projections are the extreme points on L, as shown in the following diagram.

We will say that Vi is above Vj relative to u, if its orthogonal projection onto L is further along in the u-positive direction. That is, if p: R2 ® L is the orthogonal projection to L, then the vector from p(Vj) to p(Vi) is in the same direction as u. However, this is equivalent to having the vector (Vi-Vj) make an acute angle with u, and this further corresponds to the condition that: u · (Vi-Vj) > 0. So, we use this algebraic expression to test vertices for the above-below relation in the u direction. Additionally, this test can determine if an edge ei of W is increasing or decreasing relative to u; namely, it is increasing if its edge vector evi satisfies: u · evi > 0, and decreasing when: u · evi < 0.
The straightforward brute-force way to find an extreme point is to test all points incrementally, and remember the current extreme for the points tested so far. Each new point considered only has to be compared to the current extreme one. This works for any set of n points, and is clearly an O(n) algorithm. The algorithm for finding the extreme maximum and minimum relative to u is simple, as shown in the following pseudo-code.
| Input: W
= {V0,V1,...,Vn-1}
is a set of n points u = a direction vector
Put: max = min = 0 |
For an arbitrary set of points without structure, this is the best one can do. For a set of n points, this algorithm will perform (n-1) dot product comparisons to compute the maximum of a set.
However, if the point set is the ordered set of vertices for a convex polygon, then we can find the maximum (or similarly, the minimum) relative to u more quickly, in O(log n) time, by using a binary search. But the algorithm is no longer simple. Say we know that the maximum vertex lies on the part of the polygon starting at vertex Va and ending at vertex Vb. This "chain" is given by the index range [a,b] = {a, a+1,..., a+k = b(mod n)} where k>0 is the first k for which the final equality holds. We next pick a vertex midway between Va and Vb, call it Vc. To do the binary search, we need a simple test that will restrict our search to one of the two subchains [a,c] or [c,b]. Because the polygon is convex, we can easily determine which one to pick by comparing the edge vectors eva at Va and evc at Vc. Six cases can occur, 3 each with A up and A down, as shown in the following diagram.

Then, for each case, one picks the segment, either [a,c] or [c,b], that would contain the maximum. These are the only cases that can occur because the polygon W is monotonic in the direction u. That is, it consists of exactly two monotone polyline segments, one increasing and one decreasing relative to the direction u, that join the max and min vertices. Since a convex polygon is monotonic in any direction, this is true no matter which direction u is selected. The converse is also true; namely, a polygon that is monotonic in all directions must be convex.
We can now construct the whole algorithm by starting with the chain [a=0,b=n], which is the whole polygon. Then, as long as b > a+1, we will always have that: 0 £ a < c = (a+b)/2 < b £ n and do not have to do any additional (mod n) computations. At each step, after testing Vc and determining which case we have, we either increase a to c or decrease b to c, and get a new range [a,b] that is half as large as the former one. Of course, we do have to check whether we have found a local maximum by testing each new subdivision vertex Vc. This is done by checking if Vc is greater than or equal to it's two neighboring vertices relative to u. This test works as long as the two neighbors do not both have equal heights with Vc. However, if that were true, then Vc would not be an essential proper vertex of W; instead, it would be on the interior of the segment from Vc-1 to Vc+1, and so could be removed without changing W's boundary. A convex polygon with only proper vertices is called a proper convex polygon. For any proper convex polygon, a local maximum must also be a global maximum, and we are done when we find one.
Putting this together, the complete algorithm to find a u-maximum vertex is given by the pseudo-code:
| Input: W
= {V0,V1,...,Vn-1,Vn=V0}
is a proper convex polygon u = a 2D direction vector
if V0 is a local maximum, then return 0; |
A C++ implementation is given below in the routine polyMax_2D(). As implemented, we use 1 dot product comparison to test if a vector is up or down, and two dot product tests for a local maximum. For a polygon with n vertices, the algorithm uses 2.5 such tests on average for each loop iteration, and thus has 2.5 log(n) dot product comparisons in total. For small n, this is less efficient than the brute force algorithm as long as (n-1) < 3 log(n). The break even point is at n=9, where both methods use 8 comparisons (although the binary search may terminate faster when a maximum is found early). Since many frequently occurring convex polygons have fewer than 10 vertices, it is wise to use the brute force method for these small polygons.
We can use the above algorithm to quickly find the minimum distance from a proper convex polygon W = {Vi} to a line L through two points P0 and P1. Denote this distance by d(W, L). If the polygon is not convex, we can compute and use its convex hull if we expect to compute its distance to many different lines at different orientations.
The idea here is simple. First, pick u to be perpendicular normal vector to the line L so that at least one vertex of W is on the side of L not pointed to by u. More precisely, a point P is on the side of L pointed to by u, called the u-positive side of L, if u · (P - P0) > 0. Otherwise it is on the u-backside. To start, we get a normal vector n of L by using the perp-operator (see [Hill, 1994] or Algorithm 4) to select n = (P1 - P0)^. Now, if n · (V0 - P0) £ 0, V0 is on the n-backside of L, and we put u = n. If not, we simply reverse this normal vector by putting u = -n. In either case, u is normal to L and V0 is on the u-backside.
Next, we find the maximum vertex Vmax of W relative to this normal vector u. There are two cases:
Vmax is on the u-positive side of L, and thus any polyline chain from V0 to Vmax must cross L. This means that d(W, L) = 0 since they have a common point.
Vmax is on the u-backside of L, and thus all vertices of are on L's u-backside. Then, Vmax is the closest vertex to L, and d(W, L) = d(Vmax, L) the distance from Vmax to the line (see Algorithm 2 about the Distance of a Point to a Line).
These cases are shown in the following diagram.

C++ code for this algorithm is given below in the routine dist2D_Poly_to_Line().
Here are some sample "C++" implementations of these algorithms.
// Copyright 2001, softSurfer (www.softsurfer.com)
// This code may be freely used and modified for any purpose
// providing that this copyright notice is included with it.
// SoftSurfer makes no warranty for this code, and cannot be held
// liable for any real or imagined damage resulting from its use.
// Users of this code must verify correctness for their application.
// Assume that classes are
already given for the objects:
// Point and Vector (2D) with:
// coordinates {float x, y;}
// operators for:
// == to test
equality
// != to test
inequality
// = for
assignment
// -Vector for
unary minus
// Point
= Point ± Vector
// Vector =
Point - Point
// Vector =
Vector ± Vector
// Line with defining points {Point P0, P1;}
//===================================================================
// dot product (2D) which allows vector operations in arguments
#define dot(u,v) ((u).x * (v).x + (u).y * (v).y)
// tests for vector orientation relative to a direction vector u
#define up(u,v) (dot(u,v) > 0)
#define down(u,v) (dot(u,v) < 0)
#define dr(u,Vi,Vj) (dot(u, (Vi)-(Vj))) // direction
sign of (Vi-Vj)
#define above(u,Vi,Vj) (dr(u,Vi,Vj) > 0) // true if Vi is
above Vj
#define below(u,Vi,Vj) (dr(u,Vi,Vj) < 0) // true if Vi is
below Vj
// polyMax_2D(): find a
polygon's max vertex in a specific direction
// Input: U = a 2D direction vector
// V[] = array
vertices of a proper convex polygon
// n
= number of polygon vertices, with V[n]=V[0]
// Return: index (>=0) of the maximum vertex, or
// (-1) = an
error [Note: should be impossible, but...]
int
polyMax_2D( Vector U, Point* V, int n )
{
if (n < 10) {
// use brute force search for small polygons
int max = 0;
for (int i=1; i<n; i++)
// for each point in {V1,...,Vn-1}
if (above(U,V[i],V[max]))
// if V[i] is above prior V[max]
max = i;
// new max index = i
return max;
}
// use binary search for large polygons
int a, b, c;
// indices for edge chain endpoints
Vector A, C;
// edge vectors at V[a] and V[c]
int upA, upC;
// test for "up" direction of A and C
a=0; b=n;
// start chain = [0,n] with V[n]=V[0]
A = V[1] - V[0];
upA = up(U,A);
// test if V[0] is a local maximum
if (!upA && !above(U, V[n-1], V[0])) //
V[0] is the maximum
return 0;
for(;;) {
c = (a + b) / 2;
// midpoint of [a,b], and 0<c<n
C = V[c+1] - V[c];
upC = up(U,C);
if (!upC && !above(U,V[c-1],V[c]))
// V[c] is a local maximum
return c;
// thus it is the maximum
// no max yet, so continue with the
binary search
// pick one of the two subchains [a,c]
or [c,b]
if (upA) {
// A points up
if (!upC) {
// C points down
b = c;
// select [a,c]
}
else {
// C points up
if (above(U,V[a],V[c])) { // V[a] above V[c]
b = c;
// select [a,c]
}
else {
// V[a] below V[c]
a = c;
// select [c,b]
A = C;
upA = upC;
}
}
}
else {
// A points down
if (upC) {
// C points up
a = c;
// select [c,b]
A = C;
upA = upC;
}
else {
// C points down
if (below(U,V[a],V[c])) { // V[a] below V[c]
b = c;
// select [a,c]
}
else {
// V[a] above V[c]
a = c;
// select [c,b]
A = C;
upA = upC;
}
}
}
// have a new (reduced) chain [a,b]
if (b <= a+1)
// the chain is impossibly small
return (-1);
// return an error: something's wrong
}
}
//===================================================================
//
dist2D_Poly_to_Line(): find the distance from a polygon to a line
// Input: V[] = array vertices of a proper convex
polygon
// n
= the number of polygon vertices, with V[n]=V[0]
// L
= a Line (defined by 2 points P0 and P1)
// Return: minimum distance from V[] to L
float
dist2D_Poly_to_Line( Point* V, int n, Line L )
{
Vector U, N;
int max;
// get a leftward normal N to L
N.x = -(L.P1.y - L.P0.y);
N.y = (L.P1.x - L.P0.x);
// get a normal U to L with V[0] on U-backside
if (dot(N, V[0]-L.P0) <= 0)
U = N;
else U = -N;
max = polyMax_2D( U, V, n );
// max vertex in U direction
if (dot(U, V[max]-L.P0) > 0)
// V[max] on U-positive side
return 0;
// so polygon and line intersect
else
return
dist_Point_to_Line( V[max], L); // min dist to line
}
//===================================================================
Francis Hill, "The Pleasures of 'Perp
Dot' Products" in
Graphics Gems IV (1994)
[Note: the first critical definition has a typo, and should be: a^
= (-ay, ax).]
Joseph O'Rourke, Computational Geometry in C (2nd Edition) , Sect 7.9 "Extreme Point of Convex Polygon" (1998)
|
Copyright ©
2001-2006 softSurfer. All rights reserved. |