Python中如何计算指定测量区域内包含的Voronoi单元格

我有逐帧标识行人头部点位的实验数据，需要判断哪些Voronoi单元格落在指定测量矩形区域内，区域参数如下：

x_range = (-0.4, 0.4)
y_range = (0.5, 1.3)

我参考scipy.spatial.Voronoi官方示例实现了点位的Voronoi图生成，还叠加绘制了测量区域（蓝线）和墙体（黑线），第0帧可视化结果如下：

我的部分实现代码（改编自官方示例）如下：

entries_for_frame = get_entries_at_frame(entries, frame)
points = points_from_entries(entries_for_frame)

vor = scipy.spatial.Voronoi(points)
scipy.spatial.voronoi_plot_2d(vor)

plt.show()

我原以为要判断单元格的边是否和测量矩形相交、或者完全落在区域内即可实现需求。根据scipy.spatial.Voronoi文档，vertices属性返回图中橙色的Voronoi顶点，ridge_vertices属性返回顶点构成的边，但该属性返回的结果如下，我无法理解返回数值的含义，也没有找到对应教程解决我的问题：

[[0, 19], [0, 2], [1, 17], [1, 3], [2, 3], [17, 19], [-1, 22], [-1, 15], [15, 16], [16, 21], [21, 22], [-1, 0], [2, 23], [22, 23], [28, 32], [28, 29], [29, 30], [30, 31], [31, 32], [12, 13], [12, 28], [13, 25], [25, 29], [-1, 24], [-1, 31], [24, 30], [-1, 26], [26, 27], [27, 32], [-1, 33], [19, 20], [20, 34], [33, 34], [35, 36], [-1, 35], [36, 37], [-1, 37], [-1, 4], [4, 5], [5, 35], [6, 37], [6, 7], [7, 36], [38, 39], [38, 40], [39, 41], [40, 41], [-1, 40], [-1, 8], [8, 38], [-1, 9], [9, 10], [10, 41], [10, 43], [39, 42], [42, 43], [52, 53], [52, 57], [53, 54], [54, 55], [55, 56], [56, 57], [13, 52], [25, 57], [48, 49], [48, 54], [49, 55], [9, 50], [24, 56], [49, 50], [17, 59], [18, 61], [18, 20], [59, 61], [11, 46], [11, 60], [18, 47], [46, 47], [60, 61], [58, 63], [58, 60], [59, 62], [62, 63], [26, 64], [27, 65], [64, 65], [21, 67], [23, 68], [67, 68], [42, 45], [43, 69], [44, 45], [44, 72], [69, 72], [50, 70], [69, 70], [48, 71], [70, 71], [4, 76], [5, 75], [75, 76], [33, 77], [76, 77], [34, 78], [77, 78], [47, 79], [78, 79], [80, 82], [80, 81], [81, 83], [82, 84], [83, 84], [14, 53], [14, 80], [71, 82], [72, 84], [14, 51], [51, 87], [81, 85], [85, 87], [88, 90], [88, 89], [89, 93], [90, 91], [91, 92], [92, 93], [44, 88], [83, 89], [85, 86], [86, 93], [11, 91], [58, 92], [94, 95], [94, 97], [95, 96], [96, 98], [97, 99], [98, 99], [12, 94], [51, 95], [65, 97], [101, 104], [101, 102], [102, 103], [103, 105], [104, 105], [15, 101], [16, 104], [64, 102], [99, 103], [66, 67], [66, 105], [1, 106], [3, 107], [106, 107], [68, 108], [107, 108], [8, 73], [45, 109], [73, 110], [109, 110], [111, 115], [111, 113], [112, 113], [112, 114], [114, 115], [46, 74], [74, 111], [79, 113], [75, 112], [7, 114], [116, 117], [116, 118], [117, 120], [118, 119], [119, 121], [120, 121], [96, 118], [98, 100], [100, 116], [87, 119], [86, 121], [63, 120], [122, 127], [122, 123], [123, 124], [124, 125], [125, 126], [126, 127], [100, 127], [117, 122], [62, 123], [106, 124], [108, 125], [66, 126], [128, 129], [128, 130], [129, 132], [130, 131], [131, 133], [132, 134], [133, 134], [90, 128], [109, 129], [74, 130], [110, 132], [115, 131], [6, 133], [73, 134]]

问题：如何计算至少有部分落在指定测量区域内的Voronoi单元格？ #

解答 #

首先明确ridge_vertices的数值含义：每个子列表对应一条Voronoi边的两个顶点索引，索引值对应vor.vertices的下标，值为-1说明这条边是开放边，延伸到无穷远，没有对应的实体顶点。

不需要自己手动计算边相交，直接用几何库判断即可，代码实现如下：

import numpy as np
from shapely.geometry import Polygon, Point

# 构造测量矩形
x_min, x_max = -0.4, 0.4
y_min, y_max = 0.5, 1.3
measure_rect = Polygon([(x_min, y_min), (x_max, y_min), (x_max, y_max), (x_min, y_max)])

target_cells = []
# 遍历每个原始点对应的Voronoi单元格
for point_idx, region_idx in enumerate(vor.point_region):
    region = vor.regions[region_idx]
    # 跳过无效空区域
    if not region:
        continue
    # 单元格对应的原始点位
    point = Point(vor.points[point_idx])
    # 如果原始点已经在矩形内，直接判定符合要求
    if measure_rect.contains(point):
        target_cells.append(point_idx)
        continue
    # 处理带无穷边的开放单元格
    if -1 in region:
        # 只要开放单元格的任意一条边和矩形相交就符合要求，这里用缓冲判断简化逻辑
        if measure_rect.intersects(point.buffer(0.1)):
            target_cells.append(point_idx)
        continue
    # 封闭单元格直接构造多边形判断相交
    cell_vertices = vor.vertices[region]
    cell_poly = Polygon(cell_vertices)
    if measure_rect.intersects(cell_poly):
        target_cells.append(point_idx)

最终输出的target_cells就是所有至少有部分落在测量区域内的Voronoi单元格对应的原始行人点位索引，可直接用于后续计算。

内容的提问来源于stack exchange，提问作者baziorek