SciPy空间Delaunay/ConvexHull混乱

Question

我正在尝试生成随机凸多面体。我生成一组随机 3D 坐标，然后找到它们的凸包（到目前为止还不错）。

然后我想我会使用 Delaunay 三角剖分给我一个凸包的三角剖分。这是我的基本理解开始显示的地方！

这是代码

import numpy as np
from scipy.spatial import ConvexHull
from scipy.spatial import Delaunay
import matplotlib as mpl
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# Generate random points & convex hull
points = np.random.rand(20,3)
hull = ConvexHull(points)

fig = plt.figure()
ax = fig.gca(projection = '3d')

# Plot hull's vertices
for vert in hull.vertices:    
    ax.scatter(points[vert,0], points[vert,1], zs=points[vert,2])#, 'ro')

# Calculate Delaunay triangulation & plot
tri = Delaunay(points[hull.vertices])
for simplex in tri.simplices:
    vert1 = [points[simplex[0],0], points[simplex[0],1], points[simplex[0],2]]
    vert2 = [points[simplex[1],0], points[simplex[1],1], points[simplex[1],2]]
    vert3 = [points[simplex[2],0], points[simplex[2],1], points[simplex[2],2]]
    vert4 = [points[simplex[3],0], points[simplex[3],1], points[simplex[3],2]]
    ax.plot([vert1[0], vert2[0]], [vert1[1], vert2[1]], zs = [vert1[2], vert2[2]])
    ax.plot([vert2[0], vert3[0]], [vert2[1], vert3[1]], zs = [vert2[2], vert3[2]])
    ax.plot([vert3[0], vert4[0]], [vert3[1], vert4[1]], zs = [vert3[2], vert4[2]])
    ax.plot([vert4[0], vert1[0]], [vert4[1], vert1[1]], zs = [vert4[2], vert1[2]])

plt.show()

有几件事让我担心，情节有时会漏掉船体上的一些点，这似乎是 Delaunay 四面体化，我想我不应该对此感到惊讶，但不是我想要的。

我只想对船体表面进行三角剖分，所以我想是一个包含表面刻面的单纯形？这可能吗？

谢谢

B

编辑：在 pv 的启示性 post 之后，我将代码修改如下；

import numpy as np
import pylab as pl
import scipy as sp
from scipy.spatial import ConvexHull
from scipy.spatial.distance import euclidean
import matplotlib.pyplot as plt
import mpl_toolkits.mplot3d as a3

aspect = 0
while aspect == 0:

    # Generate random points & convex hull
    points = np.random.rand(20,3)
    hull = ConvexHull(points)

    # Check aspect ratios of surface facets
    aspectRatio = []
    for simplex in hull.simplices:
        a = euclidean(points[simplex[0],:], points[simplex[1],:])
        b = euclidean(points[simplex[1],:], points[simplex[2],:])
        c = euclidean(points[simplex[2],:], points[simplex[0],:])
        circRad = (a*b*c)/(np.sqrt((a+b+c)*(b+c-a)*(c+a-b)*(a+b-c)))
        inRad = 0.5*np.sqrt(((b+c-a)*(c+a-b)*(a+b-c))/(a+b+c))
        aspectRatio.append(inRad/circRad)

    # Threshold for minium allowable aspect raio of surface facets
    if np.amin(aspectRatio) > 0.3:
        aspect = 1

ax = a3.Axes3D(pl.figure())
facetCol = sp.rand(3) #[0.0, 1.0, 0.0]

# Plot hull's vertices
#for vert in hull.vertices:    
#    ax.scatter(points[vert,0], points[vert,1], zs=points[vert,2])

# Plot surface traingulation
for simplex in hull.simplices:
    vtx = [points[simplex[0],:], points[simplex[1],:], points[simplex[2],:]]
    tri = a3.art3d.Poly3DCollection([vtx], linewidths = 2, alpha = 0.8)
    tri.set_color(facetCol)
    tri.set_edgecolor('k')
    ax.add_collection3d(tri)

plt.axis('off')
plt.show()

现在一切都如我所愿。我添加了纵横比阈值以确保更好的三角剖分。

B

Answer 1

一些事情：

您将 points[hull.vertices] 作为 Delaunay 的参数，因此 tri.simplices 中的整数是 points[hull.vertices] 中的索引，而不是 points 中的索引，因此您最终绘制错误点数
四面体有 6 个脊，但您只绘制了 4 个
如果您只需要凸包曲面的三角剖分，可以使用 hull.simplices

即

for simplex in hull.simplices:
    xs, ys, zs = points[simplex].T
    xs = np.r_[xs, xs[0]] # close polygons
    ys = np.r_[ys, ys[0]]
    zs = np.r_[zs, zs[0]]
    ax.plot(xs, ys, zs)

或者只是：

ax.plot_trisurf(points[:,0], points[:,1], points[:,2],
                triangles=hull.simplices)

SciPy空间Delaunay/ConvexHull混乱

SciPy spatial Delaunay/ConvexHull confusion

python

3d

scipy

delaunay

convex-hull