一种三维点云相似配准方法

本发明属于三维点云数据处理技术领域，具体涉及一种三维点云相似配准方法，特别适用于三维点云数据中含有大量噪声点和外点情况下的相似配准。

背景技术：

三维图像采集设备一次只能得到物体一个侧面的三维数据，为了得到物体的整体三维数据，需要从多个角度获取三维数据，并对这些三维数据进行配准。对于三维点云配准问题，待配准的点云之间不仅存在着旋转和平移变换，往往还包括了尺度变换，这就是相似配准问题。由于相似配准问题大量地存在于医学图像、卫星遥感图像等数据中，因此，相似配准是一类非常重要的配准问题。现有的相似配准方法大多存在应用局限或精度不佳的问题，如ying等人将尺度因子合并到迭代最近点(iterativeclosestpoint，icp)算法中，将配准问题转化为7d非线性空间上的约束优化问题，再采用奇异值分解(singularvaluedecomposition，svd)方法来迭代求解优化问题。然而该算法难以应用于存在大量噪声和外点的实际点云数据的相似配准中。

现有的三维点云相似配准算法多数存在配准精度不佳和应用的局限性，如对于含有大量噪声和外点的三维点云其配准精度较差。由于点云数据采集设备的物理限制和噪声等因素影响，所采集点云中存在大量噪声和外点，从而影响点云的相似配准精度，甚至获得错误的配准结果。

技术实现要素：

为了解决现有技术中存在的上述问题，本发明提供了一种三维点云相似配准方法。本发明要解决的技术问题通过以下技术方案实现：

本发明提供了一种三维点云相似配准方法，包括：

s1：获得形状点云和模型点云，建立相似配准优化模型并对参数进行初始化；

s2：利用所述相似配准优化模型获得尺寸因子初值、刚体变换中旋转矩阵的初值和平移向量的初值；

s3：利用所述尺寸因子获取多个取值，并利用每个取值进行刚体配准求解和迭代，获得每个取值对应的最优旋转矩阵和平移向量；

s4：根据每个取值对应的最优旋转矩阵和平移向量，进行尺度因子精细化迭代，获取最优尺度因子及最终的旋转矩阵和平移向量。

在本发明的一个实施例中，所述s1包括：

s11：获取形状点云和模型点云其中，np为所述形状点云中的点数，nq为所述模型点云的点数；

s12：根据所述形状点云和所述模型点云，建立基于伪huber损失函数的三维点云相似配准优化模型：

s.t.r^tr＝i3，det(r)＝1

其中，s为尺度因子，r为旋转矩阵，t为平移向量，pi为所述形状点云中第i个点的三维坐标，mc(i)为形状点云中第i个点pi在模型点云中的最优对应点，b为外点阈值，i3表示3阶的单位矩阵；

s13：设定刚体配准最大迭代次数kmax、尺度因子精细化最大迭代次数nmax、刚体变换后点云间均方误差最小值εmin，以及相似变换均方误差最小值εsmin，并且令刚体配准迭代次数初始值k＝1，尺度因子精细化迭代次数初始值n＝1。

在本发明的一个实施例中，所述s2包括：

s21：对所述模型点云m进行相似变换，获得相似变换后的形状点云p：

其中，s为尺度因子，r为旋转矩阵，t为平移向量；

s22：利用模型点云m及其相似变换后的形状点云p计算获得尺寸因子的初值：

其中，pg和mg分别为形状点云p和模型点云m的重心，pw和mw分别为形状点云p和模型点云m的质心；

s23：利用主成分分析法计算获得所述旋转矩阵r的初值r0以及所述平移向量t的初值t0。

在本发明的一个实施例中，所述s3包括：

s31：利用表达式令n＝1，根据所述尺寸因子初值s0获得多个取值

s32：对于所述多个取值进行所述形状点云p和所述模型点云m的刚性变换求解，获得每个取值对应的最优旋转矩阵和最优平移向量。

在本发明的一个实施例中，所述s32包括：

s321：根据当前取值计算第k次迭代的相似变换后的形状点云p与模型点云m的空间对应关系满足：

其中，i＝1,2,…,np，k的初始值为1，ck(i)表示第k次迭代中模型点云中的第ck(i)个点，mj为pi在模型点云m中对应的点；

s322：根据第k次迭代后获得的空间对应关系计算形状点云p与模型点云m在两次迭代之间的相对刚体变换r^*和t^*，计算公式为：

其中，r^*表示相邻两次迭代之间的相对旋转矩阵，t^*表示相邻两次迭代之间的相对平移矩阵；

s323：利用所述相对刚体变换r^*和t^*更新第k次迭代的刚体变换rk和tk：

rk＝r^*rk-1，tk＝r^*tk-1+t^*；

s324：计算第k次迭代的刚体变换rk和tk后形状点云q与模型点云m的均方误差：

s325：判断εk≤εmin或k≥kmax是否成立，若否，则令k+1并返回至步骤s321进行迭代，若是，则确定当前刚体变换矩阵rk和tk为相似变换后的形状点云与模型点云的最优刚体变换，其中，εmin为刚体变换后形状点云与模型点云的均方误差最小值，kmax为刚体配准最大迭代次数；

s326：重复步骤s321至步骤s325，获得剩余取值对应的最优旋转矩阵和最优平移向量。

在本发明的一个实施例中，所述s4包括：

s41：计算每一组相似变换的均方误差：

其中，表示求二范数的平方；

s42：在中找出获得均方误差εsl最小时的尺寸因子slmin及对应的均方误差εsl、旋转矩阵和平移向量；

s43：判断εsl≤εsmin或n≥nmax是否成立，若否，则返回至s31，令n+1，s0为步骤s42中均方误差εsl最小时的尺寸因子slmin，进行迭代，若是，则输出当前迭代所求尺度因子slmin及其对应的旋转矩阵和平移向量。

本发明的另一方面提供了一种存储介质，所述存储介质中存储有计算机程序，所述计算机程序用于执行上述实施例中任一项所述的三维点云相似配准方法的步骤。

本发明的又一方面提供了一种电子设备，包括存储器和处理器，所述存储器中存储有计算机程序，所述处理器调用所述存储器中的计算机程序时实现如上述实施例中任一项所述三维点云相似配准方法的步骤。

与现有技术相比，本发明的有益效果在于：

本发明的三维点云相似配准方法，将伪huber损失函数引入三维点云相似配准优化模型中，由于伪huber损失函数对于噪声和外点均具有抑制作用，因此能够有效提高相似配准准确性。

以下将结合附图及实施例对本发明做进一步详细说明。

附图说明

图1是本发明实施例提供的一种三维点云相似配准方法的流程框图；

图2是本发明实施例提供的一种三维点云相似配准方法的详细流程图。

具体实施方式

为了进一步阐述本发明为达成预定发明目的所采取的技术手段及功效，以下结合附图及具体实施方式，对依据本发明提出的一种三维点云相似配准方法进行详细说明。

有关本发明的前述及其他技术内容、特点及功效，在以下配合附图的具体实施方式详细说明中即可清楚地呈现。通过具体实施方式的说明，可对本发明为达成预定目的所采取的技术手段及功效进行更加深入且具体地了解，然而所附附图仅是提供参考与说明之用，并非用来对本发明的技术方案加以限制。

应当说明的是，在本文中，诸如第一和第二等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来，而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。而且，术语“包括”、“包含”或者任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的物品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素。在没有更多限制的情况下，由语句“包括一个……”限定的要素，并不排除在包括所述要素的物品或者设备中还存在另外的相同要素。

请参见图1和图2，图1是本发明实施例提供的一种三维点云相似配准方法的流程框图；图2是本发明实施例提供的一种三维点云相似配准方法的详细流程图。该三维点云相似配准方法包括：

s1：获得形状点云和模型点云，建立相似配准优化模型并对参数进行初始化；

步骤s1具体包括：

s11：获取形状点云和模型点云其中，np为形状点云p的点数，nm为模型点云m的点数，np∈n，nm∈n，n表示自然数集合；

s12：根据所述形状点云p和所述模型点云m，建立基于伪huber损失函数的三维点云相似配准优化模型。

三维点云相似配准的目的是建立两个点云之间的空间对应关系并寻找两个点云之间的最优相似变换。在求取尺度因子之后，相似配准问题转换成为刚体配准问题。通常三维点云刚体配准包括两个内容：(1)建立两个点云之间的对应关系；(2)求解两个点云之间的刚体变换。

给定形状点云和模型点云假设从形状点云p到模型点云m存在某种映射关系c:p→m，以及从形状点云p到模型点云m的相似变换关系t，则可以定义如下的相似性度量函数：

j(t(p),c(p))，

点云相似配准过程可以看成如下的优化问题：

其中，s为尺度因子，r为旋转矩阵，t为平移向量。

本实施例引入对噪声和外点不敏感的伪huber损失函数建立鲁棒的三维点云相似配准优化模型，从而来抑制噪声和外点对配准结果的影响，提高相似配准的精度。

伪huber损失函数是huber函数的平滑版本，该函数具有连续可导性，能够有效抑制噪声和外点的影响，其函数定义为：

其中，b为外点阈值，a为误差。对于较小的a值，该损失函数近似值为a/2值，而对于较大的a值，该损失函数可近似为一条斜率为b的直线，因此其对噪声和外点不敏感。

对于含噪声和外点的三维点云相似配准问题，其噪声和外点的部分可看成异常值，利用伪huber损失函数建立相似配准的目标函数能够有效降低噪声和外点对配准过程的影响，从而提高配准的效率和精度。因此，对于含有噪声和外点的两个点云p和m，可建立如下的基于伪huber损失函数的相似配准优化模型：

s.t.r^tr＝i3，det(r)＝1

其中，s为尺度因子，r为旋转矩阵，t为平移向量，pi为所述形状点云中第i个点的三维坐标，mc(i)为形状点云中第i个点pi在模型点云中的最优对应点，b为外点阈值，i3表示3阶的单位矩阵。

由此，三维点云相似配准问题转化为一个非线性目标函数的最小化优化问题。

s13：设定刚体配准最大迭代次数kmax、尺度因子精细化最大迭代次数nmax、刚体变换后点云间均方误差最小值εmin，以及相似变换均方误差最小值εsmin，并且令刚体配准迭代次数初始值k＝1，尺度因子精细化迭代次数初始值n＝1。

s2：利用所述相似配准优化模型获得尺寸因子初值、刚体变换中旋转矩阵的初值和平移向量的初值；

在本实施例中，步骤s2具体包括：

s21：对所述模型点云m进行相似变换，获得相似变换后的形状点云p：

其中，s为尺度因子，r为旋转矩阵，t为平移向量，nm为模型点云m中的点数；

s22：利用模型点云m及其相似变换后的形状点云p计算获得尺寸因子的初值。

在点云的相似变换中，对于同一个点云，其重心和质心的相对距离是成比例变化的，在不同尺度因子的相似变换下，这个相对距离比就是尺度因子。由于噪声和外点的影响，重心和质心的位置会发生变化，结果并不精确，因此本实施例先通过计算三维点云的重心和质心的距离比，来得到尺度因子的初值，再通过迭代方式对尺度因子进行精细化求解。

给出一个模型点云m，其经过相似变换后的点云为p，假设形状点云p和模型点云m重心分别为pg和mg，质心分别为pw和mw，可以得到：

其中，ri是三维点云中每个点到质心的距离。

则尺度因子s的初值为：

s23：利用主成分分析法计算获得所述旋转矩阵r的初值r0以及所述平移向量t的初值t0。

具体地，本实施例采用主成分分析(principalcomponentanalysis，pca)来求解刚体变换的初值。首先计算两个待配准三维点云的协方差矩阵和均值，使用协方差矩阵求得特征向量，然后以均值为原点，特征向量为坐标轴建立两个待配准点云的三维坐标系，通过将三维坐标系对齐来求得刚体变换中旋转矩阵r的初值r0和平移向量t的初值t0。

s3：利用所述尺寸因子获取多个取值，并利用每个取值进行刚体配准求解和迭代，获得每个取值对应的最优旋转矩阵和平移向量；

在求解三维点云相似配准优化模型时，本实施例采用由粗到精的方式，迭代求解尺度因子s。通过给定尺度因子对应的多个取值，对每一个取值进行刚体变换求解，通过误差比较来求取最优相似变换。在给定尺度因子s后，相似配准问题转换为刚体配准问题，因此相似配准优化模型求解，重点是对于刚体变换的求解。给定刚体变换初值r0和t0后，本实施例通过迭代方法求解最优刚体变换。

步骤s3具体包括：

s31：利用表达式令n＝1，根据所述尺寸因子初值s0获得多个取值

在开始尺度因子精细化迭代时，根据尺度因子s的初值s0获得一组按照等差数列设置的取值，其中第一个取值为最后一个取值为相邻取值的差值为从而获得m个取值。

例如，在第一次迭代时，n＝1，此时第一个取值为0，最后一个取值为2s0，相邻取值的差值为从而获得20个以等差数列方式排列的数值。

s32：对于所述多个取值进行所述形状点云q和所述模型点云m的刚性变换求解，获得每个取值对应的最优旋转矩阵和最优平移向量。

具体地，步骤s32包括：

s321：根据当前取值计算第k次迭代的相似变换后的形状点云p与模型点云m的空间对应关系满足：

其中，i＝1,2,…,np，k的初始值为1，ck(i)表示第k次迭代中模型点云中的第ck(i)个点，mj为pi在模型点云m中对应的点；

具体地，对于上述公式中形状点云p与模型点云m的空间对应关系的求解，本实施例采用基于delaunay三角化的最近点搜索算法来实现。将三维空间中任意分布的散乱点用直线段连接起来，形成空间上既不重叠又无间隙的紧邻的四面体集，每个四面体的顶点为原始散乱点。将形状点云p进行delaunay三角剖分，delaunay剖分必须满足两个基本准则，即空圆特性和最大化最小角特性。在此基础上，采用三角划分搜索策略找到模型点云在相似变换后的形状点云中的对应点。

s322：根据第k次迭代后获得的空间对应关系计算形状点云p与模型点云m在两次迭代之间的相对刚体变换r^*和t^*，计算公式为：

其中，r^*表示相邻两次迭代之间的相对旋转矩阵，t^*表示相邻两次迭代之间的相对平移矩阵。

lm(levenbergmarquarelt)算法是一种迭代求函数极值的算法，算法收敛速度快。本实施例采用lm算法对上式进行优化求解,利用lm算法可求解获得形状点云p与模型点云m在两次迭代之间的相对刚体变换r^*和t^*。

s323：利用所述相对刚体变换r^*和t^*更新第k次迭代的刚体变换rk和tk：

rk＝r^*rk-1，tk＝r^*tk-1+t^*；

s324：计算第k次迭代的刚体变换rk和tk后形状点云p与模型点云m的均方误差：

s325：判断εk≤εmin或k≥kmax是否成立，若否，则令k＝k+1并返回至步骤s321进行迭代，若是，则确定当前刚体变换矩阵rk和tk为相似变换后的形状点云与模型点云的最优刚体变换，其中，εmin为刚体变换后形状点云与模型点云的均方误差最小值，kmax为刚体配准最大迭代次数；

具体地，在上述步骤通过迭代求取值对应的最优旋转矩阵和最优平移向量时，迭代次数为k，初始值为1，经过k次迭代后，若εk≤εmin或k≥kmax成立，则确定当前刚体变换矩阵rk和tk为取值对应的最优旋转矩阵和最优平移向量。

s326：重复步骤s321至步骤s325，获得剩余取值对应的最优旋转矩阵和最优平移向量。

具体地，对于剩余取值按照步骤s321至步骤s325中的迭代过程，分别可获得各自对应的最优旋转矩阵和最优平移向量。

s4：根据每个取值对应的最优旋转矩阵和平移向量，进行尺度因子精细化迭代，获取最优尺度因子及对应的旋转矩阵和平移向量。

在本实施例中，步骤s4具体包括：

s41：计算每一组相似变换的均方误差：

其中，

s42：在中找出获得均方误差εsl最小时的取值sl及对应的均方误差εsl、旋转矩阵和平移向量；

s43：判断εsl≤εsmin或n≥nmax是否成立，若否，则返回至s31，令n+1，s0为步骤s42中均方误差εsl最小时的尺寸因子slmin，进行迭代，若是，则输出当前迭代所求尺度因子slmin及其对应的旋转矩阵和平移向量。

具体地，对当前尺度因子对应的所有取值进行相似变换求解后，获得每个取值对应的均方误差，找出获得均方误差εsl最小时的取值slmin及对应的均方误差εsl、旋转矩阵和平移向量，如果满足εsl≤εsmin或n≥nmax，则当前取值slmin对应的刚体变换和平移向量即为本实施例的三维点云相似配准方法所求的点云匹配的最终旋转矩阵和平移向量。

如果不满足εsl≤εsmin或n≥nmax，则返回至s31，令n＝n+1，s0替换为步骤s42中均方误差εsl最小时的尺寸因子slmin，进行迭代。例如，在尺度因子精细化迭代时，返回至s31此时，取n＝2，s0取slmin，进行迭代处理，直至满足εsl≤εsmin或n≥nmax。

综上，针对三维点云数据中含有大量噪声和外点情况下的相似配准问题，本发明实施例提出了一种基于伪huber损失函数的三维点云相似配准方法，由于伪huber损失函数对于噪声和外点均具有抑制作用，因此，将伪huber损失函数引入三维点云相似配准优化模型中，以提高三维点云相似配准的精度。此外，在建立三维点云相似配准优化模型的基础上，本发明实施例根据尺度一致性来求解初始尺度因子，将相似配准问题转换为刚体配准问题，再根据尺度因子取值区间迭代求解刚体变换矩阵和平移向量，最后给出最优刚体变换和尺度因子。该三维点云相似配准方法由于伪huber损失函数对于噪声和外点均具有抑制作用，因此能够有效提高相似配准准确性。

本发明的又一实施例提供了一种存储介质，所述存储介质中存储有计算机程序，所述计算机程序用于执行上述实施例中所述的三维点云相似配准方法步骤。本发明的再一方面提供了一种电子设备，包括存储器和处理器，所述存储器中存储有计算机程序，所述处理器调用所述存储器中的计算机程序时实现如上述实施例所述三维点云相似配准方法的步骤。具体地，上述以软件功能模块的形式实现的集成的模块，可以存储在一个计算机可读取存储介质中。上述软件功能模块存储在一个存储介质中，包括若干指令用以使得一台电子设备(可以是个人计算机，服务器，或者网络设备等)或处理器(processor)执行本发明各个实施例所述方法的部分步骤。而前述的存储介质包括：u盘、移动硬盘、只读存储器(read-onlymemory，rom)、随机存取存储器(randomaccessmemory，ram)、磁碟或者光盘等各种可以存储程序代码的介质。

以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明，不能认定本发明的具体实施只局限于这些说明。对于本发明所属技术领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干简单推演或替换，都应当视为属于本发明的保护范围。