从叶片截面散点提取Medial Axis的Voronoi算法问题咨询

我需要从叶片高度方向的2D散点中提取Medial Axis，输入散点为叶片表面全区域的无序点集。我尝试了Voronoi算法，得到的结果中蓝色点为输入点，红色点为Voronoi顶点，但效果远未达到预期。

我的两个问题 #

1. 为何当前Voronoi算法结果效果很差？我期望的效果是连贯的叶片骨架形态（注：该效果对应的输入仅为叶片轮廓边界点，与当前输入不同）。
2. 如何提取内部Voronoi顶点而非全部？内部Voronoi顶点指的是位于叶片轮廓内部区域的红色点。

我参考的算法来自论文《A New Voronoi-Based Surface Reconstruction Algorithm》，希望从无序2D散点得到理想的Medial Axis。以下是我的代码：

import numpy as np
import open3d as o3d
from scipy.spatial import Voronoi, Delaunay, voronoi_plot_2d
import matplotlib.pyplot as plt

# 获取输入点
point3d = o3d.io.read_point_cloud(r'E:\Evaluation\data\data925\model\cmm_nom\cmm925_nom457.20.ply')
pointCoordinate_3d = np.asarray(point3d.points)
pointCoordinate_2d = pointCoordinate_3d[:, :2]

# 生成Voronoi图
vor = Voronoi(pointCoordinate_2d)
plt.plot(pointCoordinate_2d[:, 0], pointCoordinate_2d[:, 1], 'b.', markersize=8, label='点云')
finite_region = np.where(np.all(np.isfinite(vor.vertices), axis=1))
plt.plot(vor.vertices[finite_region, 0], vor.vertices[finite_region, 1], 'r.', markersize=8, label='Voronoi顶点')
plt.xlim([pointCoordinate_2d[:, 0].min()-15, pointCoordinate_2d[:, 0].max()+15])
plt.ylim([pointCoordinate_2d[:, 1].min()-15, pointCoordinate_2d[:, 1].max()+15])
plt.show()

问题解答 #

问题1：当前Voronoi结果效果差的原因 #

核心原因在于输入点的分布差异：

• 你期望的流畅中轴效果，对应的输入是叶片轮廓边界点——此时Voronoi顶点会集中在轮廓的“中心骨架”区域，自然形成连贯的形态。
• 而当前输入是叶片表面全区域散点（包含大量内部点），Voronoi顶点是至少三个输入点的外接圆圆心，内部点会分割出大量细碎的局部区域，导致顶点分散在整个叶片范围内，完全偏离了中轴应有的定位，最终呈现杂乱无章的结果。

问题2：提取内部Voronoi顶点的方法 #

要筛选出叶片内部的Voronoi顶点，可按以下步骤实现：

1. 拟合叶片的精确轮廓：
   对输入的2D散点使用Alpha Shape算法（或Delaunay三角剖分提取边界），生成叶片的多边形轮廓，这是判断点是否在内部的基础。
2. 点-in-多边形筛选：
   用点-in-多边形算法逐个检查有限Voronoi顶点是否位于拟合出的轮廓内部。示例代码如下：

   from shapely.geometry import Point, Polygon
   from alphashape import alphashape  # 需要先安装alphashape库：pip install alphashape
   
   # 用Alpha Shape拟合叶片轮廓
   alpha = 0.1  # 可根据点云密度调整alpha值，越小轮廓越贴合
   contour_polygon = alphashape(pointCoordinate_2d, alpha)
   # 提取轮廓点坐标
   contour_points = np.array(list(contour_polygon.exterior.coords))
   
   # 筛选内部Voronoi顶点
   inner_voronoi_vertices = []
   finite_vertices = vor.vertices[finite_region]
   for v in finite_vertices:
       point = Point(v[0], v[1])
       if contour_polygon.contains(point):
           inner_voronoi_vertices.append(v)
   inner_voronoi_vertices = np.array(inner_voronoi_vertices)
   
   # 绘制筛选后的结果
   plt.plot(pointCoordinate_2d[:,0], pointCoordinate_2d[:,1], 'b.', markersize=8)
   plt.plot(inner_voronoi_vertices[:,0], inner_voronoi_vertices[:,1], 'r.', markersize=8)
   plt.show()
   

3. 优化中轴连贯性：
   筛选后的顶点若仍有细碎噪声，可通过DBSCAN聚类去除离群点，再用B样条插值拟合出连贯的中轴曲线。

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