Friday, February 26, 2010

Convex Hull



In computational geometry, Convex Hull or Convex Polygon is a very common term and we face many problems which require the construction of a Convex Hull to generate desired solution.

A convex hull is a convex polygon which has the minimal area to contain the set of given points. In other words, for a given set of points, a Convex Hull is such a Convex Polygon that, every point on the set is either on the Polygon or inside the Polygon.

So, a typical problem might ask, "You are given n integer co-ordinates (3 ≤ n ≤ 105), find the Convex Hull". It means it wants you to write down the co-ordinates of the picks or, find the convex area, or, find the center of gravity of it, or, sometimes, you are asked to create an arbitrary convex hull, or you may be asked, find whether a given point is inside a convex hull, and so on... Well, before getting into the harder ones, lets see how we can construct a convex hull smartly :P

There are a number of ways to compute convex hull from the given points. Like, I found those in wikipedia, (I know only three of them):
  • Direct method
  • Gift Wrapping (Packaging) method
  • Graham's scan method
  • Incremental (Sequential addition) method
  • Divide and Conquer method
  • Quick method (related to Quick Sort)
  • Inner point elimination

The graham scan method is the most widely used method that has a complexity of O( n lg n ) for n points. In this method, the algorithm starts from a extreme (topmost / leftmost / rightmost / bottommost) point and keep wrapping all the other points unless it reaches the beginning point. Here's a little animation that shows it:


For detailed reading, you can visit Wikipedia or read Introduction to Algorithm - CLRS. Here's my implementation for finding Convex Hull with the help of Graham Scan method (integer points only):
#include <algorithm>
using namespace std;

#define MAX 100009
#define i64 long long
#define sq(x) ((x)*(x))

struct point {
    i64 x, y;
} P[MAX], C[MAX], P0;

// P[] contains all points
// C[] contains convex hull ccw

inline i64 TriArea2(point a, point b, point c) {
    return (a.x*(b.y-c.y) + b.x*(c.y-a.y) + c.x*(a.y-b.y));
}

inline i64 sqDist(point a, point b) {
    return (sq(a.x-b.x) + sq(a.y-b.y));
}

bool comp(point a, point b) {
    i64 d = TriArea2(P0, a, b);
    if(d<0) return false;
    if(!d && sqDist(P0, b) > sqDist(P0, a)) return false;
    return true;
}

void ConvexHull(int np, int &nc) {
    int i, j, pos = 0;
    for(i=1; i<np; i++)
        if(P[i].y<P[pos].y || (P[i].y==P[pos].y && P[i].x>P[pos].x))
            pos = i;
    swap(P[0], P[pos]);
    P0 = P[0];
    sort(&P[1], P+np, comp);
    C[0] = P[0], C[1] = P[1], C[2] = P[2];
    for(i=j=3; i<np; i++) {
        while(TriArea2(C[j-2], C[j-1], P[i]) <= 0) j--;
        C[j++] = P[i];
    }
    nc = j;
}

Here,
TriArea2() returns the double of triangular area formed by the three points passed as its argument.
sqDist() returns the square distance between two points passed as its argument.
sort() is stl sort that sorts the P[] array according to their slope.
The first loop in ConvexHull() finds the starting point, i.e. the bottommost and rightmost point. After sorting rest of the points correctly, the second loop examines which point is to be added to the Hull. Finally nc is the number of points in convex hull updated via reference. There is no checking for n < 3 or if there are duplicates (actually if number of distinct points are ≥ 3, it doesn't matter in this method), we need to handle that according to problem requirement.

Not much tough as it seems to be ... happy coding

8 comments:

  1. Not much tough as it seems to be ...

    Because it is much more tough than it seems to be...

    ReplyDelete
  2. Replies
    1. Thanks for pointing it out, it will be &nc, not &cp, it tells you the number of nodes in convex hull, this is a rough implementation, check for a better one on my code library:
      http://zobayer2009.wordpress.com/code

      Delete
  3. do u have an implementation of divide and conquer method as well .. that will be really helpful...

    ReplyDelete
    Replies
    1. please can you code one and post it,it will definitely help a lot.

      Delete
    2. do you mean the divide and conquer method?

      Delete