改进迭代最近点算法支持下的LiDAR点云与正射影像粗-精配准方法

徐卫红，彭家真，高红旗

(浙江华东测绘与工程安全技术有限公司，310014，杭州)

0 引言

当前遥感点云数据获取最广泛的手段包括激光扫描和摄影测量，这2种方式都能获得高精度的点云数据，生成的点云已广泛用于各种生产应用[1]。LiDAR使用主动式的探测方式，基于极坐标几何点位原理快速获取空间点云信息，存在成本高、点云密度较低、缺乏地物光谱信息等缺点[2]。无人机摄影测量则是一种被动测量目标物的手段，其通过影像密集匹配可以得到比LiDAR点云密度更高、地物几何特征更加明显的空间点云。影像中丰富的光谱和纹理特征可以作为LiDAR点云数据的补充，有效地解决了单一数据源信息不足问题。通过LiDAR点云与影像融合，能产生更丰富的地物空间信息和语义信息，而2种点云数据融合的前提条件是进行点云配准。因此，对2种不同的点云数据进行配准是一个至关重要的环节。

点云配准是指通过某种配准算法计算得到2组不同点云的刚体变换参数，然后再利用变换参数将2组点云统一到同一坐标系下。国内外研究人员提出了许多点云配准算法，其中ICP算法具有简单、复杂度低等优点，是点云匹配最为常用的算法之一。然而，ICP算法存在计算效率低、鲁棒性差、无初值时收敛慢等问题。因此，针对ICP算法效率不高的问题，许多研究人员提出了大量的改进算法[3]。Li等通过动态调整刚体变换参数的因子减少ICP算法迭代次数，实验结果表明改进的算法提高了对应点搜索效率，点云配准时间相比于传统方法要更短，且精度更高[4]。孙炜等针对ICP配准时间长问题，提出了融合轮廓特征快速点云配准方法[5]。李运川等利用重心重合法提高点云重叠度，利用RANSAC算法得到点云配准较好的变换参数初始值，最后改进后的ICP算法实现点云高精度配准，且配准效率高于传统的ICP算法[6]。这些改进的算法主要是通过减少在2组点云中对应点的搜索时间来提高配准效率。

LiDAR点云与影像匹配点云的配准方法主要分为3种。第1种是将影像匹配得到的点云与LiDAR点云进行配准。刘燕萍等基于影像密集匹配得到三维点云，通过粗配准建筑物屋顶的特征角点，再根据提取屋檐特征线以实现点云的精配准[7]。李彩林等通过影像生成三维稀疏点云，以影像三维点云和激光点云拟合的最邻近曲面为约束条件，结合共线条件方程实现影像三维点云和激光点云的高精度配准[8]。Stamos等利用运动恢复结构(Structure from Motion，SfM)算法恢复物方稀疏点云，最后将LiDAR点云和摄影测量点云进行配准[9]。第2种是将LiDAR点云转化为二维图像，再利用点云配准方法进行匹配。赵中阳等先利用LiDAR点云生成的深度图像，基于航空影像和深度图像的面特征采用SIFT算子提取点特征，最后通过面特征到点特征的配准策略实现LiDAR点云与航空影像的配准[10]。徐景中等通过LiDAR点云生成深度图像，提取深度图像与对应影像的结构特征，采用粗匹配-精匹配的循环迭代匹配策略，来实现LiDAR点云与航空影像的自动配准[11]。第3种是建立影像与LiDAR点云之间的直接配准关系。Zhu等提出一种稳健有效的由粗到细的配准方法，利用空间约束和Gabor结构特征进行配准，结果表明该配准方法对几何失真和辐射变化具有较高的鲁棒性[12]。Habib等利用直线特征建立遥感影像与激光点云之间的直接映射关系，这种方法对于变化检测和系统校准应用很有用[13]。王育坚等提出一种基于保局(Principal Component Analysis，PCA)的三维点云配准方法，利用特征寻找点云之间的对应关系来实现配准[14]。然而，大多LiDAR点云与影像配准方法存在数据冗余、自动匹配精度不高或效率低下等问题。

本文基于ICP算法，采用粗-精配准匹配策略来实现LiDAR点云与正射影像的配准，使得点云与影像匹配兼顾精度与效率。首先利用体素滤波器对点云数据进行预处理，再利用PCA算法压缩数据量和获得主分量进行点云粗配准来获得点云刚体变换的初始值，最后采用ICP算法对点云进行高效地精配准。

1 LiDAR点云与正射影像配准

1.1 点云预处理

本文通过无人机采集得到影像，利用Pix4D软件经过密集匹配得到点云数据，研究区域的正射影像如图1所示。同时在激光雷达设备上输出得到同一研究区域的激光点云，激光点云高度渲染图见图2。由正射影像和LiDAR渲染图可以看出，正射影像含有丰富的纹理信息，而LiDAR点云具有高精度的空间几何信息。

图1 正射影像

图2 LiDAR点云渲染图

在进行点云配准之前，为了去除冗余点云数据，减少对应点对搜索时间，本文利用体素滤波方法对点云数据进行预处理。体素滤波在点云下采样过程中的优势在于去除点云噪声点和离群点的同时，保留了点云本身的几何结构。体素滤波的原理是根据输入的点云计算得到包围整个点云的立方体，将立方体根据设定的分辨率分割成多个小立方体，最后将小立方体内点云的质心坐标近似为小立方体内的点云。见图3，体素滤波的计算过程如下：1)计算三维点云数据坐标轴上的最大值和最小值；2)设计体素小栅格的边长；3)根据得到的坐标轴最值计算点云最小包围盒的边长；4)计算体素网格的尺寸；5)计算每个点云在体素小栅格内的索引；6)从大到小排列索引里的元素，并以各体素重心代替其栅格内的所有点，若重心不存在，则用距重心最近的数据点代替栅格内所有点。

1.2 迭代最近点算法

ICP算法是被提出用于三维点云配准[15]，该算法通过迭代搜索2个点云距离最近的对应点，以各对应点对欧氏距离平方和最小为目标函数，再利用最小二乘计算最优的2组点云配准刚体变换参数。

假定有待配准点云集P和基准点云集M，M和P分别包含Nm和Np个点云。以m=[mR|mT]来表示2个点云集间的刚体变换参数。利用矩阵矢量m在基准点云集中搜索对应于待匹配点云集中最近的点，在这个过程当中会搜寻到多个点。为了得到最优变换参数，欧氏距离平方和最小的目标函数表达式为：

图3 体素滤波过程

(1)

式中mi和pi表示点云集M和P中的各个点云。

计算待匹配点云集和基准点云集的重心：

(2)

(3)

计算2个点云集间的互协方差阵，并利用四元数法计算得到一个4×4的对称矩阵：

(4)

(5)

再利用公式(6)计算平移向量T：

mT=μm-R(mR)μp

(6)

对点集P中的每一个点pi应用变换函数E得到新的点云集Pnew，并在基准点云集M中搜索Pnew的对应点得到Mnew。

(7)

经过以上过程的多次迭代计算，直到Pnew和Mnew的误差小于给定的阈值或者达到设定的最大迭代次数时终止迭代，输出R和T。

1.3 基于PCA的粗配准

PCA是一种经典的高维数据降维方法[16]。PCA算法在点云配准过程中主要用于给精配准阶段提供刚体变换参数的初始值，提高点云配准效率。其算法原理是通过计算点云的主方向，然后根据主方向计算得到点云配准的初始值并对初始值进行校正。首先计算2个点云集的协方差阵，通过协方差阵求出点云的主轴方向，然后利用主轴方向计算点云集的旋转矩阵，并通过2个点云集的中心偏移量求出点云集的平移向量[17]。基于PCA方法进行点云粗配准的步骤如下：

计算2组点云的中心：

(8)

(9)

计算2个点云集间的协方差阵：

(10)

对CP和CM进行奇异值分解，将协方差阵转化为由2组点云的主方向UP和UM构成的矩阵：

(11)

通过点云主方向计算旋转矩阵和平移向量的初始值R0和T0：

(12)

最后，针对初始值R0和T0进行校正，由于各个主轴方向存在反向的可能性，求出的初始R0和T0并不一定能用于点云精配准。因此，还需要对R0和T0的初始值进行校正，校正过程如下。

利用得到的R0和T0和公式(14)对点云进行旋转：

pnew=P*R0+T0

(13)

在基准点云集M中搜索Pnew的最近距离对应点，计算得到对应点间的距离：

H=pnew-mi

(14)

计算初始配准误差：

erri=H*HT

(15)

(16)

式中erri为第i个对应点间的误差，error为2组点云间的平均均方误差。

将UP中3个列向量(待配准点云的主方向)与目标点云各个轴的方向进行对比，若方向相反，则把相反方向所对应的UP的列向量取反，如公式(17)：

UP(i)=-UP(i)

(17)

1.4 ICP精配准

得到预处理后的点云数后，利用迭代最近点算法对预处理后的点云进行匹配，见图4，配准过程如下。

1)体素滤波精简点云。利用体素滤波对待匹配点云P和基准点云M进行采样，计算获得精简后的待匹配点云P′和基准点云M′。

2)ICP迭代搜索最近对应点。对于待匹配点云中的每一个点云pi，在基准点云集中搜索其对应点mi。

4)实施点云变换。对待匹配点云集中的点云pi进行变换得到新的点云集P2。

5)计算变换后的点对距离。

6)迭代终止判断。若通过变换参数计算得到的Pnew和Mnew的距离误差小于给定的阈值或超过最大迭代次数时，终止迭代运算。

2 实验结果与分析

为验证本文方法的有效性，本文实验选取具有代表性的四川某一地区的LiDAR点云与无人机数字正射影像进行配准实验。该测区包含丰富的地物信息如：耕地、林地、建筑物、道路、水塘等地表信息，覆盖范围约为1 km2。

LiDAR点云数据获取采用SKY-LARK轻小型无人机LiDAR测量系统，该系统采用椭圆轨迹扫描方式，最高点频能获得高达60万点/s，最大测量距离可达500 m。采用大疆无人机航拍影像，飞行2个架次(每架次飞行时间约为18 min)，规划航线共8条，相对航高设为150 m，航向与旁向重叠度均设为80%，地面分辨率约为2 cm，整个研究区的航拍影像为870幅，在Pix4D完成立体像对匹配、相对定向、空中三角测量等处理，并生成数字正射影像图。本文方法LiDAR点云与正射影像配准后，点云三维渲染图见图5，纹理映射后结果图见图6。

图5 LiDAR三维渲染图

图6 附带纹理的三维图

此外，为对比分析本文方法针对原始的ICP算法的改进效果，在同等软硬件环境下，采用相同的数据对ICP算法进行LiDAR点云与无人机正射影像配准，并采集了10个地面控制点(Ground Control Point, GCP)进行精度检验见图7，表1为本文方法与原始ICP算法在迭代运算次数、平均误差、效率方面的对比。

图7 LiDAR三维渲染图

表1 ICP算法与本文方法对比

由对比结果可知，本文提出的方法配准结果明显优于ICP算法，能实现LiDAR点云与正射影像高精度配准。在配准效率方面，采取基于粗-精的配准策略，能去除点云噪声，并加快配准时间，本文采用体素滤波方法对原始点云数据进行预处理，并且利用PCA算法进行粗配准，以获取精配准过程中所需的刚体变换参数的初始值，相比于未经过预处理和粗配准的点云配准，处理后的点云配准速度有明显的加快，这有利于进行点云实时配准。在精度方面，本文通过PCA粗配准和ICP精配准的点云配准策略，实现了点云与正射影像高精度配准。相比原始ICP算法精度和效率分别提高了约1.8倍和3倍。

3 结论

针对传统的ICP算法难以满足LiDAR点云与正射影像高精度、高效的配准的需求，本文结合PCA粗配准与ICP精配准策略来实现LiDAR点云与无人机正射影像配准。为去除LiDAR点云数据中的离群点和噪声点，本文采用体素滤波方法对点云数据进行预处理以消除冗余数据，提高点云与影像配准效率，然后，利用PCA算法计算点云配准刚体变换参数初值，并校正后作为ICP算法的输入值，从而提高ICP算法配准的准确性和收敛效率。实验结果表明，本文的配准策略在点云与正射影像配准方面能获得了较好的配准结果，能满足包含复杂地表地物分布的LiDAR点云与影像多源数据配准需求。