Search code examples
pythonmatplotlibscipydelaunaypolyhedra

Delaunay triangularization of Polyhedron (Python)


I'm trying to get the Delaunay Triangulation of a polyhedron in python so that I can calculate the centroid. I see that there is a Delaunay function in scipy.spatial and that it works in n-dimensions. The trouble is that the documentation shows 2D use and gives me no indication of what to do with higher dimensions. Being able to decompose this object into an array would probably solve this issue for me, but I don't know how to do that.

The problem I'm running into is that I do not know how to verify that this is working correctly as it is outputting an object. I can find nothing on Google about how to graph a polyhedron or how to use this object that scipy is spitting back.

If I do

import numpy as np
from scipy.spatial import Delaunay

points = np.array([[0,0,0],[1,0,0],[1,1,0],[1,0,1],[1,1,1],[0,1,0],[0,1,1],[0,0,1]])
Delaunay(points)

I really would just like to be able to get back the coordinates of these tetrahedrons so that I can calculate the centroids of the polyhedrons. It would also be really nice if I were able to graph the tesselated polyhedron too. I saw in MATLAB that I can do this with a fuction called trimesn, and I found one from matplotlib but it seems to be really different and its documentation is not great.

from matplotlib.collections import TriMesh TriMesh.__doc__ 

u'\n      Class for the efficient drawing of a triangular mesh using\n   
Gouraud shading.\n\n    A triangular mesh is a
:class:`~matplotlib.tri.Triangulation`\n    object.\n    '

Solution

  • What tess = Delaunay(pts) returns is an object of the Delanauy class. You can check the tetrahedrons as tess.simplices. It has different attributes and methods. In 2D, for example, it can plot you triangulation, convex hull and Voronoi tesselation.

    Regarding the visualization of the final collection of tetrahedrons I didn't find a straightforward way of doing it, but I managed to get a working script. Check the code below.

    from __future__ import division
    import numpy as np
    from scipy.spatial import Delaunay
    import matplotlib.pyplot as plt
    from mpl_toolkits.mplot3d import Axes3D
    from mpl_toolkits.mplot3d.art3d import Poly3DCollection, Line3DCollection
    from itertools import combinations
    
    
    def plot_tetra(tetra, pts, color="green", alpha=0.1, lc="k", lw=1):
        combs = combinations(tetra, 3)
        for comb in combs:
            X = pts[comb, 0]
            Y = pts[comb, 1]
            Z = pts[comb, 2]
            verts = [zip(X, Y, Z)]
            triangle = Poly3DCollection(verts, facecolors=color, alpha=0.1)
            lines = Line3DCollection(verts, colors=lc, linewidths=lw)
            ax.add_collection3d(triangle)
            ax.add_collection3d(lines)
    
    pts = np.array([
                [0,0,0],
                [1,0,0],
                [1,1,0],
                [1,0,1],
                [1,1,1],
                [0,1,0],
                [0,1,1],
                [0,0,1]])
    tess = Delaunay(pts)
    
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    for k, tetra in enumerate(tess.simplices):
        color = plt.cm.Accent(k/(tess.nsimplex - 1))
        plot_tetra(tetra, pts, color=color, alpha=0.1, lw=0.5, lc="k")
    ax.scatter(pts[:, 0], pts[:, 1], pts[:, 2], c='k')
    plt.savefig("Delaunay.png", dpi=600)
    plt.show()
    

    The resulting image is

    enter image description here