Delaunay Triangulation

%matplotlib inline
from scipy.spatial import Delaunay
from mpl_toolkits import mplot3d

N=20
X,Y = np.meshgrid(np.linspace(-1,1,N), np.linspace(-1,1,N))
X += (np.random.rand(N,N)*2-1)/20
Y += (np.random.rand(N,N)*2-1)/20
points = np.vstack([X.reshape(N*N), Y.reshape(N*N)]).T
print(points[0:5])

[[-0.9977976  -0.99539644]
 [-0.89534579 -0.96342752]
 [-0.82520533 -1.04061474]
 [-0.6529773  -0.96018126]
 [-0.58339513 -1.04623718]]

fig, ax = plt.subplots(figsize=(8,8))
ax.plot(points[:,0], points[:,1], 'o')
ax.set_aspect('equal')

msk = np.linalg.norm(points, axis=1) < 1
print(msk[0:50])
fig, ax = plt.subplots(figsize=(8,8))
ax.plot(points[msk,0], points[msk,1], 'o')
ax.set_aspect('equal')

[False False False False False False False False False False  True  True
 False False False False False False False False False False False False
 False  True  True  True  True  True  True  True  True  True  True False
 False False False False False False False False  True  True  True  True
  True  True]

z = -(3*points[:,0]-0.25)**2 - (3*points[:,1]+0.25)**2
fig, ax = plt.subplots(figsize=(8,8))
p=ax.scatter(points[:,0], points[:,1], c=z)
fig.colorbar(p)
ax.set_aspect('equal')

fig = plt.figure(figsize=(8,8))
ax = plt.axes(projection='3d')
#ax.contour3D(points[:,0].reshape(N,N), points[:,1].reshape(N,N), z.reshape(N,N), 50)
ax.plot(points[:,0], points[:,1], z, 'o')

[<mpl_toolkits.mplot3d.art3d.Line3D at 0x7fed1248acf8>]

dln = Delaunay(points)
print(dln.simplices[0:10])

[[120 101 121]
 [101 102 121]
 [101 120 100]
 [ 12   4  15]
 [399 379 339]
 [378 379 399]
 [397 398 399]
 [398 378 399]
 [360 381 380]
 [300 320 280]]

fig, ax = plt.subplots(figsize=(8,8))
ax.plot(points[:,0], points[:,1], 'o')
for i, tri in enumerate(dln.simplices):
    for j in range(3):
        k = (j+1)%3
        ax.plot(points[[tri[j], tri[k]],0], points[[tri[j], tri[k]],1], '-k')
ax.set_aspect('equal')

fig = plt.figure(figsize=(8,8))
ax = plt.axes(projection='3d')
ax.plot_trisurf(points[:,0], points[:,1], z, triangles=dln.simplices, cmap='viridis')

<mpl_toolkits.mplot3d.art3d.Poly3DCollection at 0x7fed62116240>

print(dln.convex_hull)

[[  4  15]
 [397 399]
 [381 380]
 [320 280]
 [  2   0]
 [  0 100]
 [399  39]
 [397 387]
 [320 380]
 [381 382]
 [  4   2]
 [ 17  15]
 [280 100]
 [382 387]
 [ 19  18]
 [ 17  18]
 [ 19  39]]

fig, ax = plt.subplots(figsize=(8,8))
ax.plot(points[:,0], points[:,1], 'o')
for i, seg in enumerate(dln.convex_hull):
    ax.plot(points[[seg[0], seg[1]],0], points[[seg[0], seg[1]],1], '-k')
ax.set_aspect('equal')

dlnmsk = Delaunay(points[msk])
print(dlnmsk.simplices[0:10])
fig, ax = plt.subplots(figsize=(8,8))
ax.plot(points[msk,0], points[msk,1], 'o')
for i, seg in enumerate(dlnmsk.convex_hull):
    ax.plot(points[msk][[seg[0], seg[1]],0], points[msk][[seg[0], seg[1]],1], '-k')
ax.set_aspect('equal')

[[233 247 232]
 [247 261 232]
 [246 230 260]
 [212 245 228]
 [272 283 281]
 [271 272 281]
 [135 116 117]
 [274 276 282]
 [276 277 282]
 [168 150 151]]