Convex Hull
 
05_source
20161206
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
#!/usr/bin/python
import math
import matplotlib.pyplot as plt
import copy 
import random      # only needed for random points
import numpy as np # only needed for random points from normal distrib

# convenience functions to get the x,y value from a tuple, 
# and avoid using the wrong tuple index
# (sacrifice a bit of performance for legibility) 
def x(tp): 
    return tp[0]

def y(tp): 
    return tp[1]

# calculate the angle for a line between two points and the horizontal
def calc_angle(pt1, pt2): 
    dy=float( y(pt2) - y(pt1) ) 
    dx=float( x(pt2) - x(pt1) ) 
    angle=0.0
    if dx!=0:
        angle=math.atan(dy/dx)   
    else:
        if dy==0: 
            angle=0.0
        elif dy<0: 
            angle=math.atan(float('-inf'))
        else: 
            angle=math.atan(float('inf'))
    # if the calucated angle is negative, make it positive
    while (angle<0): 
        angle=angle+2*math.pi;
    return angle

# calculate the distance between two points 
def dist(pt1,pt2): 
    dx=float( x(pt2) - x(pt1) )
    dy=float( y(pt2) - y(pt1) )
    return math.sqrt( dx*dx + dy*dy)  

# does moving from point pts0 -> pt1 -> pt2 describe a counterclockwise motion ? 
def clockwise(pt1, pt2, pt3): 
    area = (x(pt2)-x(pt1))*(y(pt3)-y(pt1))-(y(pt2)-y(pt1))*(x(pt3)-x(pt1))
    rv=0               # collinear
    if area<0: rv=-1   # clockwise
    elif area>0: rv=+1 # counter clockwise
    return rv 


## MAIN ----------------------------------------------------------------------

# --- DATA ------------
points= [ (9,2), (5,6), (3,10), (11,13), (7,6), (15,18), (19,16), (13,4), (23,8), (17,12)  ]

# variation: random points
#points=[]
#for i in range(1,1000): points.append( (random.randint(1,100), random.randint(1,100)) )

# variation: sedgewick's points (approx) 
#points=[]
# points= [(373, 234), (312, 235), (313, 128), (165, 164), (185, 91), (190, 34), 
#          (452, 249), (404, 109), (237, 297), (262, 130), (124, 127), (334, 29), (459, 2)]

# variation : random points, but normally distributed
#points=list(zip(
#    [ int(r) for r in np.random.normal(100, 20, 1000)],
#    [ int(r) for r in np.random.normal(100, 20, 1000)]
#    ) )
# ---------------------

# identify the point with the lowest y coordinate and put it first in the list of points
# special case: in case of tie, decide by lowest x coordinate
points=sorted( points, key=lambda a: a[0]) # secondary sort on x  
points=sorted( points, key=lambda a: a[1]) # primary sort on y (python sort is stable) 

# enrich the tuple list with the angle between the line through (point, first point) 
# and the horizontal 
points = [ (x(pt),y(pt),calc_angle(points[0],pt)) for pt in points ] 

# order the points by angle
points=sorted(points, key=lambda a: a[2])

# after sorting we don't need the angles anymore, let's drop them
points = [ (x(pt),y(pt)) for pt in points ] 

# the core of the graham scan: 
#   iterate over the points in ascending angle order
#   check if pt1-> pt2 -> pt3 makes a counterclockwise turn
#   if NOT, then discard pt2
hull=copy.deepcopy(points)  # deep copy, because we are going to delete points 
i=0
while True: 
    if clockwise(hull[i],hull[i+1],hull[i+2])<=0: # clockwise or collinear
        del hull[i+1]
        i-=1    # take a step back 
    else: 
        i+=1    # step forward
        if (i+2)>=len(hull): # are we done?
            break

print "Points on the hull" ,hull

# to close the polygon, (for drawing, and calculating distance) append the first point 
hull.append(hull[0])

print "Perimeter distance: ", sum( [ dist( hull[i],hull[i+1] ) for i in range(len(hull)-1) ] )

# plotting the points and the hull
plt.scatter(*zip(*points),s=20)      # plot all the points
plt.plot(*zip(*hull),linestyle='-')  # plot the hull
plt.savefig("convex_hull.png")
 
Notes by Data Munging Ninja. Generated on nini:sync/20151223_datamungingninja/convexhull at 2016-12-26 08:22