在无序点集上放置网格的算法

Question

给定表示为三维笛卡尔矢量的大量无序点集（数以万计到百万级别），如何制作一个包围所有点的正方形规则网格（用户定义间距）的良好算法？一些约束条件：

1. 网格需要是正方形和规则的；
2. 我需要能够调整网格间距（一个正方形边长的长度），最好只用一个变量；
3. 我想要尽可能小的网格，即网格中的每个“块”都应该至少包含一个无序点，并且每个无序点都应该被包含在一个“块”中；
4. 算法的返回值应该是网格点的坐标列表。

为了在二维平面上说明，考虑以下点集：

对于某个间距为X的网格，算法的一种可能的返回值是这些红色点的坐标（虚线仅用于说明）：

对于间距为X/2的网格，算法的一种可能的返回值是这些红色点的坐标（虚线仅用于说明）：

我正在处理的无序点是大型蛋白质分子的原子坐标，类似于.pdb文件中所包含的内容。

Python 是首选的解决方案语言，但伪代码也可以。

编辑：我认为我对我需要的描述可能有点模糊，因此我添加了一些约束条件和图像，以澄清问题。

- tel

---

那么“最小包含正方形”是一个有效的解决方案吗？因为它至少包含一个点并且所有点都在其中？我不这么认为。 - Dan D.

他所说的“最小尺寸网格”实际上是指“具有最小单元粒度的网格”。否则，对于整个点集的边界框将是最优解。 - Chris

@DanD。这是一个起点 :) 但我想要的部分是能够任意指定网格点的间距。换句话说，我想用一个变量来改变网格的密度。 - tel

@DeepYellow 这些网格点将用作蛋白质静电场计算中的探针。这是否属于“快速接近性检查”？我更像是一位物理学家而不是计算机科学家。 - tel

@tel: 是的，它确实可以。你可能想从这里开始了解：http://en.wikipedia.org/wiki/Spatial_index#Spatial_Index。Blender已经提到了k-d树，但你应该了解其他选择。 - Codie CodeMonkey

Answer 1

我建议你制作一个 k-d树。它速度较快，简单易实现：

和维基百科代码：

class Node: pass

def kdtree(point_list, depth=0):
    if not point_list:
        return

    # Select axis based on depth so that axis cycles through all valid values
    k = len(point_list[0]) # assumes all points have the same dimension
    axis = depth % k

    # Sort point list and choose median as pivot element
    point_list.sort(key=lambda point: point[axis])
    median = len(point_list) // 2 # choose median

    # Create node and construct subtrees
    node = Node()
    node.location = point_list[median]
    node.left_child = kdtree(point_list[:median], depth + 1)
    node.right_child = kdtree(point_list[median + 1:], depth + 1)
    return node

不过，您需要稍微修改它以符合您的限制。

- Blender

Answer 2

怎么样沃罗诺伊图？可以使用福尔图算法以O（n log n）的时间生成。

我不知道它是否解决了您的问题，但沃罗诺伊图非常“自然”。它们在自然界中非常普遍。

例如（来自维基百科）：

- Michał Šrajer

---

对于任何想要将另一个实现推向世界的人：Voronoi / Delaunay的分治算法比Fortune的算法更加合理。两者都存在严重问题，会出现近似退化的浮点输入。 - the guy formerly known as d

@Jimmy：谢谢 - 已经修复了。那是从其他问题中引用的。 - Michał Šrajer

Answer 3

因为您正在寻求用户指定间距的常规正方形网格，所以似乎一个相当简单的方法应该可以解决。

首先通过数据来计算每个维度中的最小和最大坐标。计算需要覆盖最大值和最小值之间距离的用户指定间距的步数。

再次通过数据来分配每个点到网格中的一个单元格，使用具有每个坐标的最小值和指定间距的点的网格（例如X_cell = Math.floor((x_i - x_min) / spacing)）。使用字典或数组记录每个单元格中的点数。

现在打印出至少有一个点的单元格的坐标。

您有一些自由度，我尚未尝试优化：除非最小和最大坐标之间的距离是网格间距的精确倍数，否则会有一些余地，使您可以滑动网格并仍然包含所有点：目前网格从最低点的位置开始，但它可能在最高点之前结束，因此您可以在每个维度上将其向下移动一点。随着这样做，一些点将从一个单元格移动到另一个单元格，占用的单元格数量将发生变化。

如果您只考虑一次移动一个维度，那么您可以相当有效地计算出会发生什么。计算每个点与其单元格在该维度上的最大坐标之间的距离，然后对这些值进行排序。随着网格向下移动，距离其最大坐标最小的点将首先交换单元格，并且您可以按排序顺序逐个迭代这些点。如果您在执行此操作时更新单元格中的点数计数，则可以计算出哪个移位使占用的单元格数量最小。

当然，您需要考虑三个维度。您可以一次处理一个维度，直到减少单元格数量。这是一个局部最小值，但可能不是全局最小值。寻找其他局部最小值的一种方法是从随机选择的起始点重新开始。

- mcdowella

Answer 4

找到一个最小面积的正方形，将所有点都包含在内。重复将每个正方形分成4个子正方形（从1到4到16到64到...）。在其中一个正方形变为空之前停止。很容易证明，所得到的网格最多只比最优解粗糙四倍（关键见解：空正方形保证至少包含任何至少细两倍的网格中的一个正方形）。

可能通过引入随机平移可以减少该常数。

- the guy formerly known as d

Answer 5

我有2D网格聚类的经验，并在C#代码中实现了一个示例。 http://kunuk.wordpress.com/2011/09/15/clustering-grid-cluster/

这可以处理步骤1、2和4。 您需要修改代码并将其更新为3D空间。希望这能给您一些想法。

代码运行时间为O(m*n)，其中m是网格数，n是点数。

- Kunukn

Answer 6

如果您希望网格单元是正方形和规则的，那么您很可能需要一个八叉树。如果您可以放松正方形和规则的约束条件，您可以制作一个k-d树。

- James Thompson