如何确定哪些点在多边形内部，哪些不在（有大量的点）？涉及IT技术。

Question

我有一个大的数据点集（100,000+），存储在一个二维numpy数组中（第一列：x坐标，第二列：y坐标）。我还有几个一维数组，存储每个数据点的附加信息。现在我想从这些1D数组的子集中创建图形，其中仅包括给定多边形中的点。

我想出了以下解决方案，既不优雅也不快:

#XY is the 2D array.
#A is one of the 1D arrays.
#poly is a matplotlib.patches.Polygon

mask = np.array([bool(poly.get_path().contains_point(i)) for i in XY])

matplotlib.pylab.hist(A[mask], 100)
matplotlib.pylab.show()

你能帮我改进一下这段代码吗？我尝试使用np.vectorize代替列表推导式进行操作，但无法使其正常工作。

- AbuBakr

Answer 1

使用matplotlib.nxutils.points_inside_poly进行测试非常有效。

在matplotlib FAQ中，有此40年老算法的示例及进一步的解释。

更新：请注意，points_inside_poly从matplotlib版本1.2.0开始已被弃用。请改用matplotlib.path.Path.contains_points。

- Chewie

---

我知道pnpoly页面。两年前，我拿了Fortran代码并用f2py包装它......我不知道现在它作为mpl模块可用！谢谢。 - oz123

只是想指出，这个函数已经被弃用自 matplotlib 1.2.0 版本以来 - 文档告诉你应该使用 matplotlib.path.Path.contains_points 替代。 - ali_m

@ali_m 感谢你提醒，我已经更新了答案。 - Chewie

Answer 2

很抱歉，我不熟悉你所使用的库，但我认为我有一个合理的算法想法，接下来我将使用纯 Python 实现它，并且相信你可以改进并使用这些库来实现它。此外，我并不保证这是实现此目的的最佳方法，但我想尽快给出我的答复。

---

现在，这个想法来自于在寻找一组点的凸集的算法中使用两个向量的叉积，例如 Graham's Scan。假设我们有两个点p1和p2，它们定义了从原点（0,0）到（x1，y1）和（x2，y2）的点向量 p1 和 p2 。 p1 x p2 的叉积给出第三个向量 p3 ，它垂直于 p1 和 p2 ，其大小由向量所限定的平行四边形的面积确定。

一个非常有用的结果是矩阵的行列式：

/ x1, x2 \
\ y1, y2 /

...这是x1*y2 - x2*y1，它给出了向量p3的大小，并且符号表明p3是“从”平面“出来”还是“进入”平面。关键点在于，如果这个大小是正的，那么p2就在p1的“左边”，如果是负的，那么p2就在p1的“右边”。

希望这个ASCII艺术的例子能有所帮助：

    . p2(4, 5)
   /
  /
 /
/_ _ _ _ _. p1(5, 0)

x1*y2 - x2*y1 = 5*4 - 0*5 = 20，因此p2在p1的“左侧”

最后让我们看看为什么这对我们有用！如果我们有一个多边形顶点的列表和图中其他点的集合，那么对于多边形的每条边，我们都可以获得该边的向量。我们还可以获取连接起始顶点与图中所有其他点的向量，并通过测试它们是否位于边的左侧或右侧来消除每条边的一些点。在此过程结束时未被移除的所有点都在多边形内部。无论如何，下面是一些代码，以使这更加清晰！

---

获取按逆时针方向绘制多边形时访问这些顶点的顺序列表，例如某个五边形可能是：

poly = [(1, 1), (4, 2), (5, 5), (3, 8), (0, 4)]

获取包含图中所有其他点的集合，我们将逐步从该集合中删除无效点，直到最后留下的点恰好是在多边形内部的点。

points = set(['(3, 0), (10, -2), (3,3), ...])

实际上，代码本身非常紧凑，但我花了很长时间来解释它的工作原理。 to_right 接受表示向量的两个元组，并返回True如果v2位于v1的右侧。然后，循环只需遍历多边形的所有边缘，并从工作集中删除点，如果它们在任何边缘的右侧。

def to_right(v1, v2):
    return (v1[0]*v2[1] - v1[1]*v2[0]) < 0

for i in range(len(poly)):
    v1 = poly[i-1]
    v2 = poly[i]
    for p in points:
        if(to_right(v2-v1, p-v1)):
            points.remove(p)

编辑：为了澄清，如果点在多边形的左侧而不是右侧，则会被删除与指定多边形顶点的顺序有关。如果它们按顺时针顺序排列，则希望消除左侧点。目前我没有特别好的解决方案。

---

无论如何，希望我关于这些内容的理解是正确的，并且对某人有所帮助，即使不是原帖作者。该算法的渐近复杂度为O(mn)，其中n是图中点的数量，m是多边形的顶点数，在最坏情况下所有点都在多边形内部，我们必须检查每个边上的每个点，而没有任何一个点会被删除。

- Max Spencer

---

很遗憾，这并不是我想要了解的内容，因为该算法已经在matplotlib库中实现。但是阅读这篇文章非常有趣，我现在对它的工作原理有了很好的理解。感谢您的贡献！ - AbuBakr

没问题，我今年的算法课程中有一些关于这样的几何算法的材料，尝试回答你的问题实际上帮助我更好地理解了它们。 - Max Spencer