Page 7
Distance between a point and a line
Page 8
Index
See this article in english 
 

Convex Hull (Quick Hull)

Alexander Hristov

The convex hull problem in geometry tries to find the smallest convex set containing the points, or more informally, finding those 'border' points that, when united, form a convex polygon that includes all other points.

There are many approaches for handling this problem:

Here we'll talk about the Quick Hull algorithm, which is one of the easiest to implement and has a reasonable expected running time of O(n log n). (Unfortunately,it has a worst case scenario of O(n2) )

Let's assume a random set of points:

Quick Hull 0

As a start for the algorithm, we need a chord that we know as belonging to the final convex hull. One such chord can be the segment between the leftmost and the rightmost points, which are known to be on the hull:

Quick Hull

This chord will be our initial "set of hull points".

We now proceed to divide the set of the remaining points into two sets, depending on which side of the chord they fall. Let's call S1 the set of the points above the chord in our example, and S2 the set of the points below the chord. For each of the sets, we proceed recursively as follows:

1) We find the point that is furthest away from the chord. Let's call this point P. This point will belong to the hull, so we insert it in our "set of hull points" between A and B:

Quick Hull

2) We now remove from the S1 set all points that fall inside the triangle:

Quick Hull

3) Again, we form two sets : The set of points above AP and the set of points above PB (in the example above, this set is empty). We replace the original chord AB with the segments AP and PB, and apply the algorithm from step 1 recursively on the two sets. We stop when a set is empty or has only one point, in which case it's known to be on the hull

The remaining steps would be as follows:

Quick Hull Step   Quick Hull Step  Quick Hull Step

 

Quick Hull Step  Quick Hull Step  Quick Hull Step

 

Ok, if you are new to geometry algorithms, you might wonder how to check whethera point lies on one side of a segment or on the other side. The notion of vector cross product comes to help us. The cross product v x w of two vectors v and w is a vector whose length is |v||w|sin φ, (where |v| is the length of v and φ is the angle between the vectors), and who is orthogonal (perpendicular) to both v and w.  Since there are two such possible vectors, the definition arbitrarily selects the one that matches the direction in which a screw would move if rotated from v to w:

Cross Product

 

Mathematically, if the coordinates of vectors v and w are (vx,vy) and (wx,wy) respectively, the cross produt will be (vxwy - wywx). So now, if you have a segment defined by points A B and want to check on which side of AB a third point P falls, you only need to compute the cross product AB x AP and check its sign: If it's negative, it will be on the "right" side of AB (when standing on A and looking towards B). If positive, it will be on the left side:

Cross Product

This Java function performs the check for two points:

 
  public int pointLocation(Point A, Point B, Point P) {
    int cp1 = (B.x-A.x)*(P.y-A.y) - (B.y-A.y)*(P.x-A.x);
    return (cp1>0)?1:-1;
  }
  
 

The second thing you might find trouble with is computing the distance between a point and a segment. Check the related article here, and consider that when we are trying to find -in step 1 of the algorithm - the point P that is farthest away from AB, the segment AB never changes. And since we need only to know the biggest distance and not its exact value, we need not divide by the length of the AB segment, and we needn't take the square root either. Thus, we only need to be concerned by the following "pseudodistance" when looking for the farthest point:

 

 
  public int distance(Point A, Point B, Point C) {
    int ABx = B.x-A.x;
    int ABy = B.y-A.y;
    int num = ABx*(A.y-C.y)-ABy*(A.x-C.x);
    if (num < 0) num = -num;
    return num;
  }  

 
With the above two functions, the hull algorithm described above can be implemented as follows:
 
0001  public ArrayList<Point> quickHull(ArrayList<Point> points) {
0002    ArrayList<Point> convexHull = new ArrayList<Point>();
0003    if (points.size() < 3) return (ArrayList)points.clone();
0004    // find extremals
0005    int minPoint = -1, maxPoint = -1;
0006    int minX = Integer.MAX_VALUE;
0007    int maxX = Integer.MIN_VALUE;
0008    for (int i = 0; i < points.size(); i++) {
0009      if (points.get(i).x < minX) {
0010        minX = points.get(i).x;
0011        minPoint = i;
0012      } 
0013      if (points.get(i).x > maxX) {
0014        maxX = points.get(i).x;
0015        maxPoint = i;       
0016      }
0017    }
0018    Point A = points.get(minPoint);
0019    Point B = points.get(maxPoint);
0020    convexHull.add(A);
0021    convexHull.add(B);
0022    points.remove(A);
0023    points.remove(B);
0024    
0025    ArrayList<Point> leftSet = new ArrayList<Point>();
0026    ArrayList<Point> rightSet = new ArrayList<Point>();
0027    
0028    for (int i = 0; i < points.size(); i++) {
0029      Point p = points.get(i);
0030      if (pointLocation(A,B,p) == -1)
0031        leftSet.add(p);
0032      else
0033        rightSet.add(p);
0034    }
0035    hullSet(A,B,rightSet,convexHull);
0036    hullSet(B,A,leftSet,convexHull);
0037    
0038    return convexHull;
0039  }
0040  
0041  public void hullSet(Point A, Point B, ArrayList<Point> set, ArrayList<Point> hull) {
0042    int insertPosition = hull.indexOf(B);
0043    if (set.size() == 0) return;
0044    if (set.size() == 1) {
0045      Point p = set.get(0);
0046      set.remove(p);
0047      hull.add(insertPosition,p);
0048      return;
0049    }
0050    int dist = Integer.MIN_VALUE;
0051    int furthestPoint = -1;
0052    for (int i = 0; i < set.size(); i++) {
0053      Point p = set.get(i);
0054      int distance  = distance(A,B,p);
0055      if (distance > dist) {
0056        dist = distance;
0057        furthestPoint = i;
0058      }
0059    }
0060    Point P = set.get(furthestPoint);
0061    set.remove(furthestPoint);
0062    hull.add(insertPosition,P);
0063    
0064    // Determine who's to the left of AP
0065    ArrayList<Point> leftSetAP = new ArrayList<Point>();
0066    for (int i = 0; i < set.size(); i++) {
0067      Point M = set.get(i);
0068      if (pointLocation(A,P,M)==1) {
0069        //set.remove(M);
0070        leftSetAP.add(M);
0071      }
0072    }
0073    
0074    // Determine who's to the left of PB
0075    ArrayList<Point> leftSetPB = new ArrayList<Point>();
0076    for (int i = 0; i < set.size(); i++) {
0077      Point M = set.get(i);
0078      if (pointLocation(P,B,M)==1) {
0079        //set.remove(M);
0080        leftSetPB.add(M);
0081      }
0082    }
0083    hullSet(A,P,leftSetAP,hull);
0084    hullSet(P,B,leftSetPB,hull);
0085  }
0086
 

You can try this implementation interactively using the following applet. Click anywhere to drop a point. Once placed, you can drag it with the mouse.



As always, you can download below the source code of the QuickHull function alone, or the source code of the above applet if you wish to play with it (GPL license).

Source code



 

Comments

Dec 08, 2010 at 04:05 Sent by Dominique
Hello, I have a problem...try with this points set 72641 35759 -12210 19702 -110318 35567 76567 35494 75322 35912 -110085 35657 -104208 36232 53360 25281 1347 18235 -1304 17777 -3954 17827 -110031 34305 -20862 17532 -55338 18035 76255 33058 74621 33052 -111003 33098 75721 32788 -107874 27375 -109499 32128 75597 31997 -110005 32213 76913 32412 -108690 32250 -110541 32421 76652 34744 74933 35028 -110706 34922 -24823 20087 seems it doesn't work...
Aug 09, 2010 at 13:31 Sent by anonymous
how can i modifiy to show the points in trigonometrical order? infoarena.ro/job_detail/475962?action=view-source
Nov 23, 2009 at 20:33 Sent by Petar
I think you have a mistake in the textual explanation: "If it's negative, it will be on the "right" side of AB (when standing on A and looking towards B). If positive, it will be on the left side ". However, in the code: 0030 if (pointLocation(A,B,p) == -1) 0031 leftSet.add(p); 0032 else 0033 rightSet.add(p); 0034 } Pozdravi
Nov 23, 2009 at 14:25 Sent by Juliane
first of all: thanks for your great work second: your calculation for distance(...) and pointlocation(...) is basically the same thing (the cross product)- stating this assumtion is correct, that means, while sorting points to right and left side of AB, you could determine the most distant point P in this content as well - no need to go two times through all points - or am I wrong?
Nov 11, 2009 at 15:55 Sent by anonymous
Please note this small correction: the cross produt will be (vxwy - vywx).
Jun 17, 2009 at 23:12 Sent by anonymous
Lighter than qhull, but powerful even so ! Thanks a lot.
Apr 22, 2009 at 06:16 Sent by Diego
i translate to c++, and i had some problems some values duplicate... i use vector...
Dec 28, 2008 at 19:40 Sent by Carlos
The line set.remove(p); when size==1 can be safely removed.
Dec 28, 2008 at 19:06 Sent by Carlos
Thanks for this clarifying piece of code!
Oct 01, 2008 at 22:46 Sent by Samuel Félix
Wow, that is a good work, contratulation =D

 

Add a Comment

Name (optional)
EMail (optional, will not be displayed)

Text