/
pyhull_demo.py
126 lines (96 loc) · 3.3 KB
/
pyhull_demo.py
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
113
114
115
116
117
118
119
120
121
122
123
124
125
126
#!/usr/bin/env python
"""
Generates the Delaunay triangulation, Voronoi tessellation and convex hull for
a random set of 2D points and plot them. For testing purposes.
Red dots - random points.
Green lines - Delaunay triangulation.
Blue lines - convex hull.
Black lines - Voronoi tesselation.
Dash black lines - Voronoi tesselation with points at infinity.
"""
from __future__ import division
__author__ = "Shyue Ping Ong"
__version__ = "0.1"
__maintainer__ = "Shyue Ping Ong"
__email__ = "shyuep@gmail.com"
__status__ = "Beta"
__date__ = "11/19/12"
import itertools
import argparse
import numpy as np
from pyhull.delaunay import DelaunayTri
from pyhull.convex_hull import ConvexHull
from pyhull.voronoi import VoronoiTess
import matplotlib.pyplot as plt
parser = argparse.ArgumentParser(description="Demo script for pyhull.",
epilog="Author: Shyue Ping Ong")
parser.add_argument("-d", "--dimension", dest="dim", type=int,
choices=[2, 3],
default=2,
help="Dimension of points.")
parser.add_argument("-n", "--npoints", dest="npts", type=int,
default=30,
help="Number of points.")
args = parser.parse_args()
npts = args.npts
dim = args.dim
points = np.random.randn(npts, dim)
legkeys = []
legitems = []
if dim == 3:
from mpl_toolkits.mplot3d import Axes3D
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
for pt in points:
p, = plt.plot(pt[0], pt[1], 'ro') if dim == 2 \
else plt.plot([pt[0]], [pt[1]], [pt[2]], 'ro')
legkeys.append(p)
legitems.append("Points")
d = DelaunayTri(points)
for s in d.simplices:
for data in itertools.combinations(s.coords, dim):
data = np.array(data)
p, = plt.plot(data[:,0], data[:,1], 'g-') if dim == 2 else\
ax.plot(data[:,0], data[:,1], data[:,2], 'g-')
legkeys.append(p)
legitems.append("Delaunay tri")
d = ConvexHull(points)
for s in d.simplices:
for data in itertools.combinations(s.coords, dim):
data = np.array(data)
p, = plt.plot(data[:,0], data[:,1], 'b-') if dim == 2 else\
ax.plot(data[:,0], data[:,1], data[:,2], 'b-')
legkeys.append(p)
legitems.append("Convex hull")
if dim == 2:
d = VoronoiTess(points)
points = d.points
avg = np.average(points, 0)
vertices = d.vertices
for nn, vind in d.ridges.items():
(i1, i2) = sorted(vind)
if i1 == 0:
c1 = np.array(vertices[i2])
midpt = 0.5 * (np.array(points[nn[0]]) + np.array(points[nn[1]]))
if np.dot(avg - midpt, c1 - midpt)> 0:
c2 = c1 + 10 * (midpt-c1)
else:
c2 = c1 - 10 * (midpt-c1)
p1, = plt.plot([c1[0],c2[0]], [c1[1], c2[1]], 'k--')
else:
c1 = vertices[i1]
c2 = vertices[i2]
p, = plt.plot([c1[0],c2[0]], [c1[1], c2[1]], 'k-')
legkeys.append(p)
legitems.append("Voronoi tess.")
legkeys.append(p1)
legitems.append("Voronoi tess. - inf vert.")
#Set some sensible limits
maxc = np.amax(points, 0)
minc = np.amin(points, 0)
plt.xlim((minc[0]-0.1, maxc[0]+0.1))
plt.ylim((minc[1]-0.1, maxc[1]+0.1))
font = {'family' : 'Times New Roman',
'size' : 18}
plt.rc('font', **font)
plt.show()