2

我试图找出一个点是否在 3D 多边形中。我使用了另一个我在网上找到的脚本来处理许多使用光线投射的 2D 问题。我想知道如何将其更改为适用于 3D 多边形。我不会看有很多凹面或孔或任何东西的非常奇怪的多边形。这是python中的2D实现:

def point_inside_polygon(x,y,poly):

    n = len(poly)
    inside =False

    p1x,p1y = poly[0]
    for i in range(n+1):
        p2x,p2y = poly[i % n]
        if y > min(p1y,p2y):
            if y <= max(p1y,p2y):
                if x <= max(p1x,p2x):
                    if p1y != p2y:
                        xinters = (y-p1y)*(p2x-p1x)/(p2y-p1y)+p1x
                    if p1x == p2x or x <= xinters:
                        inside = not inside
        p1x,p1y = p2x,p2y

    return inside

任何帮助将不胜感激!谢谢你。

4

3 回答 3

5

我检查了 QHull 版本(从上面)和线性规划解决方案(例如看到这个问题)。到目前为止,使用 QHull 似乎是最好的选择,尽管我可能会错过一些scipy.spatialLP 的优化。

import numpy
import numpy.random
from numpy import zeros, ones, arange, asarray, concatenate
from scipy.optimize import linprog

from scipy.spatial import ConvexHull

def pnt_in_cvex_hull_1(hull, pnt):
    '''
    Checks if `pnt` is inside the convex hull.
    `hull` -- a QHull ConvexHull object
    `pnt` -- point array of shape (3,)
    '''
    new_hull = ConvexHull(concatenate((hull.points, [pnt])))
    if numpy.array_equal(new_hull.vertices, hull.vertices): 
        return True
    return False


def pnt_in_cvex_hull_2(hull_points, pnt):
    '''
    Given a set of points that defines a convex hull, uses simplex LP to determine
    whether point lies within hull.
    `hull_points` -- (N, 3) array of points defining the hull
    `pnt` -- point array of shape (3,)
    '''
    N = hull_points.shape[0]
    c = ones(N)
    A_eq = concatenate((hull_points, ones((N,1))), 1).T   # rows are x, y, z, 1
    b_eq = concatenate((pnt, (1,)))
    result = linprog(c, A_eq=A_eq, b_eq=b_eq)
    if result.success and c.dot(result.x) == 1.:
        return True
    return False


points = numpy.random.rand(8, 3)
hull = ConvexHull(points, incremental=True)
hull_points = hull.points[hull.vertices, :]
new_points = 1. * numpy.random.rand(1000, 3)

在哪里

%%time
in_hull_1 = asarray([pnt_in_cvex_hull_1(hull, pnt) for pnt in new_points], dtype=bool)

产生:

CPU times: user 268 ms, sys: 4 ms, total: 272 ms
Wall time: 268 ms

%%time
in_hull_2 = asarray([pnt_in_cvex_hull_2(hull_points, pnt) for pnt in new_points], dtype=bool)

生产

CPU times: user 3.83 s, sys: 16 ms, total: 3.85 s
Wall time: 3.85 s
于 2016-10-29T23:59:13.267 回答
4

这里提出了一个类似的问题,但重点是效率

@Brian@fatalaccidentsscipy.spatial.ConvexHull在此建议的方法有效,但是如果您需要检查多个点,则速度会非常慢。

好吧,最有效的解决方案也来自scipy.spatial,但使用了Delaunay镶嵌:

from scipy.spatial import Delaunay

Delaunay(poly).find_simplex(point) >= 0  # True if point lies within poly

这是有效的,因为如果该点不在任何单纯形中(即在三角剖分之外),-1则返回。.find_simplex(point)(注意:它适用于 N 维,而不仅仅是 2/3D。)


性能比较

首先一点

import numpy
from scipy.spatial import ConvexHull, Delaunay

def in_poly_hull_single(poly, point):
    hull = ConvexHull(poly)
    new_hull = ConvexHull(np.concatenate((poly, [point])))
    return np.array_equal(new_hull.vertices, hull.vertices)

poly = np.random.rand(65, 3)
point = np.random.rand(3)

%timeit in_poly_hull_single(poly, point)
%timeit Delaunay(poly).find_simplex(point) >= 0

结果:

2.63 ms ± 280 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
1.49 ms ± 153 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

所以这个Delaunay方法更快。但这取决于多边形的大小!我发现对于一个由超过 65 个点组成的多边形,Delaunay方法变得越来越慢,而ConvexHull方法的速度几乎保持不变。

对于多个点

def in_poly_hull_multi(poly, points):
    hull = ConvexHull(poly)
    res = []
    for p in points:
        new_hull = ConvexHull(np.concatenate((poly, [p])))
        res.append(np.array_equal(new_hull.vertices, hull.vertices))
    return res

points = np.random.rand(10000, 3)

%timeit in_poly_hull_multi(poly, points)
%timeit Delaunay(poly).find_simplex(points) >= 0

结果:

155 ms ± 9.42 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)
1.81 ms ± 106 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

所以Delaunay给了极大的速度增加;更不用说我们必须等待多长时间才能获得 10'000 分或更多。在这种情况下,多边形大小不再有太大的影响。


综上所述,Delaunay不仅速度快得多,而且代码也非常简洁。

于 2020-03-13T14:28:42.820 回答
1

感谢所有评论过的人。对于任何为此寻找答案的人,我发现了一个适用于某些情况(但不是复杂情况)的方法。

我正在做的是像 shongololo 建议的那样使用 scipy.spatial.ConvexHull ,但略有不同。我正在制作点云的 3D 凸包,然后将要检查的点添加到“新”点云中并制作新的 3D 凸包。如果它们相同,那么我假设它必须在凸包内。如果有人有更强大的方法来做到这一点,我仍然会很感激,因为我认为这有点骇人听闻。代码如下所示:

from scipy.spatial import ConvexHull

def pnt_in_pointcloud(points, new_pt):
    hull = ConvexHull(points)
    new_pts = points + new_pt
    new_hull = ConvexHull(new_pts)
    if hull == new_hull: 
        return True
    else:
        return False

希望这可以帮助将来寻找答案的人!谢谢!

于 2015-03-31T18:59:24.033 回答