Algorithm 13
[Home] [Overview] [History] [Algorithms] [Books] [Web Sites] [Gift Shop]


Extreme Points of Convex Polygons and
Distance of a Polygon to a Line

by Dan Sunday

Polygon Extreme Point Algorithms
Brute-Force Search
Binary Search of a Convex Polygon
Pseudo-Code
Distance of a Polygon to a Line
Implementations
polyMax_2D()
dist2D_Poly_to_Line()
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.


Polygon Extreme Point Algorithms

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.

Brute-Force Search

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 each point Vi in {V1,...,V
n-1} (i=1,n-1)
{
    if (u · (Vi - Vmax) > 0) { // Vi is above the prior max
        max = i;               // new max index = i
        continue;              // a max can't be a min
    }
    if (u · (Vi - Vmin) < 0)   // Vi is below the prior min
        min = i;               // new min index = i
}

return max = index of max vertex

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.

Binary Search

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.

Pseudo-Code: Binary Search

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;
Put a=0; b=n;
Put A = the edge vector at V0;
forever {
    Put c = the midpoint of [a,b] = (a+b)/2;
    if Vc is a local maximum, then it is the maximum
        return c;

    // no max yet, so continue with the binary search
    // pick one of the two subchains [a,c] or [c,b]
    Put C = the edge vector at Vc;
    if (A points up) {
        if (C points down) {
            b = c;                    select [a,c]
        }
        else C points up {
            if Va is above Vc {
                b = c;                select [a,c]
            }
            else Va is below Vc {
                a = c;                select [c,b]
                A = C;
            }
        }
    }
    else A points down {
        if (C points up) {
            a = c;                    select [c,b]
            A = C;
        }
        else C points down {
            if Va is below Vc {
                b = c;                select [a,c]
            }
            else Va is above Vc {
                a = c;                select [c,b]
                A = C;
            }
        }
    }
    if (b <= a+1) then the chain is impossibly small
        return an error; since something's wrong
}

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.


Distance of a Polygon to a Line

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:

  1. 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.

  2. 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().


Implementations

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,...,V
n-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
}
//===================================================================


References

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.
Email comments and suggestions to
feedback@softsurfer.com