刘德儿,刘 鹏,肖 健
(1.江西理工大学建筑与测绘工程学院,江西 赣州 341000;2.江西理工大学外语外贸学院,江西 赣州 341000)
随着城市化建设的不断推进,构建智慧城市已经成为提升城市竞争力的主要方法,因此对建筑物的精细化建模也提出更高的要求[1]。无人机倾斜摄影能够获取建筑物俯视角度及立面丰富的表面纹理信息,而后使用高效的重建算法将多视角的影像解算出空间结构信息,从而得到三维模型[2]。可如今城市建筑的密集化对建筑立面造成的遮档和光线的削弱,使得建筑立面信息的采集不够完整,因此模型精度也将受到影响,容易出现空洞和拉花。因此,使用跨源数据构建完整的精细化模型将是未来的发展趋势[3-4]。
融合配准[5]是三维重建流程中的核心步骤,如今的主流方法为粗配准和精确配准算法相结合[6-7]。如刘剑等[8]人在计算快速点特征直方图(Fast Point Feature Histograms,FPFH)的局部特征信息描述符后,与三角网相结合,使特征提取问题简洁化,因此能够提高特征点的粗匹配的效率;田青华等[9]人提出了将方向直方图签名描述子(Signature of Histograms of OrienTation,SHOT)和空间位置相结合,再采用最小方差法寻找正确的对应点对估算转换关系;秦楠楠等[10]人利用主成分分析(Principal Component Analysis,PCA)得到点云数据的特征向量和质心,从而对点云数据进行空间变换,得到初始匹配结果;刘江[11]使用数据中心重合法的粗配准方法得到初始位置。精确配准算法应用最多的是迭代最近点(Iterative Closest Point,ICP)算法,由Besl[12]等人提出,是一个拼接效果极佳的算法。针对ICP算法的优缺点,衍生出许多改进算法。文献[13]使用最优节点优先算法来提升K-D树对邻近点的搜索速度,降低计算复杂度,从而提高迭代计算的效率;Chetverikov[14]等人提出TrICP(TrimmedIterative Closest Point,TrICP)算法,计算速度快,对于重叠区域低的数据具有较高的鲁棒性。以上配准方法大多是针对数据完整、重叠率高的条件下进行的,但是对于缺乏对称性和低重叠率的数据,则不能保证其适用性。
为获取完整的三维模型,本文将融合倾斜影像和激光点云数据用于建筑物的无缝三维重建。其中三维激光点云和航拍影像属于异维异质数据,故先将倾斜影像几何重构生成影像点云再与激光点云融合建模。本文提出使用PCA方法确定跨源点云的姿态方向对齐,然后利用几何中心重合法完成数据矫正,最后运用TrICP算法进行精确配准,获得高精度的转换关系,最终得到数据完整的无缝三维模型。
本文实验数据为江西理工大学校本部图书馆。通过图1可以看出受高大树木及其他建筑的遮挡,使得数据采集有着不小的难度。使用独立的数据采集方法,会不可避免地出现数据缺失,所以此建筑的数据可以很好地验证本文的三维重建方法。
图1 建筑物实景
激光点云数据由Riegl公司型号VZ1000的三维激光扫描仪采集获得,由于高大数木对建筑立面遮挡严重且建筑物形状复杂,故扫描站次较多。
将建筑点云源数据进行相应的预处理[15-16],首先对数据的离群点去噪,而后拼接配准,并手动分割出图书馆点云数据。运用文献[16]的方法去除噪声,该方法主要运用二维的Gamma分布,将邻域均值以及邻域距离的变化斜率作为两个限制条件剔除噪声点。通过手动选取同名点粗配准和ICP精确配准。手动配准时根据两个站点扫描场景的重叠区域挑选出至少4对同名点,随之计算出变换矩阵,将第2站点数据转换到第1站点坐标系中,然后使用ICP算法实现两个站点的精确拼接,以此逐步拼接。为避免累积误差,依次将第3,4,…,12站的数据转换到第1站点云数据的坐标系中,即所有站点的数据都处于统一的全局坐标系中。完整的激光点云数据如图2所示。
图2 激光点云预处理拼接
校园的倾斜影像通过搭载五镜头的无人机航拍采集获取,采集前应预先制订好航行路线,设置相应的重叠度参数,选定空旷的起降区域。为了更好的验证本文三维重建方法,从拍摄的海量影像中挑选出含有此建筑的影像179张,经过现场及后期检查,获取的影像清晰容易辨别,亮暗适中,可用于后期航片解算。
倾斜摄影有着4个倾斜角度镜头,因此影像的同名点数量增多,虽然计算复杂度增加,但计算精度得到提升。CotextCapture软件具备高效的匹配算法以及区域切块处理的技术,可以轻松地完成高计算复杂度的任务。将挑选出的影像加载进CotextCapture软件的新建工程中,剔除不合格影像和匀光匀色处理。为了保证影像点云与激光点云的尺度一致,在进行空中三角测量计算之前加入尺度约束,从而降低尺度问题引起的配准难度。为方便后续三维重建处理,对影像点云手动分割修整去除目标建筑物以外的物体。高密度点云如图3所示,可以看出建筑立面的数据缺失严重。
图3 影像处理结果
3.1.1 PCA坐标系求解
PCA作为一种经常使用的数据分析方法,是通过将高维的原始数据降成各维度线性无关的低维数据,从而有利于提取数据的特征分量,将方差较小的数据特征忽略掉,达到降维处理的目的,突出数据最主要的方面,甚至可以实现去噪的效果。本文对于点云数据进行降维分析只需找到相互垂直的三个主要特征分量,即对方差贡献最大的三个向量,用以对应点云的分布方向以及表达点的空间位置。
(1)
(2)
(3)
对Cov(X,Y,Z)进行奇异值分解,得到特征值λi对应的特征向量wi,并对特征值λi进行降序排列得到对角矩阵∑=diag[λi,λ2,λ3],和两组正交矩阵U(3×3)、V(3×3)。
Cov(X,Y,Z)=U∑VT
(4)
其中,正交矩阵U由∑对应的特征向量组成,表示点云PCA坐标系三个主轴的方向向量。U中三个向量分别表示点云三个主元方向,并且三个向量相互垂直,如图4所示,其中第一个向量表示点在空间分布趋势最大的方向。在提取出主成分后,将激光点云和影像点云转换到同一个坐标系下,其中转换关系就是协方差矩阵的对应关系,对两个点云经过全局进行转换,从而达到点云粗配准的效果。
图4 跨源点云的主元方向
3.1.2 矫正PCA坐标系
通过航空影像获取的影像点云在建筑体表面分布均匀,密度较高,且呈现单层次分布;而三维激光扫描仪工作时受到站点视角以及障碍物遮挡等情况的影响,使得激光点云在空间中的分布并不均匀,有的区域密度很高,有的区域点云分布稀少,以及建筑物顶部数据的缺失。所以通过这两种方式获取的点云数据具有不一样的特性,从而导致两个点云数据PCA配准的结果并不是很准确。其中坐标轴的方向大致平行,但是偶尔出现反向的情况,以及数据质心的偏差。
首先通过点云的包围盒来矫正PCA坐标系的原点。点云的外包盒一般有球形、方向包围盒和轴向包围盒AABB(Axis-Aligned Bounding Boxes,AABB)三种,其中AABB包围盒是平行于坐标轴,即每一个面都有与其相互垂直的坐标轴。由于PCA配准后两个点云的主方向大致确定,不再将其方向二次旋转,且AABB包围盒对方向很敏感,所以选择AABB包围盒进行二次匹配。
1)使用源点云A和目标点云B的包围盒顶点的极值坐标计算包围盒几何中心μa和μb,其极值为点云包围盒的顶点坐标的极大值AMax、BMax和极小值AMin、BMin。
(5)
(6)
2)得到几何中心μa和μb后,计算它们的差值,也就是源点云的二次平移向量T′。
(7)
3)利用新计算的二次平移向量T′对源点云进行平移。
(8)
式中,A′为源点云平移后的数据点集。
对于三个轴向的正反方向情形将引起23种对应关系,本文将设定一个距离阈值d,使用K-D树加速搜索对应点对,小于该阈值d的对应点将被标记,然后计算被标记的点云数量所占目标点云PM的重叠率γ=Ir/PM,比较这8种对应关系的重叠率γ大小,将重叠率最大的对应关系作为主轴的方向。
(9)
将两类点云的几何中心匹配和矫正坐标轴方向,可以减小两类点云的欧式距离,缩小粗配准的偏差,此方法也能够得到两类点云更大区域的重叠,较优的粗配准效果,为精配准提供更佳的位姿状态。
经典ICP配准算法主要思想是通过源点云向目标点云寻找欧氏距离最近的点,使每一个源点都有着对应目标点,解算出变换矩阵完成转换;然后进行第二次寻找最近点,再解算出新的变换矩阵完成转换;如此往复迭代,直到满足迭代终止条件。
TrICP算法则是使用重叠率Υ来自适应的匹配算法,将非匹配点去除以改进ICP算法。从算法结果来看,与原算法主要的区别在于可以通过不同的重叠率来影响配准结果;从算法原理来看,改变了MSE的统计方法,即算法的全过程都使用截断最小二乘(The Least Trimmed Squares,LTS)计算误差平方和。如式(10)所示,当h=n时,把所有的残差考虑在内,等价于原始算法。
(10)
qcl(i,R,T)=argmin||q-pi(R,T)||
(11)
(12)
(13)
(14)
式中,Sumdi代表Npo个重叠点对的残差平方和,本次迭代的MSE=Sumdi/Npo,随之判断是否满足终止条件,选择是否继续迭代计算;
(3)通过奇异值分解计算出变换关系(R,T),使得满足Npo点对的Sumdi最小;
(4)将源点云按照变换关系(R,T)向目标点云方向运动,重复步骤(1)。
为了验证本文算法对于跨源数据三维重建的有效性,将在Inteli5 CPU、Nvidia1060显卡以及16GB内存的PC设备上进行实验,并使用PCL(Point Cloud Library)库编程实验进行分析。使用经过了预处理的激光点云和倾斜影像数据,如图3和图4所示,并通过体素网格下采样方法对点云数据进行抽稀,最终用于验证本文算法的激光点云数据的数目为998545,影像点云数目为1428361。
跨源点云经过预处理后,首先进行各自的PCA分析,计算出两类点云各自质心以及三个最大的特征向量,进行点云旋转和平移,完成PCA变换。其中激光点云质心O1(9.9412,47.6741,6.2562),影像点云质心O2(-20.5584,-47.7046,9.5665)。通过PCA变换,如图4,质心O1和O2将作为新坐标系的原点,并完成几何中心的匹配。比较了另外两种常见算法,以及相应的精度比较。由于跨源点云的差异性,配准后也难以使两个点云完美的贴合在一起,且不一样的精度统计方法,得到的精度结果也大不一样。考虑到本文实验数据重叠率低的特点,即源点云不是目标点云的子集,且各自有较大面积的非对应区域。所以本文在精度统计时,将非对应点剔除,统计对应点的距离,然后使用计算出的均方根值(RMS)作为精度表达。
对实验数据进行多种常规的点云粗配准方法实验后,如图5所示,展示点云变换后的正视、侧视和俯视配准效果,(a)为建筑正视图;(b)为建筑侧视图;(c)为YOZ平面切片。其中正态分布变换(NDT)对于跨源点云的粗配准效果不佳,远达不到精配准的初始位姿要求。基于FPFH的采样一致性粗配算法对跨源点云具有一定的适用性,只是对于本实验数据需要设置较大的计算参数。本文的粗配方法为经过PCA变换后矫正点云中心和轴向,使得跨源点云较为紧密的拼接在一起,从图中观察可知效果明显优于另外两种方法,并由表1可知,初始位姿误差也远远低于另外两种方法,达到了0.32 m,将其作为后续精配准工作的初始位姿。
图5 粗配准实验结果对比
表1 初始配准精度比较
在激光点云得到初始位姿后,使用参数为60 %重叠度和最大30次迭代次数的TrICP算法对其进行精确配准,并与参数为最大30次迭代次数的经典ICP算法、GICP算法和FPFH+NDT算法[17]结果进行比较。如图6所示,在本文粗配准方法获得的初始位姿条件下,比较了ICP、GICP和TrICP精配结果的局部差异。可发现三种算法对俯视角的配准效果基本一致,但是多个立面的顶部区域有着细微的差异,可以看出在经典ICP算法结果中,激光点云整体上移,陷入了局部最优。
在表2的误差分析中,统计了各个立面之间的独立误差。其中文献[17]使用的FPFH+NDT算法对于本文实验效果较差,不适用于此类跨源点云的配准工作。对NDT算法分别进行粗配准和精配准的实验可以看出,NDT算法同样需要一个初值才能获得更好的运行效果。在三个立面的误差分析中,正立面和侧立面的误差,TrICP要明显小于ICP和GICP;背立面的误差,TrICP却大于ICP和GICP。在整体精度中,TrICP的精度相对于传统ICP,总体提高了53.5 %,达到了0.0811 m。
图6 精配准的细节对比
表2 精配准误差分析(包括不同立面的单独误差,总误差)
图7所示为最终跨源点云的配准结果。配准误差较大的部分点云主要集中在靠近地面区域,误差偏小的点大多分布在远离地面的建筑立面上,因此说明倾斜影像获取立面信息的准确度不一致。影像点云的数据缺失在建筑两侧的立面区域最为严重,图7中显示激光点云能够对其予以修补。总的来说,配准误差由建筑立面从上往下逐渐增大。考虑到激光雷达的高精度扫描,影像点云越接近于地面时,获取的立面信息偏差也增大,跨源点云的配准精度也随之降低。从实验结果看来,本文通过实验完成了异维异质的跨源数据配准,对于这类联合多源数据的三维重建问题具有一定的指导意义。
图7 精配准结果
本文提出一种跨源数据融合配准的三维重建方法。该方法通过PCA变换后的点云,并用中心匹配来获取一个较好的初始位姿状态,然后使用TrICP算法完成精配准工作。针对激光点云数据获取时受到视角的影响难以获取建筑顶部信息,以及倾斜影像对立面信息获取的不完整性,使用两种数据融合完成三维重建,能够弥补激光点云视角的限制和倾斜影像在建筑物立面的空洞、纹理扭曲等缺陷,提高建筑物的三维模型的信息完整性。如何消除建筑立面上小下大的误差分布特点,将是后续的主要研究工作。