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")
|