I need to do Delaunay Triangulation for a set of 3D points. I wrote an script for it (below), but it seems to that the output has no tetrahedrons in them. Please give me some inputs/ideas. I am using Python3. Thank you very much.
from scipy.spatial import Delaunay
import matplotlib.pyplot as plt
import numpy as np
points= np.array([[1,2,2],[1,3,6],[4,3,4],[5,3,2]])
tri= Delaunay(points)
fig= plt.figure()
ax= fig.gca(projection= '3d')
ax.plot_trisurf(points[:,0],points[:,1],points[:,2],triangles= tri.simplices)
plt.plot(points[:,0],points[:,1],points[:,2],'+')
plt.show()
The tetrahedrons are given in the tri.simplices member, which holds an n x 4 array of indices (n being the number of tetrahedrons). The tetrahedron is given as a set of four indices, which correspond to the indices of the four points of the tetrahedron in the points array.
For example the following code will plot the wireframe of the first tetrahedron:
tr = tri.simplices[0] # indices of first tetrahedron
pts = points[tr, :] # pts is a 4x3 array of the tetrahedron coordinates
# plotting the six edges of the tetrahedron
for ij in [[0, 1], [0, 2], [0, 3], [1, 2], [1, 3], [2, 3]]:
ax.plot3D(pts[ij, 0], pts[ij, 1], pts[ij, 2])
See my previous answers here, here and here for further example code.