Объем выпуклой оболочки с QHull от SciPy

Я пытаюсь получить объем выпуклой оболочки набора точек, используя оболочка 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__']

python scipy convex-hull qhull
person Diana    schedule 14.07.2014    source источник
comment
Попробуйте обновить свой вопрос реальным вопросом (не здесь то, что я получаю). Мне потребовалось некоторое время, чтобы понять, что нигде нельзя найти общую площадь и объем, несмотря на то, что вы указали правильный вариант.   -  person Hugues Fontenelle    schedule 14.07.2014
comment
Мое дикое предположение заключается в том, что SciPy не оборачивает этот конкретный флаг опции.   -  person Hugues Fontenelle    schedule 14.07.2014
comment
Сложный способ — это реализовать: wiki.scipy.org/Cookbook/Finding_Convex_Hull   -  person Hugues Fontenelle    schedule 14.07.2014
comment
Одна вещь, которая могла бы помочь, была бы полной pts. Таким образом, мы могли бы попробовать это сами.   -  person DrV    schedule 14.07.2014
comment
Это не реализовано в оболочках Scipy Qhull. Его можно легко добавить, если есть необходимость.   -  person pv.    schedule 14.07.2014

Ответы (2)


arrow_upward
15
arrow_downward

Кажется, не существует какого-либо очевидного способа прямого получения результатов, которые вам нужны, независимо от того, какие параметры вы передаете. Не должно быть слишком сложно вычислить себя, если вместо 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]))

РЕДАКТИРОВАТЬ Согласно комментариям, следующие более быстрые способы получения выпуклого объема корпуса:

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
person Jaime    schedule 14.07.2014
comment
Это здорово. Я искал информацию об объеме из qhull, чтобы сэкономить время вычислений, но я не думаю, что это можно сделать на Python быстрее, чем ваша версия. Итак, чтобы подтвердить, Делоне правильно ли разбивает выпуклую оболочку? Кроме того, эта einsum хорошенькая штучка... Не знал этого. :) - person Diana; 15.07.2014
comment
Делоне оказался очень медленным, поэтому вместо этого я использовал ConvexHull, получил крайние точки (те, что на корпусе) и запустил Делоне только по этим точкам. На моих данных это на 1-2 порядка быстрее. - person Diana; 15.07.2014
comment
Это действительно умно... Это может не сильно улучшиться, но вы можете попробовать вообще пропустить вызов Delaunay и построить триангуляцию вашей выпуклой оболочки, выбрав точку на оболочке, а затем вычислив объем всех тетраэдров, которые содержат эту точку и точки на каждой из симплициальных граней выпуклой оболочки (т. е. атрибут .simplices ConvexHull объектов). Он будет производить более вытянутые тетраэдры, поэтому он может быть менее численно стабильным. Но это должно быть быстрее. - person Jaime; 15.07.2014
comment
Я тоже так пробовал, но разница в результатах была слишком большой. Может быть, я сделал что-то не так. Численная нестабильность не может объяснить разницу в два порядка в значении объема, учитывая, что это уже выпуклая оболочка, верно? Я, пожалуй, еще раз посмотрю на это. - person Diana; 15.07.2014
comment
Взгляните на последнее редактирование, для случайных данных я получаю 2-кратное ускорение и точность до 15 знаков после запятой... - person Jaime; 16.07.2014
comment
Я только что попробовал это и в своем коде и понял, что я сделал неправильно в прошлый раз (перепутал значения вершин и индексы вершин). Мой ConvexHull не имеет атрибута vertices, но вместо этого я просто использовал первую точку в своем списке точек. Вероятно, это также помогает с числовой устойчивостью, если она оказывается внутренней точкой. Я думаю, что это настолько быстро, насколько это возможно. Спасибо! - person Diana; 16.07.2014
comment
Хотя это очень хорошо и поучительно, это не последний ответ на вопрос - person gota; 22.07.2019

arrow_upward
11
arrow_downward

Хотя этот вопрос отпраздновал свой второй день рождения, я хотел бы отметить, что теперь оболочка scipy автоматически сообщает объем (и площадь), вычисленный Qhull.

person Urukann    schedule 24.08.2016
comment
Обратите внимание, что для этого вам нужна как минимум версия scipy 0.17.0. - person awakenting; 26.10.2017