CN110458174A - 一种无序点云关键特征点精确提取方法 - Google Patents

Abstract

一种无序点云关键特征点精确提取方法，其包含有如下步骤：1）以每个点云点为中心，计算一环k邻域并构建三角化网格；2）计算每个点云点的离散曲率、离散法向量和邻域点距离用来构建窗口函数，计算梯度矩阵；3）对梯度矩阵计算角点响应值分类特征点。本发明使用梯度变化识别点云数据中角点、折痕点以及平面点，保留了局部细节特征，识别结果准确、可靠；可广泛应用于航空、汽车等零部件的三维激光扫描检测，具有很强的适用性。

Description

一种无序点云关键特征点精确提取方法

技术领域

本发明涉及计算机图形学和三维点云特征提取领域，具体涉及一种无序点云关键特征点精确提取方法。

背景技术

飞机零部件具有大尺寸，弱刚性及结构复杂的特点，不仅需要测量单一零件的空间尺寸，而且复杂部件的装配特征也需要进行测量。这些结构件的实测点云量大且无序，准确提取关键特征点是精简数据、提高数据分析处理效率的必要手段。特征提取是测量信息处理中最基本、最关键的部分，逆向工程、图像处理的精度严重依赖于特征提取的精度。

国内外在点云关键特征提取领域已经有了一定的研究成果，主要采用三种方式：一是基于法向量的特征提取，通过结合聚类算法判断法矢夹角超过阈值的点为特征点，这种方法适用于简单规则形状的零件，但对于局部细微特征的识别精度较低。二是基于曲率的特征提取，依据特征点处曲率信息发生突变进行识别。这种方法虽然对噪声点具有较强的抗性，但提取的特征点里面会存在伪特征点，需要额外的滤波算法进行精简。三是基于某种数学模型或者算子改进，这种方法通过设计特定的求解流程来提高特征点识别精度，但仍然存在错误识别以及特征点遗漏的问题。

发明内容

本发明的目的是针对现有的点云关键特征提取存在的局部细微特征的识别精度较低且存在伪特征点以及特征点遗漏，需要额外的滤波算法进行精简的问题，提供一种无序点云关键特征点精确提取方法，目的在于避免真实特征点2环邻域中的非特征点被错误识别的问题，使用梯度变化准确识别角点、平面点和折痕点并且保证局部细节特征。

本发明的技术方案是：

一种无序点云关键特征精确识别与提取方法，其特征是首先对每个点云点进行邻域搜索，计算邻域距离和夹角对邻域点进行逆时针排序；然后对一环k邻域构建三角网格，计算顶点的离散法向，利用泰勒展开计算离散曲率；最后，对邻域距离、离散法向和离散曲率进行加权构造梯度矩阵并计算角点响应值，通过设定阈值识别特征点。具体步骤如下所述：

步骤1、计算每个点云点P_i的一环k邻域P_i-k。

(1)设定邻域k值取值范围[10,20]，每个数据点找到距离最近的k+5个邻域点，对邻域点沿着逆时针方向进行排序，计算邻近两邻域点与点云点P_i组成的向量夹角θ；

(2)计算空间距离最小的邻域点为起始点位，计算阈值角度θ₀＝360/k，沿逆时针方向当两邻域点的夹角θ小于阈值角度θ₀，删除模长较大的邻域点，另一个点标记为合格点进行保留，合格点作为新的起始点位。循环(2)直到邻域点减少到k个时结束。

步骤2、基于步骤1获得的邻域关系，构造三角化网格来计算每个点云点P_i的离散法向量，使用泰勒展开计算离散曲率，具体为：

(1)如图1所示，计算所有邻域三角形边长a/b/c，计算任一邻域三角形的形状权重Δ_j，公式如下：

(2)根据邻域三角形顶点坐标值计算法向量N_fj，根据(1)获得的形状权重，计算点云点Pi的离散法向N_p。

(3)如图2所示，根据(2)获得的法向量计算过点云点P_i的切平面F，对P_i与P_i-k所在的空间区域进行参数曲面S(u,v)拟合。构建P_i作为原点的局部坐标系xyz，曲面坐标系u/v投影到F面上对应x/y轴。计算邻域点P_i-k投影到切平面F的投影点Q_j-k，计算任一投影点Q_j-k与x轴夹角σ_j，Q_j-k与y轴夹角η_j。

(4)对曲面S进行泰勒展开，计算公式如下：

R_n是一个高阶无穷小，后续计算忽略此值。将上式写成矩阵形式，计算每个数据点P_i点的离散一阶{S_u,S_v}与二阶{S_uu,S_uv,S_vv}导数，计算公式如下：

其中，D_j＝|P_i-k-P_i|，Δu_j＝|Q_j-k-P_i|cosσ_j，Δv_j＝|Q_j-k-P_i|cosη_j。

(5)根据(4)获得的离散一阶导数，计算外恩加滕变换矩阵W，公式如下：

式中：E、F、G为曲面的第一类基本量，L、N、M₁、M₂为曲面的第二类基本量。

(6)计算外恩加滕变换矩阵W的特征值k₁,k₂，获得点云点P_i的离散高斯曲率κ，公式如下：

κ＝k₁·k₂

步骤3、对邻域距离、离散法向和离散曲率进行加权构造梯度矩阵并计算特征值，通过设定阈值识别特征点，具体如下：

(1)计算邻域点间距d，使用步骤2获得的离散法向和离散高斯曲率计算窗口函数w，公式如下：

其中，σ_k,σ_n,和σ_c是带宽参数，σ_k＝max(d)，σ_n＝180，σ_c＝max(κ_pi-κ_pi-k)。

(2)通过(1)获得的窗口函数，计算P_i沿着任一P_i-k点向量方向的梯度变化，公式如下

(3)将P_i的一阶和二阶离散导数代入上式，计算P_i在一环邻域内的整体梯度E，公式如下：

将上式进行改写，计算邻域点P_i-k在局部坐标系xyz下的方向余弦cosα，cosβ，cosγ，获得P_i点梯度的协方差矩阵M(P_i)，计算公式如下：

(4)如图3所示，计算梯度矩阵M(P_i)的三个特征值并进行排序，得到λ₁≥λ₂≥λ₃，设定阈值Ω，根据大小关系进行如下判断：

1)角点：λ₁≥λ₂≥λ₃≥Ω。

2)棱边点：λ₁≥λ₂，λ₁≥Ω，λ₃≠0；或者λ₂≥λ₁，λ₂≥Ω，λ₃≠0。

3)平面点：λ₁≤Ω，λ₂≤Ω，λ₃≤Ω。

所述的曲面S(u,v)需要进行最小二乘求解获得参数方程，具体过程如下：

(1)曲面方程一般表达式记为：

z＝a₁x²+a₂xy+a₃y²+a₄x+a₅y+a₆

(2)将P_i点以及任意5个邻域点P_i-k分别带入上式中，获得六元一次方程组进行求解，得到系数集合A＝{a1,...a6}。重复步骤(2)进行不超过N＝C_k ⁵次排列组合，获得N组系数集。

(3)计算所有邻域点P_i-k到每组系数拟合曲面的Z向投影距离之和，选择最小值所对应的曲面系数作为拟合结果。

通过计算角点响应函数R来简化特征值大小关系的判断过程，具体过程为：

(1)对梯度矩阵E(P_i)进行求导，获得角点响应值R，计算公式如下：

(2)根据Ω的取值范围为[0.01，0.1]，判断角点相应值R与阈值Ω的大小关系如下：

1)当R>Ω：P_i点角点。

2)当R<-Ω：P_i点为折痕点。

3)当|R|<Ω：P_i点为平面点。

本发明的有益效果是：

一、本发明使用梯度的变化进行特征点的提取，如图4所示，相比较其他方法在局部细节上具有更高的识别精度，局部细节特征点识别率在97％以上。

二、本发明对于阶差、下陷以及倒角类特征点(特征表现为多条轮廓线，且距离间隔最小可达0.5mm)的识别精度在0.15mm(实测点云平均点间距的1.5倍)以下，能够比较准确剔除邻近两条轮廓之间的非特征点，无需额外的滤波算法。

附图说明

图1一环邻域三角化网格构建示意图。

图2一环邻域局部坐标系构建及邻域点投影示意图。

图3梯度矩阵提取特征点示意图。

图4本发明方法与其他特征提取算法性能比较图。

图5本发明方法操作流程图。

具体实施方式

下面结构附图和实施例对本发明作进一步的说明。

如图1-5所示。

本发明所适用的邻域三角网格拓扑请参照图1所示，局部坐标系的构建过程及一环邻域点投影方向请参照图2，计算梯度矩阵获得的特征值来提取特征点的过程请参照图3。

一种无序点云关键特征精确识别与提取方法，如图5所示，它包括以下步骤。

步骤1、计算每个点云点P_i的一环k邻域P_i-k。

(1)设定邻域k值取值范围[10,20]，每个数据点找到距离最近的k+5个邻域点，对邻域点沿着逆时针方向进行排序，计算邻近两邻域点与点云点P_i组成的向量夹角θ；

(2)计算空间距离最小的邻域点为起始点位，计算阈值角度θ₀＝360/k，沿逆时针方向当两邻域点的夹角θ小于阈值角度θ₀，删除模长较大的邻域点，另一个点标记为合格点进行保留，合格点作为新的起始点位。循环(2)直到邻域点减少到k个时结束。

步骤2、基于步骤1获得的邻域关系，构造三角化网格来计算每个点云点P_i的离散法向量，使用泰勒展开计算离散曲率，具体为：

(1)如图1所示，计算所有邻域三角形边长a/b/c，计算任一邻域三角形的形状权重Δ_j，公式如下：

(2)根据邻域三角形顶点坐标值计算法向量N_fj，根据(1)获得的形状权重，计算点云点Pi的离散法向N_p。

(3)如图2所示，根据(2)获得的法向量计算过点云点P_i的切平面F，对P_i与P_i-k所在的空间区域进行参数曲面S(u,v)拟合。

所述的曲面S(u,v)可通过最小二乘求解获得参数方程，具体过程如下：

(1)曲面方程一般表达式记为：

z＝a₁x²+a₂xy+a₃y²+a₄x+a₅y+a₆

(2)将P_i点以及任意5个邻域点P_i-k分别带入上式中，获得六元一次方程组进行求解，得到系数集合A＝{a1,...a6}。重复步骤(2)进行不超过N＝C_k ⁵次排列组合，获得N组系数集。

(3)计算所有邻域点P_i-k到每组系数拟合曲面的Z向投影距离之和，选择最小值所对应的曲面系数作为拟合结果。

构建P_i作为原点的局部坐标系xyz，曲面坐标系u/v投影到F面上对应x/y轴。计算邻域点P_i-k投影到切平面F的投影点Q_j-k，计算任一投影点Q_j-k与x轴夹角σ_j，Q_j-k与y轴夹角η_j。

(4)对曲面S进行泰勒展开，计算公式如下：

R_n是一个高阶无穷小，后续计算忽略此值。将上式写成矩阵形式，计算每个数据点P_i点的离散一阶{S_u,S_v}与二阶{S_uu,S_uv,S_vv}导数，计算公式如下：

其中，D_j＝|P_i-k-P_i|，Δu_j＝|Q_j-k-P_i|cosσ_j，Δv_j＝|Q_j-k-P_i|cosη_j。

(5)根据(4)获得的离散一阶导数，计算外恩加滕变换矩阵W，公式如下：

(6)计算外恩加滕变换矩阵W的特征值k₁,k₂，获得点云点P_i的离散高斯曲率κ，公式如下：

κ＝k₁·k₂

步骤3、对邻域距离、离散法向和离散曲率进行加权构造梯度矩阵并计算特征值，通过设定阈值识别特征点，具体如下：

(1)计算邻域点间距d，使用步骤2获得的离散法向和离散高斯曲率计算窗口函数w，公式如下：

其中，σ_k,σ_n,和σ_c是带宽参数，σ_k＝max(d)，σ_n＝180，σ_c＝max(κ_pi-κ_pi-k)。

(2)通过(1)获得的窗口函数，计算P_i沿着任一P_i-k点向量方向的梯度变化，公式如下

(3)将P_i的一阶和二阶离散导数代入上式，计算P_i在一环邻域内的整体梯度E，公式如下：

将上式进行改写，计算邻域点P_i-k在局部坐标系xyz下的方向余弦cosα，cosβ，cosγ，获得P_i点梯度的协方差矩阵M(P_i)，计算公式如下：

(4)如图3所示，计算梯度矩阵M(P_i)的三个特征值并进行排序，得到λ₁≥λ₂≥λ₃，设定阈值Ω，根据大小关系进行如下判断：

1)角点：λ₁≥λ₂≥λ₃≥Ω。

2)棱边点：λ₁≥λ₂，λ₁≥Ω，λ₃≠0；或者λ₂≥λ₁，λ₂≥Ω，λ₃≠0。

3)平面点：λ₁≤Ω，λ₂≤Ω，λ₃≤Ω。

具体实施时还可通过计算角点响应函数R来简化特征值大小关系的判断过程，具体过程为：

(1)对梯度矩阵E(P_i)进行求导，获得角点响应值R，计算公式如下：

(2)根据Ω的取值范围为[0.01，0.1]，判断角点相应值R与阈值Ω的大小关系如下：

1)当R>Ω：P_i点角点。

2)当R<-Ω：P_i点为折痕点。

3)当|R|<Ω：P_i点为平面点。

本发明未涉及部分与现有技术相同或可采用现有技术加以实现。

Claims (4)

1.一种无序点云关键特征精确识别与提取方法，其特征是首先对每个点云点进行邻域搜索，同时计算邻域距离和夹角对邻域点进行逆时针排序；然后对一环k邻域构建三角网格，计算顶点的离散法向，利用泰勒展开计算离散曲率；最后，对邻域距离、离散法向和离散曲率进行加权构造梯度矩阵并计算角点响应值，通过设定阈值识别特征点。

2.根据权利要求1所述的方法，其特征在于它包括以下步骤：

步骤1、计算每个点云点P_i的一环k邻域P_i-k。

(1)设定邻域k值取值范围[10,20]，每个数据点找到距离最近的k+5个邻域点，对邻域点沿着逆时针方向进行排序，计算邻近两邻域点与点云点P_i组成的向量夹角θ；

(2)计算空间距离最小的邻域点为起始点位，计算阈值角度θ₀＝360/k，沿逆时针方向当两邻域点的夹角θ小于阈值角度θ₀，删除模长较大的邻域点，另一个点标记为合格点进行保留，合格点作为新的起始点位；循环(2)直到邻域点减少到k个时结束；

步骤2、基于步骤1获得的邻域关系，构造三角化网格来计算每个点云点P_i的离散法向量，使用泰勒展开计算离散曲率，具体为：

(1)计算所有邻域三角形边长a/b/c，计算任一邻域三角形的形状权重Δ_j，公式如下：

(2)根据(1)获得的形状权重，计算点云点P_i的离散法向N_p；

式中：N_fj为邻域三角形的法向量，计算邻域三角形任意两顶点构成的向量叉乘获得。

(3)根据(2)获得的法向量计算过点云点P_i的切平面F，对P_i与P_i-k所在的空间区域进行参数曲面S(u,v)拟合，构建P_i作为原点的局部坐标系xyz，x轴的方向选择为最短模长|P_i-P_i-k|的向量方向。

(4)对曲面S进行泰勒展开，计算公式如下：

式中u/v为曲面的曲纹坐标方向，投影到F面上对应x/y轴。R_n是一个高阶无穷小，后续计算忽略此值。将上式写成矩阵形式，计算每个数据点P_i点的离散一阶{S_u,S_v}与二阶{S_uu,S_uv,S_vv}导数，计算公式如下：

其中，D_j＝|P_i-k-P_i|，Δu_j＝|Q_j-k-P_i|cosσ_j，Δv_j＝|Q_j-k-P_i|cosη_j；式中Q_j-k为邻域点P_i-k投影到切平面F的投影点，σ_j为任一投影点Q_j-k与x轴夹角，η_j为任一投影点Q_j-k与y轴夹角。

(5)根据(4)获得的离散一阶导数，计算外恩加滕变换矩阵W，公式如下：

式中：E、F、G为曲面的第一类基本量，L、N、M₁、M₂为曲面的第二类基本量。

(6)计算外恩加滕变换矩阵W的特征值k₁,k₂，获得点云点P_i的离散高斯曲率κ，公式如下：

κ＝k₁·k₂

步骤3、对邻域距离、离散法向和离散曲率进行加权构造梯度矩阵并计算特征值，通过设定阈值识别特征点，具体如下：

(1)计算邻域点间距d，使用步骤2获得的离散法向和离散高斯曲率计算窗口函数w，公式如下：

其中，σ_k,σ_n,和σ_c是带宽参数，σ_k＝max(d)，σ_n＝180，σ_c＝max(κ_pi-κ_pi-k)。

(2)通过(1)获得的窗口函数，计算P_i沿着任一P_i-k点向量方向的梯度变化，公式如下

(3)将P_i的一阶和二阶离散导数代入上式，计算P_i在一环邻域内的整体梯度E，公式如下：

将上式进行改写，获得P_i点梯度的协方差矩阵M(P_i)，计算公式如下：

式中cosα，cosβ，cosγ分别对应邻域点P_i-k在局部坐标系x/y/z下的方向余弦,Sx，Sy,Sz为曲面S的一阶偏导数。

(4)计算梯度矩阵M(P_i)的三个特征值并进行排序，得到λ₁≥λ₂≥λ₃，设定阈值Ω，根据大小关系进行如下判断：

1)角点：λ₁≥λ₂≥λ₃≥Ω；

2)棱边点：λ₁≥λ₂，λ₁≥Ω，λ₃≠0；或者λ₂≥λ₁，λ₂≥Ω，λ₃≠0；

3)平面点：λ₁≤Ω，λ₂≤Ω，λ₃≤Ω。

3.根据权利要求2所述的方法，其特征是所述的曲面S(u,v)需要进行最小二乘求解获得参数方程：

(1)曲面方程一般表达式记为：

z＝a₁x²+a₂xy+a₃y²+a₄x+a₅y+a₆

(2)将P_i点以及任意5个邻域点P_i-k分别带入上式中，获得六元一次方程组进行求解，得到系数集合A＝{a1,...a6}；重复步骤(2)进行不超过N＝C_k ⁵次排列组合，获得N组系数集；

(3)计算所有邻域点P_i-k到每组系数拟合曲面的Z向投影距离之和，选择最小值所对应的曲面系数作为拟合结果。

4.根据权利要求2所述的方法，其特征是通过计算角点响应函数R来简化特征值大小关系的判断过程：

(1)对梯度矩阵E(P_i)进行求导，获得角点响应值R，计算公式如下：

(2)根据Ω的取值范围为[0.01，0.1]，判断角点相应值R与阈值Ω的大小关系如下：

1)当R>Ω：P_i点角点。

2)当R<-Ω：P_i点为折痕点。

3)当|R|<Ω：P_i点为平面点。