使用SciPy中的QHull计算凸包体积

Question

使用SciPy中的QHull计算凸包体积

12

我正在尝试使用SciPy对QHull的包装来获取一组点的凸壳体积。

根据QHull的文档，我应该传递"FA"选项以获取总表面积和体积。

这是我的结果。我做错了什么吗？

> pts
     [(494.0, 95.0, 0.0), (494.0, 95.0, 1.0) ... (494.0, 100.0, 4.0), (494.0, 100.0, 5.0)]


> hull = spatial.ConvexHull(pts, qhull_options="FA")

> dir(hull)

     ['__class__', '__del__', '__delattr__', '__dict__', '__doc__', '__format__', '__getattribute__', '__hash__', '__init__', '__module__', '__new__', '__reduce__', '__reduce_ex__', '__repr__', '__setattr__', '__sizeof__', '__str__', '__subclasshook__', '__weakref__', '_qhull', '_update', 'add_points', 'close', 'coplanar', 'equations', 'max_bound', 'min_bound', 'ndim', 'neighbors', 'npoints', 'nsimplex', 'points', 'simplices']

 > dir(hull._qhull)
     ['__class__', '__delattr__', '__doc__', '__format__', '__getattribute__', '__hash__', '__init__', '__new__', '__reduce__', '__reduce_ex__', '__repr__', '__setattr__', '__sizeof__', '__str__', '__subclasshook__']

- Diana

1

尝试用一个真正的问题来更新你的问题（而不是“这是我得到的”）。花了我一些时间才想明白，尽管你提供了正确的选项，但总面积和体积都找不到。 - Hugues Fontenelle

我猜测SciPy没有包装那个特定的选项标志。 - Hugues Fontenelle

实现它的困难方式：http://wiki.scipy.org/Cookbook/Finding_Convex_Hull - Hugues Fontenelle

有一个有帮助的事情是完整的“pts”。这样我们就可以自己尝试。 - DrV

它在Scipy Qhull包装器中没有实现。如果有需要，可以很容易地添加。 - pv.

2个回答

13

尽管这个问题已经两岁了，但我想指出的是现在scipy包装器会自动报告由Qhull计算得出的体积（和面积）。

- Urukann

1

请注意，您需要至少使用scipy版本0.17.0才能使其正常工作。 - awakenting

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Jaime · Accepted Answer

无论您传递什么参数，似乎没有直接获得所需结果的明显方法。如果使用Delaunay（还提供大多数凸包相关信息），自己计算应该不难。

def tetrahedron_volume(a, b, c, d):
    return np.abs(np.einsum('ij,ij->i', a-d, np.cross(b-d, c-d))) / 6

from scipy.spatial import Delaunay

pts = np.random.rand(10, 3)
dt = Delaunay(pts)
tets = dt.points[dt.simplices]
vol = np.sum(tetrahedron_volume(tets[:, 0], tets[:, 1], 
                                tets[:, 2], tets[:, 3]))

编辑根据评论，以下是获取凸包体积更快的方法：

def convex_hull_volume(pts):
    ch = ConvexHull(pts)
    dt = Delaunay(pts[ch.vertices])
    tets = dt.points[dt.simplices]
    return np.sum(tetrahedron_volume(tets[:, 0], tets[:, 1],
                                     tets[:, 2], tets[:, 3]))

def convex_hull_volume_bis(pts):
    ch = ConvexHull(pts)

    simplices = np.column_stack((np.repeat(ch.vertices[0], ch.nsimplex),
                                 ch.simplices))
    tets = ch.points[simplices]
    return np.sum(tetrahedron_volume(tets[:, 0], tets[:, 1],
                                     tets[:, 2], tets[:, 3]))

使用一些虚构的数据，第二种方法似乎比第一种方法快大约2倍，并且数值精度非常好（15位小数！！！），尽管可能存在一些更为异常的情况：

pts = np.random.rand(1000, 3)

In [26]: convex_hull_volume(pts)
Out[26]: 0.93522518081853867

In [27]: convex_hull_volume_bis(pts)
Out[27]: 0.93522518081853845

In [28]: %timeit convex_hull_volume(pts)
1000 loops, best of 3: 2.08 ms per loop

In [29]: %timeit convex_hull_volume_bis(pts)
1000 loops, best of 3: 1.08 ms per loop