Skip to content

Instantly share code, notes, and snippets.

@arthur-e
Forked from tixxit/hull.py
Last active April 7, 2024 17:19
Show Gist options
  • Save arthur-e/5cf52962341310f438e96c1f3c3398b8 to your computer and use it in GitHub Desktop.
Save arthur-e/5cf52962341310f438e96c1f3c3398b8 to your computer and use it in GitHub Desktop.
Graham's scan convex hull algorithm, updated for Python 3.x
def convex_hull_graham(points):
'''
Returns points on convex hull in CCW order according to Graham's scan algorithm.
By Tom Switzer <thomas.switzer@gmail.com>.
'''
TURN_LEFT, TURN_RIGHT, TURN_NONE = (1, -1, 0)
def cmp(a, b):
return (a > b) - (a < b)
def turn(p, q, r):
return cmp((q[0] - p[0])*(r[1] - p[1]) - (r[0] - p[0])*(q[1] - p[1]), 0)
def _keep_left(hull, r):
while len(hull) > 1 and turn(hull[-2], hull[-1], r) != TURN_LEFT:
hull.pop()
if not len(hull) or hull[-1] != r:
hull.append(r)
return hull
points = sorted(points)
l = reduce(_keep_left, points, [])
u = reduce(_keep_left, reversed(points), [])
return l.extend(u[i] for i in range(1, len(u) - 1)) or l
@SalahAdDin
Copy link

Hey man, thank you very much!

@james-langbein
Copy link

Thank you so much for posting this code, I've been banging my head against a brick wall all day trying to get scipy.spatial.ConvexHull to work, and this code worked perfectly first time!

@jxi24
Copy link

jxi24 commented Apr 2, 2019

I think I found a bug in your code. For certain sets of points, example given the points:
[[181, 864],
[182, 859],
[182, 860],
[182, 861],
[182, 862],
[182, 863],
[182, 864]]

The above algorithm finds the hull to just be [181, 864],[182, 864]]. This is obviously incorrect. If you change line 15 to be:
while len(hull) > 1 and turn(hull[-2], hull[-1], r) == TURN_RIGHT

The hull is given by: [[181, 864],[182, 859], [182, 864]].

This allows the hull to contain points that have no turns which occurs for topologies in which most of the points occur on a line with a few not on the line.

@10avinash
Copy link

I don't know but it was showing an error stating that subtraction was deprecated.
Hence, I changed part of code to this,

def cmp(a, b):
    return ((a>b).astype(np.float32) - (a<b).astype(np.float32)).astype(np.bool)

@samykbs
Copy link

samykbs commented Jun 23, 2020

great code, do you recommend some lectures to get more familiar with graham algorithm thanks

@Vanman989
Copy link

Awesome. Is it possible to adjust it for tolerance?

@charlie5610
Copy link

Thank You!

@SinghShreya05
Copy link

I too faced this error :

return (a > b) - (a < b)
TypeError: numpy boolean subtract, the `-` operator, is not supported, use the bitwise_xor, the ^ operator, or the logical_xor function instead.

To fix this, I typecasted it to float :

def cmp(a, b):
        return float(a > b) - float(a < b)

@aclark-aquaveo
Copy link

Nice! Thanks.

@ferrine
Copy link

ferrine commented Oct 26, 2022

Numba version for the graham_hull

import numba
import numpy as np

@numba.njit
def cmp(a, b):
    return (a > b) - (a < b)

@numba.njit
def turn(p, q, r):
    return cmp((q[..., 0] - p[..., 0])*(r[..., 1] - p[..., 1]) - (r[..., 0] - p[..., 0])*(q[..., 1] - p[..., 1]), 0)


@numba.njit
def reduce_keepleft(points):
    left = []
    for r in points:
        while len(left) > 1 and turn(left[-2], left[-1], r) != 1:
            left.pop()
        if not len(left) or not np.allclose(left[-1], r, 1e-5, 1e-8, False):
            left.append(r)
    return left

@numba.njit
def sort(points):
    for i in range(points.shape[-1]-1, -1, -1):
        idx = np.argsort(points[:, i], kind="mergesort")
        points = points[idx]
    return points

@numba.njit
def stack(alist):
    out = np.empty((len(alist), *alist[0].shape))
    for i, r in enumerate(alist):
        out[i] = r
    return out

@numba.njit
def convex_hull_graham(points):
    points = sort(points)
    l = reduce_keepleft(points)
    u = reduce_keepleft(points[::-1])
    hull = l + u[1:-1]
    return stack(hull)

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment