一种基于散乱点云的铆钉孔孔位提取方法*

田清廉，熊天辰，黄 翔，李泷杲，郝 龙

（1.南京航空航天大学机电学院，南京 210016；2.江西航天海虹测控技术有限责任公司，南昌 330024）

我国航空事业不断发展，对零件上铆钉孔的制孔精度要求也越来越高，而飞机结构与系统中使用的大量零件包含铆钉孔特征[1]。通常对孔特征检测时，使用三坐标测量仪测出各个孔的坐标尺寸，但这种测量方式效率低，且不适合检测壁板等大型零件；使用专用孔检测工具 （如孔检测样板）对孔的尺寸和形位进行检测，但制造这种高精度检具耗时长且成本高[2]。随着非接触式数字化测量技术的提高，非接触式数字化测量设备也开始更多地应用在飞机零件铆钉孔检测上，这种检测方式具有精度高、效率高、柔性好等特点[3–4]。

目前，在非接触式圆孔测量方面已有了很多的研究与成果，主要有基于视觉图像的测量方法和基于三维点云的测量方法。基于视觉图像的测量方法主要是通过图像处理中的Canny 算子等边界识别算法提取照相测量图像中的边界特征进行圆孔参数拟合[5–6]，但该方法存在空间圆透视投影畸变的问题[7]，影响检测精度。基于三维点云的测量方法是使用三维激光扫描等测量设备获取零件点云，再提取点云中的孔边界特征进行圆孔参数拟合。针对点云中边界特征提取，Jenke 等[8]提出了一种将点云网格化再通过网格的拓扑结构识别边界特征的方法，这种方法可以有效提取点云孔洞边界特征，但需要网格化整个点云，数据量大时算法效率低。因此出现了直接对三维点云提取孔边界特征的方法，该方法通过建立点云空间拓扑关系，并设定某一特定判断机制用于识别边界特征点。Milroy[9]和Yang[10]等通过求解曲率极值提取点云边界，这种方法可以有效地提取比较光滑平坦的数据模型的边界点，不太适用于曲率变化大的数据模型。张杰等[11]提出了一种基于T–Scan测量的薄壁钣金件孔特征重构方法，通过扫描线点云结构中的孔截断线特征识别平面点和孔径点，拟合求解平面和空间圆，但该方法仅限于扫描线点云类型且提取精度受圆孔直径影响。刘增艺等[12]提出了一种曲面孔位视觉测量技术，使用手持式双目视觉线结构光测量曲面得到曲面点云，然后通过交互式方法提取散乱点云圆孔边缘点，进行圆孔拟合，由于该方法在圆孔识别过程中需要人工预先分割出圆孔所在区域，提取效率低且不能满足多孔自动提取需求。

基于以上分析，为解决目前面向铆钉孔检测领域的测量点云中的孔位与孔径参数提取方法限制多和人工交互较多的问题，本文提出了一种基于散乱点云的铆钉孔孔位与孔径提取方法，该方法可以自动识别散乱点云中的全部铆钉孔特征，并精确提取孔位和孔径参数。采用本文方法对试验件散乱点云进行数据处理，验证该方法的实用性与精度。

1 散乱点云预处理

各种数字化测量设备获得的点云根据其点的分布特征分为4 种类型：散乱点云、扫描线点云、网格化点云和多边形点云[13]。散乱点云相对于其他点云类型，其测量点无明显的几何分布特征，呈散乱无序状态，测量点是随机生成的。由于散乱点云无序和随机的分布特性，使用非接触式测量设备测量铆钉孔得到的点云分布情况如图1所示，大部分测量点分布在孔周围平面上，而孔边界测量点可能分布在孔内径或在孔周围平面上，有极少的点正好在孔边界上。

图1 铆钉孔边界测量点云分布Fig.1 Distribution of measurement point cloud at rivet hole boundary

根据铆钉孔边界测量点云分布情况，本文通过点的k邻域点几何分布来判断是否为边界点。若某一点为边界点，则其k邻域点的分布将偏向于一侧 （图2（a））；若某一点不是边界点，则k邻域点比较均匀地分布在该点周围 （图2（b））。因此，对于铆钉孔的散乱点云，需要进行边界点提取预处理。

图2 基于k 邻域分布的边界点识别Fig.2 Boundary point recognition based on k neighborhood distribution

由于散乱点云在数据存储结构中每个点是无规律排布的，在这样的数据结构中查找每个点的k近邻必须循环每个点，效率低下，因此必须建立点云数据的空间拓扑结构。常见的k近邻搜索算法有八叉树法、空间栅格法和kd树法[14–16]。本文采用kd树法建立散乱点云的kd树结构，从而查询每个点Pi的k邻域点集Nj（j=0，1，2，…，k–1）。

以点Pi和其k邻域点Nj（j=0，1，2，…，k–1）构成的点集X来最小二乘拟合得到Pi点的局部切平面。在切平面上建立平面坐标系，如图3所示，以Pi到切平面的投影点Pi′作为坐标原点，Pi′点和N0点的投影点N0′构成的向量Pi′N0′作为坐标系x轴，平面法向与向量Pi′N0′的叉乘×Pi′N0′作为坐标系y轴。将点集X中每个点的三维坐标都转换到该平面坐标系上，得到点集X的二维坐标点集合X′（u，V）。

图3 建立局部平面坐标系Fig.3 Establishing local plane coordinate system

以X′（u，V）中的Pi″点作为向量起点，Nj′点作为向量终点，得到平面向量集计算中每个向量到局部坐标系的x轴的夹角αj和与y轴的夹角βj。若βj＞π/2，则αj=αj+π。然后将αj以升序进行排列得到角度序列ηj，计算ηj相邻角度之间的夹角θj=ηj–ηj–1，其中，j∈[1，k–1]，如图4所示。

图4 计算相邻向量夹角Fig.4 Calculating angle between adjacent vectors

当相邻角度序列θj中的最大夹角θmax超过最大夹角阈值εθ时，则点Pi为边界点，否则点Pi为非边界点，如图5所示。阈值εθ的大小要视点云分布情况而定，一般将阈值εθ设置为π /2 可得到良好的识别结果。

图5 边界点识别Fig.5 Boundary point recognition

2 铆钉孔边界提取

点云预处理后得到的边界点中不仅包含有多个铆钉孔边界特征，还包含有点云外边界特征以及点云漏洞边界特征，如图6所示。边界点云在数据结构上仍为散乱点云，需要对边界点云建立点云拓扑结构，并对边界点云进行分割，得到属于不同边界特征的边界点云块，最后根据点云块的椭圆度和圆度提取其中的铆钉孔边界。

图6 点云中的边界特征Fig.6 Boundary features in point clouds

2.1 边界点欧式聚类分割

边界点云中的不同边界特征具有同边界特征内相邻点距离连续的特点，即同边界特征内的相邻点距离间距较小，不同边界特征的最近点距离较大。基于此分布特点，本文采用欧式距离聚类分割算法[17]对边界点云中的边界特征进行分割。其边界点云聚类分割步骤如下。

步骤1。对边界点云X建立kd树点云结构，方便后续点邻域搜索。

步骤2。创建空的聚类集C和点集Q。

步骤3。对任意点pi∈X，执行以下操作：（1）把pi加入Q中 （2）对每个点pj∈Q，通过边界点kd树的k近邻搜索算法找到pj的k邻域，放入点集Pjk；对每个pjk∈Pjk，若pjk不在Q中，且pjk与pj的欧式距离rj

步骤4。当边界点云X中每个点都执行过步骤3 后，就得到了聚类集C。

步骤5。删除聚类集C中的点数量少于nmin（nmin为最小聚类点数量阈值）的点云块，得到最终聚类集C。

2.2 椭圆拟合结果提取铆钉孔

边界点云聚类分割后得到属于不同边界特征的点云块，对点云块分别拟合平面并把边界点投影到各自拟合平面上。由于铆钉孔边界点的拟合平面与理论铆钉孔所在圆柱相交成一椭圆形状，边界点在拟合平面上的投影点分布在该椭圆周围，如图7所示。因此本文对所有的边界点聚类分割点云块拟合平面椭圆，通过椭圆拟合结果判断是否为铆钉孔边界点。

图7 铆钉孔边界点椭圆拟合Fig.7 Ellipse fitting of rivet hole boundary points

首先对分割点云块ci∈C最小二乘拟合平面，构建局部平面二维坐标系：以平面上的任一向量⇀ 作为坐标系x轴，平面的法向⇀与⇀ 的叉乘作为坐标系y轴，点云块ci中的任一点p0在平面上的投影点p0′作为坐标系原点。将点云块ci中所有点的三维坐标转换到构建的局部平面二维坐标系下，得到二维坐标点集合ci′{pj′（xj′，yj′），j=1，…，mj}。

然后对平面二维坐标点集合ci′最小二乘拟合平面椭圆。设椭圆方程为Ax2+Bxy+Cy2+Dx+Ey+F=0，并满足椭圆参数约束4AC–B2> 0，拟合优化目标函数为

令W=[A，B，C，D，E，F]T

则优化目标为

式中，H=为椭圆参数约束。

由于||WTX||2时，W存在缩放因子使得所有W′=aW也满足优化目标，因此可令WTHW=1，于是优化目标就变为

构造拉格朗日函数，即

对式（4）进行求解。令S=XXT，则SW=λHW，通过求解广义逆矩阵得到6 个可能的解W。由于S为正定矩阵，所以λ> 0，因此可以用WTHW=1 和λ> 0 这两个条件来筛选得到最终合格的解W=[A，B，C，D，E，F]T。

根据椭圆参数W=[A，B，C，D，E，F]T计算椭圆圆心Oe（x0，y0）和长短轴半径ae、be：

计算椭圆拟合误差Error 与长短轴之比Ratio：

其中，F（p′j）=（Ax′j2+Bxj′yj′+Cy′j2+Dxj′+Eyj′+F），n为 点云块ci∈C中点的数量。

若Error >εe（εe为拟合误差阈值），则为该分割点云块椭圆度低，为非椭圆，即非铆钉孔；若Error <εe，则比较长短轴之比。

若1–εr< Ratio < 1+εr（εr为长短轴比阈值），则该分割点云块圆度高，即为铆钉孔边界，否则为非铆钉孔边界。

若对铆钉孔半径大小有要求，可设置最大半径阈值rmax和最小半径阈值rmin，将长短轴的均值与rmax和rmin比较来提取满足半径要求的铆钉孔。

点云中可能有其他不是铆钉孔的圆孔，需要对其进行剔除。计算非铆钉孔理论圆心与该算法拟合椭圆圆心的空间距离dc，若dc小于一定阈值εd，则剔除该提取孔。

3 铆钉孔参数提取

铆钉孔孔位孔径参数提取采用先投影孔边界点到孔定位面然后最小二乘拟合端面圆的方法。由于孔的边界点中部分点在孔内壁上，使用孔边界点直径求解孔定位面会降低参数提取精度，所以需要寻找孔周围一定范围内的点来拟合孔定位面。然后根据边界点拟合椭圆和孔定位面求得初始铆钉孔参数，根据边界点到初始铆钉孔参数的距离偏差剔除边界点中的噪声点，最后再将孔边界点投影到孔定位面上拟合圆孔。

3.1 孔定位面拟合

在铆钉孔边界提取过程中求得了所有边界点云块的椭圆拟合方程、椭圆圆心Oe和长短轴ae、be，以拟合椭圆的圆心Oe作为搜索中心，长短轴均值的ω倍作为搜索半径Rs，通过建立的散乱点云kd树结构进行半径邻域搜索找出铆钉孔周围Rs范围内的点集Pr，如图8所示。

图8 铆钉孔邻域点云搜索Fig.8 Point cloud search of rivet hole neighborhood

将点集Pr精确拟合平面。设孔定位面方程为ax+by+cz+d=0，构建拟合目标函数：

用特征值法[18]求解式（8），得到平面参数a、b、c、d。然后计算点集中每个点pi∈Pr到平面的距离di，当di>εd（εd为平面拟合误差阈值）时，则认为该点是噪点并剔除，反之保留，最后使用保留的点重新拟合孔定位面。

3.2 端面圆孔拟合

在孔定位面上建立平面二维坐标系，将铆钉孔边界点的三维坐标转换成平面二维坐标点pj（xj，yj）。然后根据拟合椭圆和定位面计算初始铆钉孔参数，将椭圆圆心投影到孔定位面上，投影点即为铆钉孔圆心，椭圆短轴长度即为铆钉孔半径。

计算平面二维坐标点pj（xj，yj）到初始铆钉孔的距离δj，当δj>εc（εc为圆孔拟合误差阈值）时，则认为该点为噪点并剔除，反之保留，最后使用保留的点重新拟合端面圆。

设端面圆孔方程为（x–A）2+（y–B）2=R2，令a=–2A，b=–2B，c=A2+B2–R2，则方程变为x2+y2+ax+by+c=0，构建拟合目标函数：

式（9）中分别对a、b、c求偏导，令偏导等于0，可得：

解式（10）可得端面圆参数a、b、c，则端面圆的二维圆心坐标为（x0，y0）=（a/–2，b/–2），半径

把二维圆心（x0，y0）坐标逆变换到三维空间中，得到三维空间圆心坐标（x′0，y′0，z′0），三维空间圆半径R等于平面圆半径r，孔法向（i，j，k）等于拟合平面的法向，铆钉孔参数提取完成。

4 试验验证

4.1 提取效果验证

应用本文方法对散乱点云中的铆钉孔特征进行提取，验证其提取效果，如图9所示。图9（a）为结构光扫描仪测得的某壁板零件的散乱点云数据，所含点数为90688；用本文边界点提取算法提取点云中的边界点，k=50、εθ=π/2 时如图9（b）所示，提取边界点数为3392；用本文边界点云欧氏距离聚类分割算法分割边界点，dth=2mm、nmin=10 时分割结果如图9（c）所示，共分割出14 个边界特征；依据椭圆拟合结果提取铆钉孔边界特征，εe=1、εr=0.05 时提取结果如图9（d）所示，14 个边界特征中的9 个铆钉孔特征全部提取成功。可以看出，本文方法可以有效提取散乱点云中的铆钉孔边界特征。

图9 铆钉孔特征提取效果Fig.9 Results of rivet hole feature extraction

4.2 提取精度验证

为验证本文提出算法的铆钉孔参数提取精度，采用天远Robot-Scan 扫描仪对铆钉孔试验件进行测量和数据处理，Robot-Scan 相机像素为1000 万像素，平均单张点间距为0.16mm，单张测量精度为15μm。

如图10所示，使用Robot-Scan 扫描仪对试验件上的10 个理论直径d=8mm 的铆钉孔进行测量得到试验件散乱点云数据。应用本文方法对数据进行处理，得到铆钉孔参数。然后将三坐标测量机测量试验件得到的参数作为真值，与本文方法提取结果进行比较得出本文方法铆钉孔参数提取偏差。三坐标测量机测量结果和本文参数提取结果如表1所示，铆钉提取偏差如图11所示。可以看出，本文提出方法对试验件上铆钉孔的孔位提取偏差均值为0.017mm，最大值为0.029mm；铆钉孔孔径提取偏差均值为0.0797mm，最大值为0.095mm。

图10 使用Robot–Scan 对铆钉孔试验件测量Fig.10 Measurement of rivet hole test pieces with Robot–Scan

表1 铆钉孔参数提取试验结果Table 1 Results of rivet hole parameter extraction experiment

图11 铆钉孔提取偏差Fig.11 Extraction deviation of rivet holes

从提取结果分析可以看出该方法孔位提取精度较高，可达0.029mm；铆钉孔孔径提取误差相比较孔位提取误差较大，这是由于测量设备测量得到的散乱点云具有一定的点间距，孔边界测量点不一定全部落在孔内径上，从而使得该方法对铆钉孔的直径提取精度相对于孔位提取精度较低。

5 结论

（1）本文提出了一种基于散乱点云的铆钉孔孔位参数提取方法。用该方法对试验件测量散乱点云数据进行了处理。试验表明，该方法能够精确提取散乱点云中的全部铆钉孔特征，且孔位提取精度高，提取精度可达0.029mm。

（2）使用本文方法需要注意，由于测量点不能够保证全部落在铆钉孔径上，导致了铆钉孔孔径提取误差相对较大。因此，为保证有较高的孔径提取精度，需要测量设备精度和点测量密度高，且尽量使测量点能够较多地落在孔内径上。