CN1684105A - 大规模三维场景多视点激光扫描数据自动配准方法 - Google Patents

Abstract

本发明涉及的大规模三维场景多视点激光扫描数据自动配准方法，其步骤包括：获取数据，在被测场景前放置激光测距仪，调节激光测距仪使其Z轴垂直地面，X和Y轴平行于地面；对被测场景进行逐行逐列扫描，获得每一个视点的三维数据；相邻视点扫描数据要求保持10％～20％的重叠；提取被测物结构并对其数据进行特征提取；计算虚拟特征与构建特征单元，并进行粗配准精配准和全局配准，最后建立配准模型本发明可不借助任何外协设备和人工参与的方法，可解决大规模三维场景多视点激光扫描数据自动配准；具有普遍适用性，即适用于现代结构的场景又适用于复杂的大兴文物和古建。

Description

大规模三维场景多视点激光扫描数据自动配准方法

技术领域

本发明属于计算机信息处理领域的三维数据配准方法，特别涉及一种大规模三维场景多视点激光扫描数据自动配准方法。

背景技术

三维数据配准是数字化和反求工程的关键技术。特别是随着中远距离三维激光扫描技术的出现，简化了三维数据获取过程，促进了数字城市、数字导航、虚拟现实、城市规划、文物数字化保护等相关领域的发展，基于地面的三维激光扫描数据备受重视。为了获得大规模三维场景完整几何描述，需要在几十个甚至几百个视点位置采集激光三维扫描数据，故而存在将不同视点的激光三维扫描数据转换在同一坐标系的配准问题。

传统的多视点三维数据配准的研究多针对单一物体，最有代表性的是最近点迭代算法，即ICP算法该方法通过迭代在两视点激光扫描点集内寻找对应点，计算两视点的位置关系。不足之处在于需人工设定初始转换估计，易陷入局部最小，迭代时间长，不适合大规模配准。对于大规模场景常用配准方法有两种：靶标控制点法和GPS等外协设备法。靶标控制点法要求在相邻视点重叠区域摆放靶标，利用激光对靶标特殊的反射特性自动寻找靶标控制点，只要相邻视点的激光扫描数据中有三个以上的同名控制点，即可求出相邻视点的坐标转换关系。如果场景中存在与靶标具有相同反射特性的物质，会误识靶标，配准失败，靶标控制点法需要大量人工干预，且有些环境无法放置靶标。GPS等外协设备法主要是借助外协设备获得激光扫描仪在各视点的位姿参数，不但设备价格昂贵，而且位姿参数精度不高，配准结果与实际情况有差异。因此，急需一种不借助任何外协设备直接将各视点的三维激光扫描数据无缝拼合在同一坐标系的自动大规模三维场景多视点三维数据配准方法。经文献检索分析，发现I.Stamos在2003年发表的名为“基于特征的城市大规模场景距离数据自动配准(AutomatedFeature-Based Range Registration of Urban Scenes of Large Scale)”的文章中提出了采用特征直线实现相邻视点三维激光数据自动配准方法，并称直线特征减少了几何复杂性，适合大规模配准。该方法不足之处在于只采用线性特征不具有普遍性，对于线性特征较少的大型文物、古建等三维场景配准常常失败。此外，没能解决因遮挡等实际因素影响重叠区域缺乏有效线性特征的配准问题。

上述方法意在解决有重叠区域的两视点三维数据配准，至于几十个甚至几百个视点，必需采取有效的配准策略方可将其三维数据无锋拼合在同一坐标系。常用的配准策略有序列配准和同步配准。序列配准指一对视点配准后，另一对视点开始配准，且须包含前一对视点中的一个视点，环环相扣直到遍历所有视点为止。序列配准策略的优点在于同一时间只有一对视点参加配准，占用内存少；缺点在于存在较大的累积误差，特别是闭环形场景，首末视点间出现缝隙。同步配准则是所有视点同时配准，不会出现累积误差，配准精度高，但计算量大，对计算设备要求高。序列配准和同步配准均不适合大规模场景。

发明内容

本发明的目的是提供一种不借助任何外借设备直接将大规模三维场景多视点激光扫描数据自动无缝拼合在同一坐标系的配准方法。

本发明的技术方案如下：

本发明提供的大规模三维场景多视点激光扫描数据自动配准方法，包括如下步骤：

1)数据获取

由于大规模场景尤其室外场景，不仅包含主要建筑结构，而且含有树木、行人车辆等遮挡物，以及激光穿过透明物体形成的散乱点。遮挡物位于三维激光扫描仪与主要建筑结构之间，而激光穿过透明物体形成散乱点位于主要建筑结构之后。在测量时，在被测场景前放置激光测距仪，调整三维激光扫描仪，使其Z轴垂直地面，X和Y轴与地面平行；对被测场景进行逐行逐列扫描，获得每一个视点的三维数据；相邻视点扫描数据保持10％～20％的重叠；

2)提取主要建筑结构

输入每一视点对应的激光三维扫描数据，计算各垂直扫描线上采样点与激光测距仪的水平距离，并根据水平距离，利用哈夫变换检测垂直线段；以其所含的点数计算垂直线的长度，并选择最长垂直线段对应的水平距离为基准，基准之前δ₁米和基准之后δ₂米范围内的数据点为主要建筑结构；所述δ₁为1米，δ₂为0.5米；

因为激光测距仪对场景是逐行逐列采样，对于每一垂直扫描线，其上属于同一垂直面的点到激光测距仪的水平距离相等，且主要建筑结构一般与地面垂直本发明逐一分析每一垂直扫描线，根据采样点到激光测距仪的水平距离，采用统计分布的方法确定主要建筑结构的区域范围，将主要建筑结构提取出来；

3)对保留的主要建筑结构的数据点进行特征提取

三维特征反映了场景内在的结构特征，简化了场景的几何复杂度，是多视点三维数据配准的有效基元。本发明首先采用加权局部曲线拟合分别处理水平扫描线和垂直扫描线，计算扫描线各点的一次微分；接着根据点的上下两邻域点的微分符号和点与其上下两邻域点的距离检测场景边缘点；然后利用链码跟踪技术形成不同的边缘链，并将同一边缘链的点拟合成首尾相连的一组线段，称为特征线段。相邻特征线段的交点为特征点；

a)扫描线上各点微分

对于每条扫描线上的采样点用(1)式进行局部二次拟合，并以该点沿扫描线方向上下k个邻近点为拟合区间[-K，K]，K为3-5；

$\hat{d}\left(i\right)=a_0+a_{1}i+a_2{i}^2---\left(1\right)$

式中：i-扫描线上采样点排序，i＝0，1，…，n；

第i点到激光测距仪的计算水平距离；

a₀，a₁，a₂-拟合方程(1)式的拟合系数；

采用自适应加权拟合方法，求取拟合系数a₀，a₁，a₂，且加权拟合误差函数为

$E^{\′}\left(i\right)=\mathrm{\Σw}\left(i\right)\{d\left(i\right)-\hat{d}\left(i\right){\}}^2---\left(2\right)$

式中：d(i)-第i点到激光测距仪的实际水平距离；

w(i)-扫描线上第i点的权值，代表点的可靠性；采用迭代权值修正法，逐步消除局部点的影响，使拟合仅在可靠点内进行；具体方法如下：假设已知第m次迭代中系数为a₀ ^m，a₁ ^m，a₂ ^m，求第k+1次迭代中系数为a₀ ^m+1，a₁ ^m+1，a₂ ^m+1；首先计算第m次迭代的拟合残余量和第m次迭代的权值；

$e^k\left(i\right)={(d\left(i\right)-{\hat{d}}^m\left(i\right))}^2---\left(3\right)$

$w^{m+1}\left(i\right)=\left\{\begin{array}{cc}1& e^m\left(i\right)\≤b_s\\ \frac{b_s}{e^m\left(i\right)}& e^m\left(i\right)>b_s\end{array}\right.---\left(4\right)$

式中：b_s-误差允许阈值；

初始权值w⁰(i)＝1，用(3)、(4)、(5)式不断修正拟合系数和权值，直至拟合曲线满足误差要求且包含最多点数为止；标记参与拟合的点；

然后，利用(1)式计算扫描线上有标记各点的一次微分；

$\frac{\∂\hat{d}\left(i\right)}{\∂i}=a_1+2a_{2}i---\left(5\right)$

b)检测边缘点

对于扫描线上有标记的各点，考察其沿扫描线方向上下有标记两点的一次微分符号，选择异号的点为法向不连续点，并标记为0；

对于扫描线上有标记的各点，计算其与沿扫描线方向上下有标记两点距离，考察其中最大距离，选择最大距离大于阈值ε为深度不连续点，并标记为1；

c)特征点线提取

首先采用链码跟踪方法将边缘点区分成不同的边缘链，并将每一边缘链上的点按顺时针排序，然后通过分裂合并的方法将每一边缘链拟合成一组首尾相接的特征线段，具体方法如下：

·初始化

将边缘链分为等长的若干段L₀，L₁，…，L_n

·分裂

利用最小二乘法对当前每一段L_i拟合成直线选择离该拟合直线最远且最小残差E_i大于允许阈值T的一点为基准，将本段边缘链断开为两段，并不断分裂至满足拟合误差为止；依次处理其它各段；

·合并

选择拟合残差小于阈值T的相邻两段合并成一段为特征线段，特征线段的端点为特征点；

4)计算虚特征与构造特征单元

特征线对应场景的边缘信息，特征点对应场景的角点或尖点，它们构成了场景基本结构框架。然而，由于真实场景不可避免地存在遮挡现象，如果遮挡严重，在两视点的重叠区内缺乏有效特征，导致配准失败。本发明提出虚特征概念，将实际不相交的特征线段延长相交形成虚特征线段和虚特征点，并加入上述提取的特征线段和特征点中，虚特征说明了特征之间的几何位置关系；

本发明将特征点连同其相邻的两特征线段嵌套成特征单元PU{P|L₁，L₂：<n₁，n₂，n₃，p>}。特征单元包含丰富的几何信息，即三个方向矢量和一个位置矢量：两个特征线段L₁和L₂的方向矢量(n₁，n₂)，矢量方向由特征点指向线段，第三个矢量为前两矢量的合成矢量，反映了相邻两特征线段构成平面的法向n₃＝n₁×n₂，×为矢量叉乘位置矢量为特征点p；

5)配准

局部配准指带有重叠区域的相邻视点两两配准。本发明采用由粗到精的配准方式，使两视点数据完全融合；

先进行粗配准：利用特征单元将两视点的三维数据自动拼合在同一坐标系一对特征单元便可构造一个空间转换关系：从对应的特征线的方向及其合成法向计算旋转矩阵；从对应的特征点对计算平移向量。给定一对分别属于视点S₁和视点S₂的特征单元PU{P|L₁，L₂：<n₁，n₂，n₃，p>}∈ S₁和PU′{P′|L₁′，L₂′：<n₁′，n₂′，n₃′，p′>}∈S₂，旋转矩阵由n_i和n_i′(i＝1，2，3)计算，即使下式误差函数最小，并采用四元素法求出旋转矩阵；

$\mathrm{Err}=\underset{i=1}{\overset{3}{\mathrm{\&Sigma;}}}{\left|\right|n}_i^{\&prime;}-R{n}_i{\left|\right|}^2---\left(6\right)$

设p和p′在各自局部坐标系的坐标为(x，y，z)与(x′，y′，z′)，计算向量；

$t=R\left[\begin{array}{c}x\\ y\\ z\end{array}\right]-\left[\begin{array}{c}{x}^{\&prime;}\\ y^{\&prime;}\\ z^{\&prime;}\end{array}\right]---\left(7\right)$

遍历所有特征单元构造N个候选转换估计，采用匹配程度确定出最佳配准转换估计。匹配程度以得到匹配的特征单元数目计算。特征单元匹配的约束条件如下：

·位置约束

p′是视点S₂中与视点S₁中p点最近的特征点，且p和p′间的距离小于阈值Δ₁；

·方向约束

L_i与L_i′两线段所在直线间的距离和方向夹角分别小于Δ₂和Δ₃；

再进行精配准：经过粗配准之后，带有重叠区域的两视点的三维数据基本拼合在同一坐标系内，只在某些细节部分配准质量稍差些，需进一步精配准。具体方法：1)在粗配准的基础上，确定重叠区。2)将重叠区划分成N块，每块含u×v个点。3)计算每块对应点的平均差异，找出差异较大的M块作为适配区。在适配区内利用ICP算法进行精配准，可大大减少运算时间。4)利用了激光采样点的激光反射率和点位信息，结合马氏距离寻找对应点。

对于视点超过2个以上的多视点三维场景，所述步骤还包括全局配准：

全局配准指将所有视点三维数据拼合在同一坐标系的策略。首先所有视点两两配准，计算匹配程度；再以视点为结点，以匹配程度为连结权，生成无向带权图；然后利用最小生成树原则建立全局配准模型图。由全局配准模型图确定一个固定视点S_a，并使其它视点的数据向其转化。其中建立全局配准模型的具体方法如下：

设第i对有重叠区的两视点计算出的转换矩阵为T_i＝[R_i，t_i]和对应的局部匹配程度为g(T_i)；以视点为结点，以g(T_i)为结点连接权重，生成一个包含所有视点在内的无向带权连通图G＝<V，E>，V为结点，E为边，且e_i＝g(T_i)∈E，并采用逐步短接法求出最小生成树，建立全局配准模型。

本发明的技术效果为：

本发明不借助任何外协设备和人工参与的方法，解决了大规模三维场景多视点激光扫描数据自动配准。其优势在于：1)构造特征单元，描述了不同特征之间的关系，简化了几何复杂性，更适于大规模配准。2)引入虚特征解决了因遮挡等因素影响相邻视点重叠区缺少有效特征的配准问题。3)采用最小生成树原则建立了全局配准模型图，配准精度可靠。4)本发明提供的方法具有普遍适用性，既适用于现代结构的场景又适用于复杂的大型文物和古建。

附图说明

图1为本发明对场景被测物进行激光扫描的示意图；

图2(a)为垂直扫描面内激光采样点分布情况；

图2(b)为垂直扫描面内激光采样点按水平投影值的统计情况；

图3(a)和3(b)为主要建筑结构提取的实例；

图4(a)、4(b)、4(c)、4(d)和4(e)为在遮挡严重情况下两视点三维数据配准实例；

图5为全局配准模型生成过程；

图6(a)、6(b)为全局配准实例；

具体实施方式

以下结合附图及具体的实施例对本发明的技术方案作进一步描述：

本发明要求输入的是中远距离的三维激光测距仪获取的激光扫描数据。三维激光测距仪的作用距离一般大于200m，沿激光束方向测量精度是5cm/200m。本发明的具体实施步骤如下：

实施例1

1)数据获取

激光测距仪放在被测场景前20-50米，调节测距仪使Z轴垂直地面，X和Y轴平行于地面，对被测场景进行逐行逐列扫描(如图1)，获得对应视点的三维数据。另外，相邻视点扫描数据保持10％～20％的重叠，以便使各视点三维数据顺利配准。

2)提取主要建筑结构

输入每一视点对应的激光扫描数据，计算各垂直扫描面各点与激光测距仪的水平距离，并根据水平距离，利用哈夫变换检测垂直线段【如图2(a)和图2(b)】，以其所含的点数计算垂直线的长度，选择最长垂直线段对应的水平距离为基准d_mi，基准之前δ₁米和基准之后δ₂米范围内的数据为主要建筑结构保留。图2(a)某一垂直扫描面激光采样点的分布图，由图可以看出，保留区段内为主要建筑结构，其保留区段之前和之后的景物和散乱点删去；图2(b)为图2(a)中该垂直扫描面激光采样点水平距离统计分布示意图，由图可以看出，主要建筑墙面上的点在同一垂直线上，即水平距离相等。

本实施例δ₁＝1.0m，δ₂＝0.5m可提取主要建筑结构，图3(a)和3(b)为主要建筑结构提取实例。图3(a)为对被测物进行逐行逐列激光扫描的全景实测图，其中A为主要建筑结构(墙面)，B为地面，C为主要建筑结构的遮挡物；图3(b)为使用上述方法提取主要建筑结构的结果图，图中A为提取出的主要建筑结构(墙面)。

3)特征提取

a)计算扫描线上各点的一次微分

对于每条扫描线上采样点用(1)式进行局部二次拟合，并以该点沿扫描线方向上下k个邻近点为拟合区间[-k，k]，K为3-5；本实施例k＝3(当然K也可为4或5；可根据扫描点密度而定，至于究竟应该取多少，属本技术领域技术人员应知应会)；

式中：i-扫描线上采样点排序，i＝0，1，…，n；

第i点到激光测距仪的计算水平距离；

a₀，a₁，a₂-拟合方程(1)式的拟合系数

传统的方法根据(2)式的拟合误差函数E(i)，利用最小二乘法计算拟合系数a₀，a₁，a₂。

$E\left(i\right)=\mathrm{\&Sigma;}{\{d\left(i\right)-\hat{d}\left(i\right)\}}^2---\left(2\right)$

式中：d(i)-第i点到激光测距仪的实际水平距离；

但最小二乘法对不连续点和局外点不敏感，故采用自适应加权曲线拟合方法求取拟合系数a₀，a₁，a₂，式(3)为加权拟合误差函数。

$E^{\&prime;}\left(i\right)=\mathrm{\&Sigma;w}\left(i\right)\{d\left(i\right)-\hat{d}\left(i\right){\}}^2---\left(3\right)$

式中：w(i)-扫描线上第i点的权值，代表点的可靠性。采用迭代权值修正法，逐步消除局部点的影响，使拟合仅在可靠点内进行。具体方法如下：

假设已知第m次迭代中拟合系数为a₀ ^m，a₁ ^m，a₂ ^m，求第k+1次迭代中拟合系数a₀ ^m+1，a₁ ^m+1，a₂ ^m+1。首先计算第m次迭代的拟合残余量和第m次迭代的权值

$e^k\left(i\right)={(d\left(i\right)-{\hat{d}}^m\left(i\right))}^2---\left(4\right)$

$w^{m+1}\left(i\right)=\left\{\begin{array}{cc}1& e^m\left(i\right)\&le;b_s\\ \frac{b_s}{e^m\left(i\right)}& e^m\left(i\right)>b_s\end{array}\right.---\left(5\right)$

式中：b_s-误差允许阈值

初始权值w⁰(i)＝1，利用(3)、(4)、(5)式不断修正拟合系数和权值，至拟合曲线满足误差要求且包含最多点数为止。标记参与拟合的点。

然后，利用(1)式计算扫描线上有标记各点的一次微分。

b)检测边缘点

对于扫描线上任意一个有标记的点p，及其上下有标记的邻域点(p₁和p₂)，检测邻域点(p₁和p₂)一次微分的符号，若异号则p点为法向不连续，并标记为0。

计算p与p₁和p₂的距离，如果dis＝max(dis(p，p₁)，dis(p，p₂)＞ε(ε为阈值)，则p点为深度不连续，标记为l。

c)特征点线提取

首先采用链码跟踪方法将边缘点区分成不同的边缘链，并将每一边缘链上的点按顺时针排序，然后通过分裂合并的方法将每一边缘链拟合成一组首尾相接的特征线段，具体方法如下：

·初始化

将边缘链分为等长的若干段L₀，L₁，…，L_n

·分裂

利用最小二乘法对当前每一段L_i拟合成直线。选择离该拟合直线最远且最小残差E_i大于允许阈值T的一点为基准，将本段边缘链断开为两段，并不断分裂至满足拟合误差为止。依次处理其它各段。

·合并

选择拟合残差小于阈值T的相邻两段合并成一段为特征线段，特征线段的端点为特征点。

4)计算虚特征与构造特征单元

选择实际不相交的特征线段延长相交形成虚特征线段和虚特征点加入上述提取的特征线段和特征点中。将特征点连同其相邻的两特征线段嵌套成特征单元PU{P|L₁，L₂：<n₁，n₂，n₃，p>}。n₁和n₂为特征线段L₁和L₂的方向矢量，n₃＝n₁×n₂为前两矢量合成的第三矢量，×为矢量叉乘p为特征点。

5)配准

a)粗配准：给定一对分别属于两视点S₁和S₂的特征单元PU{P|L₁，L₂：<n₁，n₂，n₃，p>}∈S₁和PU′{P′|L₁′，L₂′：<n₁′，n₂′，n₃′，p′>}∈S₂，那么旋转矩阵会由n_i和n_i′(i＝1，2，3)得到，即使下式误差函数最小，并采用四元素法求出旋转矩阵。

$\mathrm{Err}=\underset{i=1}{\overset{3}{\mathrm{\&Sigma;}}}{\text{||}{n}_i^{\&prime;}\text{-R}{n}_i\text{||}}^2---\left(7\right)$

设p和p′在各自局部坐标系的坐标为(x，y，z)与(x′，y′，z′)，计算向量。

$t=R\left[\begin{array}{c}x\\ y\\ z\end{array}\right]-\left[\begin{array}{c}{x}^{\&prime;}\\ y^{\&prime;}\\ z^{\&prime;}\end{array}\right]---\left(8\right)$

每对对应特征单元得到一个对应的配准转换估计，选择匹配程度最大的为最佳配准转换估计。特征单元匹配的约束条件如下：

·位置约束

p′是视点S₂中与视点S₁中p点最近的特征点，且p和p′间的距离小于阈值Δ₁。

·方向约束

L_i与L_i′两线段所在直线间的距离和方向夹角分别小于Δ₂和Δ₃。

b)精配准：按粗配准获得的转换估计，将两视点的三维数据拼合在同一坐标系，确定出重叠区域。将两视点对应的重叠区分别等分成10块，计算每块内对应点的平均差异。选择差异最大的3块作为适配区。在适配区内利用ICP算法进行精配准为增加对应点搜索的鲁棒性，利用激光采样点的反射率和三维信息，以马氏距离衡量两点p∈S1和p′∈S2的差异，即

$d^2(p,p^{\&prime;})=\underset{p^{\&prime;}\&Element;S_2}{\min }{(p-p^{\&prime;})}^T{C}_{p^{\&prime;}}^{-1}(p-p^{\&prime;})---\left(9\right)$

式中：d²(p，p′)-p和p′点的马氏距离；

C_p′-p′点及其k×k个邻域点的协方差，本实施例k＝3。

图4(a)、4(b)、4(c)、4(d)和4(e)为招待所两视点三维数据配准结果；图4(a)为在S1视点获取的激光采样点云；图4(b)为在S2视点获取的激光采样点云；图4(c)和图4(d)为利用上述方法的粗配准结果，其中图4(c)没采用虚特征，图4(d)采用了虚特征；图4(e)为利用上述方法进行精配准的结果。

实施例2

为使用本发明的方法对学生餐厅进行大规模三维场景8个视点的激光扫描数据自动配准，其步骤1)-5)同实施例1，所不同的是：由于视点数目超过2个，需要进行全局配准：

首先所有(有重叠区)视点进行两两配准，利用最小生成树原则建立全局配准模型，然后根据全局配准模型选择一个视点S_a固定，其它视点的激光扫描数据转换在S_a所在的坐标系内。

具体地讲，每一对有重叠区的两视点都可计算出转换矩阵T_i＝[R_i，t_i]和局部匹配程度g(T_i)(以得到配准的点数计)。以各视点为结点，以g(T_i)为结点的连接权重，则会生成一个包含所有视点在内的无向带权连通图G＝<V，E>，V为结点，E为边(权重)，且e_i＝g(T_i)∈E。根据局部匹配程度g(T_i)，在无向连通带权图G＝<V，E>中，采用逐步短接法可求最小生成树，建立全局配准模型。

图5为其全局配准模型生成过程示意图，图中的结点1、2、3、4、5和6表示了6个视点；第1步，计算所有视点两两配准的局部匹配程度；第2步为将匹配程度最大化的两视点进行短接；第3-5步为将匹配程度次之的两视点进行短接，直至所有视点短接完成；第6步为将分枝最多的视点4固定为Sa。

图6(a)为本实施例《学生餐厅》8个视点全局配准结果图；图6(b)为图6(a)中窗户部位的局部放大图。

Claims (2)

1一种大规模三维场景多视点激光扫描数据自动配准方法，包括如下步骤：

1)获取数据

在被测场景前放置激光测距仪，调节激光测距仪使其Z轴垂直地面，X和Y轴平行于地面，对被测场景进行逐行逐列扫描，获得每一个视点的三维数据；相邻视点扫描数据保持10％～20％的重叠；

2)提取主要建筑结构

输入每一视点对应的激光三维扫描数据，计算各垂直扫描线上采样点与激光测距仪的水平距离，并根据水平距离，利用哈夫变换检测垂直线段；以其所含的点数计算垂直线的长度，并选择最长垂直线段对应的水平距离为基准，基准之前δ₁米和基准之后δ₂米范围内的数据点为主要建筑结构；所述δ₁为1米，δ₂为0.5米；

3)对保留的主要建筑结构的数据点进行特征提取

a)扫描线上各点微分

对于每条扫描线上的采样点用(1)式进行局部二次拟合，并以该点沿扫描线方向上下k个邻近点为拟合区间，其拟合区间[-K，K]为[-3，3]；

$\hat{d}\left(i\right)=a_0+a_{1}i+a_2{i}^2---\left(1\right)$

式中：i-扫描线上采样点排序，i＝0，1，…，n；

-第i点到激光测距仪的计算水平距离；

      a₀，a₁，a₂-拟合方程(1)式的拟合系数；

采用自适应加权拟合方法，求取拟合系数a₀，a₁，a₂，且加权拟合误差函数为

$E^{\&prime;}\left(i\right)=\mathrm{\&Sigma;w}\left(i\right){\{d\left(i\right)-\hat{d}\left(i\right)\}}^2---\left(2\right)$

式中：d(i)-第i点到激光测距仪的实际水平距离；

      w(i)-扫描线上第i点的权值，代表点的可靠性；采用迭代权值修正法，逐步消除局部点的影响，使拟合仅在可靠点内进行；具体方法如下：假设已知第m次迭代中系数为a₀ ^m，a₁ ^m，a₂ ^m，求第k+1次迭代中系数为a₀ ^m+1，a₁ ^m+1，a₂ ^m+1；首先计算第m次迭代的拟合残余量和第m次迭代的权值；

$e^k\left(i\right)={(d\left(i\right)-{\hat{d}}^m\left(i\right))}^2---\left(3\right)$

$w^{m+1}\left(i\right)=\left\{\begin{array}{cc}1& e^m\left(i\right)\&le;b_s\\ \frac{b_s}{e^m\left(i\right)}& e^m\left(i\right)>b_s\end{array}\right.---\left(4\right)$

式中：b_s-误差允许阈值；

初始权值w⁰(i)＝1，用(3)、(4)、(5)式不断修正拟合系数和权值，直至拟合曲线满足误差要求且包含最多点数为止；标记参与拟合的点；

然后，利用(1)式计算扫描线上有标记各点的一次微分；

$\frac{\&PartialD;\hat{d}\left(i\right)}{\&PartialD;i}=a_1+2a_{2}i---\left(5\right)$

b)检测边缘点

对于扫描线上有标记的各点，考察其沿扫描线方向上下有标记两点的一次微分符号，选择异号的点为法向不连续点，并标记为0；

对于扫描线上有标记的各点，计算其与沿扫描线方向上下有标记两点距离，考察其中最大距离，选择最大距离大于阈值ε为深度不连续点，并标记为1；

c)特征点线提取

首先采用链码跟踪方法将边缘点区分成不同的边缘链，并将每一边缘链上的点按顺时针排序，然后通过分裂合并的方法将每一边缘链拟合成一组首尾相接的特征线段，具体方法如下：

·初始化

将边缘链分为等长的若干段L₀，L₁，…，L_n

·分裂

利用最小二乘法对当前每一段L_i拟合成直线选择离该拟合直线最远且最小残差E_i大于允许阈值T的一点为基准，将本段边缘链断开为两段，并不断分裂至满足拟合误差为止；依次处理其它各段；

合并

选择拟合残差小于阈值T的相邻两段合并成一段为特征线段，特征线段的端点为特征点；

4)计算虚特征与构造特征单元

选择实际不相交的特征线段延长相交形成虚特征线段和虚特征点加入上述提取的特征线段和特征点中；将特征点连同其相邻的两特征线段嵌套成特征单元PU{P|L₁，L₂：<n₁，n₂，n₃，p>}；n₁和n₂为特征线段L₁和L₂的方向矢量，n₃＝n₁×n₂为前两矢量合成的第三矢量，×为矢量叉乘p为特征点；

5)配准

a)粗配准

给定一对分别属于两视点S₁和S₂的特征单元PU{P|L₁，L₂：<n₁，n₂，n₃，p>}∈S₁和PU′{P′|L₁′，L₂′：<n₁′，n₂′，n₃′，p′>}∈S₂，那么旋转矩阵会由n_i和n_i′(i＝1，2，3)得到，即使下式误差函数最小，并采用四元素法求出旋转矩阵；

$\mathrm{Err}=\underset{i=1}{\overset{3}{\mathrm{\&Sigma;}}}{\left|\right|n_i^{\&prime;}-{\mathrm{Rn}}_i\left|\right|}^2---\left(6\right)$

设p和p′在各自局部坐标系的坐标为(x，y，z)与(x′，y′，z′)，计算向量；

$t=R\left[\begin{array}{c}x\\ y\\ z\end{array}\right]-\left[\begin{array}{c}{x}^{\&prime;}\\ y^{\&prime;}\\ z^{\&prime;}\end{array}\right]---\left(7\right)$

每对对应特征单元得到一个对应的配准转换估计，选择匹配程度最大的为最佳配准转换估计。特征单元匹配的约束条件如下：

·位置约束

p′是视点S₂中与视点S₁p点最近的特征点，且p和p′间的距离小于阈值Δ_l；

·方向约束

L_i与L_i′两线段所在直线间的距离和方向夹角分别小于Δ₂和Δ₃；

b)精配准

按粗配准获得的转换估计，将两视点的三维数据拼合在同一坐标系，确定出重叠区域；将两视点对应的重叠区分别等分成N块，计算每块内对应点的平均差异；选择差异最大的M块作为适配区，在适配区内利用ICP算法进行精配准为增加对应点搜索的鲁棒性，利用激光采样点的反射率和三维信息，以马氏距离衡量两点p∈S1和p′∈S2的差异，即

$d^2(p,p^{\&prime;})=\underset{p^{\&prime;}\&Element;S_2}{\min }{(p-p^{\&prime;})}^T{C}_{p^{\&prime;}}^{-1}(p-p^{\&prime;})---\left(8\right)$

式中：d²(p，p′)-p和p′点的马氏距离；Cp′-p′点及其k×k个邻域点的协方差，k＝3。

2、按权利要求1所述的大规模三维场景多视点激光扫描数据自动配准方法，其特征在于，对于视点超过2个以上的多视点三维场景，所述步骤还包括全局配准：

首先所有有重叠区的视点两两进行配准，利用最小生成树原则建立全局配准模型，然后根据全局配准模型选择一个视点S_a固定，其它视点的激光扫描数据转换在S_a所在的坐标系内，其中建立全局配准模型的具体方法如下：

设第i对有重叠区的两视点计算出的转换矩阵为T_i＝[R_i，t_i]和对应的局部匹配程度为g(T_i)；以视点为结点，以g(T_i)为结点连接权重，生成一个包含所有视点在内的无向带权连通图G＝<V，E>，V为结点，E为边，且e_i＝g(T_i)∈E，并采用逐步短接法求出最小生成树，建立全局配准模型。

Images