12

QHullの SciPy ラッパーを使用して、一連の点の凸包のボリュームを取得しようとしています。

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__']
4

2 に答える 2

15

ConvexHull渡すパラメータに関係なく、求めている結果を直接取得する明確な方法はないようです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]))

EDITコメントによると、凸包の体積を取得するより高速な方法は次のとおりです。

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 番目の方法の方が約 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
于 2014-07-14T10:34:07.717 に答える
11

この質問は 2 回目の誕生日を迎えましたが、scipy ラッパーは Qhull によって計算された体積 (および面積) を自動的に報告することを指摘したいと思います。

于 2016-08-24T08:16:36.350 に答える