王帅帅, 柏艳红, 王 银, 孙志毅
(太原科技大学电子信息工程学院, 太原 030024)
三维激光扫描技术广泛应用于医学研究[1]、文物保护[2]和三维重建[3]等领域.在三维数据采集过程中, 由于受激光扫描仪测量角度或外界障碍物遮挡等因素的影响, 被扫描物体的部分点云信息易于遗漏.为了获取目标物体表面完整的点云数据, 须采用点云配准[4-5]技术将多视角下的点云数据转换到同一坐标系下, 为后续的三维重建奠定基础.由于点云配准是三维重建的关键技术, 故快速且精确的点云配准技术具有重要的应用价值.目前,应用最广泛且最经典的点云配准方法是Besl等[6]提出的迭代最近点(iterative closest point, ICP)算法, 该方法以四元数运算模型为基础, 不断迭代获取使得两片点云间欧式距离最小的变换矩阵, 但对两片点云的初始位姿依赖性较强,在迭代时容易陷入局部极值, 导致配准结果较差.针对ICP算法的缺陷, 研究者大多通过有效的粗配准获得良好的初始配准位姿, 再采用ICP算法完成精配准[7].Lei等[8]针对重叠小、噪声大且位姿差异大的点云, 提出了一种基于均匀采样点的全局点云配准方法, 利用由距离误差表示的质量函数识别最优变换并实现点云的粗配准, 进而采用改进的ICP算法完成精配准, 显著提升了配准效率和精度; 荆路等[9]基于尺度不变特征变换(scale invariant feature transform, SIFT)关键点结合快速点特征直方图(fast point feature histogram, FPFH)描述子, 利用采样一致性初始配准 (sample consensus initial alignment, SAC-IA)算法完成初始配准, 该算法虽提高了配准精度, 但配准时间较长; Zaganidis等[10]采用点云的语义信息作为搜索对应关系的先验信息, 辅助正态分布变换(normal distribution transform, NDT)配准点云, 改善了配准精度和效率,但适用场景有限; 李宇翔等[11]基于内部形态描述子(intrinsic shape signatures, ISS)关键点[12]利用余弦相似度匹配对应点对获得良好的初始位姿, 然后采用点到平面的ICP算法进行精配准, 该算法对低重叠率的点云具有较好的鲁棒性, 但方向直方图(signature of histograms of orientations, SHOT)描述子匹配效率较低; Li等[13]提出了一种基于对应关系的灵活配准方法, 该方法配准精度高、速度快且鲁棒性较强, 但其配准精度在一定程度上依赖于关键点检测器的定位精度; Li等[14]采用高斯混合模型拟合点云模型, 配准精度虽得以提高, 但配准效率较低.本文拟提出一种结合内部形态描述子关键点及二进制方向直方图描述子(binary signature of histograms of orientations, BSHOT)的点云配准方法, 通过粗配准求解初始变换矩阵, 利用该矩阵将源点云转换至目标点云坐标系, 再采用ICP算法完成精配准.
本文配准过程由粗配准和精配准两部分组成: 1) 粗配准.由于BSHOT特征描述子结合了点云的密度和噪声等几何分布信息, 所以首先计算点云的分辨率, 利用ISS算法提取两片点云的关键点并在关键点的邻域构建局部坐标系, 统计邻域点和查询点法向量的夹角余弦值并计入直方图中; 然后将高维直方图特征转化为BSHOT描述子, 利用汉明距离匹配对应点对,再采用随机采样一致性(random sampling consistency, RSC)算法去除误匹配点对, 得到点云配准的初始变换矩阵; 2) 精配准.通过点对点的ICP算法估算出最优变换矩阵.具体配准流程如图1所示.
图1 本文配准算法流程Fig.1 Registration algorithm process of this paper
当传感器视角发生稍许变化时,点云中的关键点依然保持不变,其稳定性强、辨识度高且信息量丰富,关键点的数量远小于原始点云数量,故本文采用ISS算法提取关键点.一个内部形态描述子ISS由一个局部参考坐标系和一个编码三维形态特征的高分辨率特征向量组成.
ISS算法通过计算各个查询点与邻域范围内所有点的协方差矩阵的特征值, 并设定特征值间的比例阈值选取关键点.设点云P中共有m个点, 任意一点pi的原始坐标为(xi,yi,zi),i=1,2,…m, ISS关键点提取过程如下:
1) 对点云数据P中的任一查询点pi构建局部参考坐标系, 并将rf视为搜索半径;
2) 统计点云P中每个查询点pi在半径为rf的球形邻域内的所有点pj, 并计算其权值
(1)
3) 利用rf范围内的所有邻域点pj计算查询点pi的加权协方差矩阵
(2)
5) 设定两个小于1的阈值ε1,ε2,若满足条件:
(3)
即视查询点pi为ISS关键点.
为使特征描述子保持旋转、平移和尺度不变性,同时减少描述子匹配的时间, 受文献[15]启发, 采用BSHOT算法对所提取关键点的邻域特征进行描述.由于BSHOT描述子是对SHOT描述子的二进制编码,故须计算SHOT描述子, 步骤如下:
1) 求解邻域协方差矩阵
(4)
式中r为邻域半径,di为邻域点qi到关键点q的欧氏距离, 并根据所求特征向量构建局部参考坐标系;
2) 求解协方差矩阵按降序排列的特征值(λ1>λ2>λ3)及相互正交的单位特征向量(v1,v2,v3),分别对应于局部参考坐标系的x+,y+,z+方向.以局部参考坐标系为基准建立关键点的球形邻域,沿球的径向、经度、纬度方向分别划分为2, 8, 2个部分,则整个球形邻域被划分为如图2所示的32个子空间;
图2 BSHOT空间结构Fig.2 Space structure of BSHOT
3) 统计子空间内各邻域点qi与关键点q法向量间的夹角余弦值cosθ并计入11维直方图, 再组合32个子空间的直方图, 构成352维的SHOT高维直方图特征.
一个352维的SHOT描述子需要1 408 B的存储空间, 将其编码为352位(即44 B)的二进制描述子,不仅可降低描述子所需存储空间,而且能提高特征匹配的效率.对于一个SHOT描述子Si,i={0,1,2,…,351},Si中的每个值都是0和1之间的十进制数.假设以Bi表示BSHOT描述子,i={0,1,2,…,351}, 且Bi的每个值都为二进制数0或1.从Si的初始端连续取4个值{S0,S1,S2,S3},对应位被编码为BSHOT描述子{B0,B1,B2,B3}.
设Ssum=S0+S1+S2+S3, 若Si的4个值都为0, 则{B0,B1,B2,B3}为{0,0,0,0}; 若有一个值Si∈{S0,S1,S2,S3}超过Ssum的90%, 则对应的Bi编码为1,如S2值大于Ssum的90%, 则{B0,B1,B2,B3}编码为{0,0,1,0}.若上述情况不成立,再检查任意两个值的和是否大于Ssum的90%, 如若S0+S1超过Ssum的90%, 则将{B0,B1,B2,B3}编码为{1,1,0,0}.以此类推, BSHOT描述子通过上述迭代方式将SHOT描述子由十进制编码为二进制.
经过计算点云分辨率、提取ISS关键点及计算BSHOT描述子后,可得到源点云和目标点云中各片点云的关键点及其特征向量.由于点云密度不一、拓扑结构复杂且存在噪声点,所以不同点云中的同一个关键点的描述子会出现偏差,当利用双向汉明距离进行对应关键点匹配时易产生错误的匹配点对,这将会影响后续初始变换矩阵的求解.
为了解决可能出现的误匹配问题,现采用RSC算法剔除误匹配点对,并通过随机取样剔除局外点,建立一个仅由局内点数据组成的点集.下面对输入点云数据随机选取点集计算预设的模型参数.假设模型的初始状态至少含3个局内点,局内点在数据中的占比
ε=ni/N,
(5)
式中ni为局内点数,N为局内点与局外点之和.若以(1-ε)k表示经k次迭代计算后模型至少采样到一个局外点的概率, 则局内点的采样概率P=1-(1-ε)k,P越接近1表示模型中局内点数量越多.在保证一定置信度下迭代结束后估计出最优的模型参数.
剔除错误的对应点对后,利用剩余的正确对应点关系, 以目标点云作为参考, 求解两片点云间的初始变换矩阵M, 将源点云转换到以目标点云为标准的坐标系中.两片点云间的变换关系为
(6)
式中(x′i,y′i,z′i)为目标点云上的点,(xi,yi,zi)为源点云上与目标点云对应的点.由于变换矩阵M中含有12个与旋转或平移相关的参数,所以至少需要6组对应点通过求解线性方程组得出变换矩阵M.
硬件环境为Inter (R) Core i5-4200H@2.8 GHz处理器,8.0 GB运行内存, 操作系统为Windows10 64位.软件环境为cmake-gui, Visual Studio 2017, 开源点云库PCL1.8.1.实验数据采用斯坦福大学不同视角下的dragon和bunny点云模型, dragon000与dragon024的点云数量分别为41 841, 34 836, bunny000与bunny045的点云数量分别为40 256,40 097.首先加载源点云和目标点云,由于BSHOT特征结合了点云的几何分布信息,故须在计算点云分辨率后提取两片点云中的ISS关键点, 然后根据关键点的邻域特征计算BSHOT描述子.相同关键点的特征描述子也近似相同,本文算法剔除误匹配后的对应点如图3所示.
图3 对应点匹配结果Fig.3 Corresponding point matching results
在相同软硬件条件下, 将本文ISS-BSHOT+点到点的ICP算法与传统点到点的ICP算法[6]、SIFT-FPFH+点到点ICP算法[9]、ISS-SHOT+点到面的ICP算法[11]进行对比实验.为了验证本文算法的稳健性,将目标点云沿x轴正方向平移0.25 m后进行配准实验.配准前源点云和目标点云的初始位姿如图3所示.由图3可见,配准前源点云与目标点云初始位姿相差较大且无重叠.本文采用均方根误差(root mean square error, RMSE)作为配准精度的评价指标,衡量两片点云在空间上的接近程度,结果如表1所示.图4为不同算法的粗配准结果,其中红色点云为目标点云,蓝色点云为经过初始变换矩阵变换后的源点云.
由图4和表1可知:本文粗配准算法的精度略低于文献[11],较文献[9]高47.4%;虽然本文粗配准算法在将SHOT特征描述子转换为BSHOT特性描述子时会丢失部分信息,使得利用BSHOT特征描述子查找到的对应点数量低于SHOT描述子,进而导致粗配准精度略低,但是仍能为后续精配准提供良好的初始位姿.
为兼顾配准效率和精度,本文采用点到点的ICP算法进行精配准,设定ICP最大迭代次数为10,几种不同算法的精配准结果如图5所示,其中红色点云为目标点云,蓝色点云为经过ICP所求变换矩阵变换后的源点云.由图5和表1可知:传统ICP配准算法的配准精度最低且用时较长;本文配准算法用时最少,配准效率较文献[9]和文献[11]分别高32.5%,51.6%,配准精度较文献[9]高15.4%,但略低于文献[11].其原因是采用二进制的BSHOT特征描述子进行关键点的特征匹配时速度较快,但粗配准精度略低于文献[11],进而导致精配准精度下降.