在python中查找点是否在3D多边形中

Question

我试图找出一个点是否在 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

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

Answer 1

这里提出了类似的问题，但重点是效率。

\n\n

@Brian和@fatalaccidentsscipy.spatial.ConvexHull此处建议的方法有效，但速度非常慢，但如果您需要检查多个点，则

\n\n

嗯，最有效的解决方案，也来自scipy.spatial，但是利用Delaunay曲面细分：

\n\n

from scipy.spatial import Delaunay\n\nDelaunay(poly).find_simplex(point) >= 0  # True if point lies within poly\n

\n\n

这有效，因为-1返回的是.find_simplex(point)如果该点不在任何单纯形中（即在三角剖分之外），则

\n\n\n\n

性能对比

\n\n

首先是为了一点：

\n\n

import numpy\nfrom scipy.spatial import ConvexHull, Delaunay\n\ndef in_poly_hull_single(poly, point):\n    hull = ConvexHull(poly)\n    new_hull = ConvexHull(np.concatenate((poly, [point])))\n    return np.array_equal(new_hull.vertices, hull.vertices)\n\npoly = np.random.rand(65, 3)\npoint = np.random.rand(3)\n\n%timeit in_poly_hull_single(poly, point)\n%timeit Delaunay(poly).find_simplex(point) >= 0\n

\n\n

结果：

\n\n

2.63 ms \xc2\xb1 280 \xc2\xb5s per loop (mean \xc2\xb1 std. dev. of 7 runs, 100 loops each)\n1.49 ms \xc2\xb1 153 \xc2\xb5s per loop (mean \xc2\xb1 std. dev. of 7 runs, 100 loops each)\n

\n\n

因此该Delaunay方法速度更快。但这取决于多边形的大小！我发现对于由超过 65 个点组成的多边形，该Delaunay方法变得越来越慢，而ConvexHull接近速度几乎保持不变。

\n\n

对于多个点：

\n\n

def in_poly_hull_multi(poly, points):\n    hull = ConvexHull(poly)\n    res = []\n    for p in points:\n        new_hull = ConvexHull(np.concatenate((poly, [p])))\n        res.append(np.array_equal(new_hull.vertices, hull.vertices))\n    return res\n\npoints = np.random.rand(10000, 3)\n\n%timeit in_poly_hull_multi(poly, points)\n%timeit Delaunay(poly).find_simplex(points) >= 0\n

\n\n

结果：

\n\n

155 ms \xc2\xb1 9.42 ms per loop (mean \xc2\xb1 std. dev. of 7 runs, 10 loops each)\n1.81 ms \xc2\xb1 106 \xc2\xb5s per loop (mean \xc2\xb1 std. dev. of 7 runs, 1000 loops each)\n

\n\n

所以Delaunay给出了极大的速度提升；更不用说要等多久才能达到10'000点甚至更多。在这种情况下，多边形大小不再有太大的影响。

\n\n\n\n

总之，Delaunay不仅速度快很多，而且代码也非常简洁。

\n

Answer 2

我检查了 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