[frames] | no frames]

Source Code for Module gmisclib.convex_hull2d

1  #!/usr/bin/env python
2
3  """convexhull.py
4
5  Calculate the convex hull of a set of n 2D-points in O(n log n) time.
6  Taken from Berg et al., Computational Geometry, Springer-Verlag, 1997.
7  Prints output as EPS file.
8
9  When run from the command line it generates a random set of points
10  inside a square of given length and finds the convex hull for those,
11  printing the result as an EPS file.
12
13  Usage: convexhull.py <numPoints> <squareLength> <outFile>
14
15  Dinu C. Gherman
16
17  Small Bug: Only works with a list of UNIQUE points, Evan Jones, 2005/05/18
18  If the list of points passed to this function is not unique, it will raise an assertion.
19  To fix this, remove these lines from the beginning of the convexHull function:
20
21  Taken from http://aspn.activestate.com/ASPN/Cookbook/Python/Recipe/66527
22  and modified to work with complex numbers.
23
24  """
25
26
27  import sys, random
28
29 -class DuplicatePoints(ValueError):
30 - def __init__(self, s):
31 ValueError.__init__(self, s)
32 33 34 ###################################################################### 35 # Helpers 36 ###################################################################### 37
38 -def _myDet(p, q, r):
39 """Calc. determinant of a special matrix with three 2D points. 40 41 The sign, "-" or "+", determines the side, right or left, 42 respectivly, on which the point r lies, when measured against 43 a directed vector from p to q. 44 """ 45 46 # We use Sarrus' Rule to calculate the determinant. 47 # (could also use the Numeric package...) 48 # sum1 = q[0]*r[1] + p[0]*q[1] + r[0]*p[1] 49 # sum1 = q.real*r.imag + p.real*q.imag + r.real*p.imag 50 # sum2 = q[0]*p[1] + r[0]*q[1] + p[0]*r[1] 51 # sum2 = q.real*p.imag + r.real*q.imag + p.real*r.imag 52 sum = q.real*(r.imag-p.imag) + p.real*(q.imag-r.imag) + r.real*(p.imag-q.imag) 53 54 # return sum1 - sum2 55 return sum
56 57
58 -def _isRightTurn((p, q, r)):
59 "Do the vectors pq:qr form a right turn, or not?" 60 if p==q or p==r or q==r: 61 if p==q: 62 dup = p 63 elif p==r: 64 dup = p 65 elif q==r: 66 dup = q 67 raise DuplicatePoints, 'Two points are identical at %s' % str(dup) 68 return _myDet(p, q, r) < 0
69 70
71 -def _isPointInPolygon(r, P):
72 "Is point r inside a given polygon P?" 73 74 # We assume the polygon is a list of points, listed clockwise! 75 for i in xrange(len(P[:-1])): 76 p, q = P[i], P[i+1] 77 if not _isRightTurn((p, q, r)): 78 return False # Out! 79 80 return True # It's within!
81 82