Convex hulls in Python: the Graham scan algorithm

convex_hull.png

The boundary of the smallest convex polygon that encloses all of the points in a set makes up the convex hull.

This is the Graham scan algorithm in action, which is one common algorithm for computing the convex hull in 2 dimensions.

The animation was created with Matplotlib.

Computing the convex hull is a preprocessing step to many geometric algorithms and is the most important elementary problem in computational geometry, according to Steven Skiena in the Algorithm Design Manual.

Graham scan algorithm

This algorithm is pretty straightforward to learn.

You may have heard that you can use sorting to find a convex hull and wondered how and where sorting would come into play.

The basic idea

  • Initialize an empty stack that will contain the convex hull points.
  • Pick a starting point and add it to the stack.
  • Sort the rest of the points in counterclockwise order around the starting point.
  • Sweep through the sorted points.
  • Initially add each new point to the stack and then check to make sure that the hull is still convex with the new point.
  • Delete any points that create concave angles - these points lie inside of the hull.

Python implementation of the Graham scan algorithm

The full code can be found here.

We will compute the convex hull of a set of 50 random points in a 100 x 100 grid.

points = [(random.randint(0,100),random.randint(0,100)) for i in range(50)]

Initialize an empty stack - I'm using a Python list for the stack.

hull = []

Pick a starting point that will definitely be part of the convex hull.

Since the convex hull is made up of the outer boundary of the set of points, it makes sense that a minimum or maximum point would be part of it.

In a 2-D set of points you could use the point with the smallest x value.

If more than one point has the smallest x value, take the one with the smallest y value.

points.sort(key=lambda x:[x[0],x[1]])
start = points.pop(0)

Next, sort the rest of the points in counterclockwise order, which we can do by sorting them by their slope.

def get_slope(p1, p2):
    if p1[0] == p2[0]:
        return float('inf')
    else:
        return 1.0*(p1[1]-p2[1])/(p1[0]-p2[0])

points.sort(key=lambda p: (get_slope(p,start), -p[1],p[0]))

Now we can sweep through the points in a counterclockwise direction.

Add each point to the hull stack initially, and then we check to make sure the three points making up each new corner of the polygon create a convex angle.

Once there are at least 3 points in the stack, we start looking at each triplet of points from the current point to the previous two points in the stack hull[-3:].

Call them p1, p2 and p3, with p3 being the current point we're looking at.

The points create two vectors p1 -> p2 and p2 -> p3.

Since we are sweeping in a counterclockwise direction, we want to rotate left to get to the current point from the previous two points.

convexangle.png

Rotating left means that this corner of the convex hull polygon we are forming is, indeed, convex, and as we know, all of the corners of the convex hull need to be convex.

We can find out the rotation direction by computing the cross product of the vectors created by p1, p2 and p3.

def get_cross_product(p1,p2,p3):
    return ((p2[0] - p1[0])*(p3[1] - p1[1])) - ((p2[1] - p1[1])*(p3[0] - p1[0]))

If the result is negative, then the three points are rotating right in a clockwise direction, which would add a concave angle to the polygon, so we want to get rid of the second point, p2, because it lies inside the convex hull.

concaveangle.png

After removing P2 from the stack, check to make sure the new triplet in hull[-3:] rotates left - if it still rotates right, keep removing the middle point until the triplet rotates left.

while len(hull) > 2 and get_cross_product(hull[-3],hull[-2],hull[-1]) < 0:
    hull.pop(-2)

If the cross product is zero, meaning the points are in a straight line or collinear, it's up to you and your project requirements whether or not you want to keep or drop the point.

Then move on to the next point in the sweep.

When the sweep is done, the points that remain in hull are the points that form the convex hull.

Computing the convex hull in higher dimensions

The gift wrapping algorithm is typically used for finding the convex hull in a higher dimensional space.

In the 2-D case, this algorithm is known as the Jarvis march.

Python libraries

Here are a few options for computing convex hulls in your projects.

Let me know of any other libraries you know of!

A few creative applications of convex hulls

  • Face-swapping apps. This face-swap tutorial for OpenCV, shows how you can use the convex hull of various facial landmarks in images to determine the face boundary.
  • In video game simulations where you are dealing with object collisions, you might represent a concave object by its convex hull. Here is an interesting post on that.
  • Cooking - in this ingredient polyhedron, the ingredient proportions within the convex hull are proportions that will work for a recipe.

Thanks for reading!

As always, let me know if you have questions or comments by writing them below or reaching out to me on Twitter @LVNGD.

blog comments powered by Disqus

Recent Posts

cloud grid
Deploying a Flask app on AWS Lambda with Zappa
July 10, 2021

In this post I wanted to document the steps I took to deploy a Flask app to AWS Lambda using Zappa. Configuring everything properly in AWS was the most complicated part, so I hope this post can help other AWS noobs who might be struggling!

Read More
D3 Data Join
Data Joins in D3
March 10, 2021

Data joins in D3 can be a tricky thing to wrap your mind around, but once you do, you can take your visualizations to the next level with animations. Data joins are a core concept in D3, so it's a good idea to get acquainted with them.

Read More
rectCollideBlogCover.png
Rectangular Collision Detection in D3 Force Layouts
Feb. 4, 2021

D3 has a built-in force to detect circle collisions in force layouts, but what if you're working with rectangles? In this post we will go over how to detect and resolve collisions, and then adapt D3's built-in forceCollide to work on rectangles.

Read More
Get the latest posts as soon as they come out!