朱广堂,叶珉吕
(佛山市测绘地理信息研究院,广东 佛山 528000)
三维激光扫描仪虽然能以较快的速度及较高的精度获取数据,但由于仪器测距、测角误差,物体反射率,人为操作,遮挡及环境变化等原因,扫描获取的点云数据中存在大量的噪声点。现有的去噪算法思想大多是来自目前比较成熟的网格模型和数字图像去噪算法,如经典的拉普拉斯(Laplace)算法[1]、双边滤波算法[2]、高斯滤波算法、平均滤波算法及中值滤波算法等。国内外学者对上述多种算法进行了结合或改进,提出了许多更为有效的算法[3-9],但此类算法本质上都是通过离散点在一定方向上的移动,达到模型光滑的效果。虽然模型拥有“视觉满意”的效果,但这也改变了点位的真实位置,因此现有的方法并不适用于对精度、数据真实性等方面要求较高的测绘领域。同时,目前点云去噪结果一般通过去噪前后的图像对比定性地进行评价,去噪结果无法定量分析,评价结果具有较大的主观性。针对上述问题,本文提出基于移动最小二乘(moving least squares, MLS)曲率特征的点云去噪及定量评价方法,并通过实例验证方法的可行性。
一个好的点云去噪算法应同时具有以下两个特点:①有效去除噪声的同时能够保持原始点云数据的特征信息;②算法复杂程度低,执行效率高[10]。点云数据的曲率、法向量等几何信息能够较好地反映点云的表面特征,为了使去噪结果能够较好地保留点云特征,首先需要计算点云的曲率信息。因为点云是离散的点数据,曲率信息需要在离散点集拓扑结构的基础上通过曲面拟合计算得到。离散点集的拓扑结构一般通过KD-Tree、四叉树、八叉树等方法搜索各点K邻域建立,KD-Tree数据结构针对海量、空间分布不均匀的离散点集拥有较好的空间搜索效率,因此本文采用KD-Tree进行K邻域搜索,构建点集拓扑结构。搜索任意一点P的邻域的流程如下[11]:
(1) 如遇到根节点,则从根结点开始搜索。
(2) 如遇到叶子结点,则采用递归的方法返回所有与P邻近的点,并计算最小距离值。
(3) 通过判断以P为球心、最小距离值为半径的球体与分割平面是否相交来确定分割面另一边上是否有点更接近点P。如果相交,分割面另一边上可能有邻近点,则向上层区域回溯去搜寻更近点;如果不相交,则继续沿着当前分支运行。
(4) 直至完成所有节点的搜索,P的邻域点集建立完成。
所有点的邻域点集构建完成后,可在此基础上计算曲率、法向量等几何信息。移动最小二乘法是对最小二乘法的改进,最早由文献[12]提出用于构建原始三维表面的投影表面,是目前最为精确的点云数据法向量计算方法。计算任意一点Pi的法向量基本流程如下[13]。
设点Pi处的切平面方程为ax+by+cz+d=0,其中a2+b2+c2=1,则平面法线即点Pi的法向量为ni=[abc]T。假定点Pi及其k个邻近点构成点集{Pi,Nk(Pi)},由点到切平面的距离构建点集的观测方程为
BY=D
(1)
其中
式中,(xi,yi,zi)为点Pi的坐标;(xi,j,yi,j,zi,j)(j=1,2,…,k)为点Pi的k个邻近点的坐标;di,0为点Pi到切平面的距离;di,1,di,2,…,di,k分别为点Pi第k个邻近点到切平面的距离。
采用高斯函数确定任意点Pj∈{Pi,Nk(Pi)}到切平面距离来确定该点的权Pi,j
(2)
式中,‖Pj-Pi‖表示Pj到Pi的距离;h是与邻点之间距离相关的一个常数,较大的h表示更大范围内的逼近。
由约束准则
得
BTPBY=0
(3)
由a2+b2+c2=1及式(3)可解得Y,即a、b、c、d的值,进而得到点Pi的法向量ni=[abc]T。
最终根据法向量与曲率的关系,任意一点Pi的曲率可由下式计算
(4)
式中,nNk(j,Pi)表示点Pi第k个邻近点的法向量。
目前还没有成熟、通用的点云数据去噪方法,日常数据处理更多是通过仪器配套软件的去噪功能或在软件中人机交互手动去除,此类方法必然会出现误删、多删、漏删的情况。本文根据不同物体曲率不同或曲率突变的原理,结合上述点云法向量、曲率计算方法,提出基于曲率信息的点云数据去噪方法,其具体步骤如下:
(1) 导入点云数据,利用移动最小二乘法计算每个点的法向量。
(2) 遍历每个点,根据某点及其k个邻近点的法向量,由式(4)计算曲率。
(3) 设定阈值,通过判断每点曲率与阈值的关系,去除噪声点。
算法中的参数k及阈值可根据点云数据的扫描分辨率及物体的形状等因素灵活设定,以达到最优的去噪效果。
使用三维激光扫描仪获取某规则构筑物的点云数据,如图1所示,构筑物表面遮挡物所引起的噪声会干扰构筑物几何特征信息的提取,不利于数据的后续处理与应用,因此必须在尽可能保持原始数据完整性的前提下去除此类噪声。使用本文方法对点云进行去噪,结果如图2所示。
由图2可以看出,本文方法能够在保持原始数据完整性的前提下,去除遮挡等原因引起的噪声,未出现误删、多删、漏删等情况。
现实生活中更多的是具有非平面、不规则特征的物体。使用三维激光扫描仪获取某大楼的点云数据,如图3所示,该大楼表面是一任意曲面,假设贴附于大楼表面的两个大字作为遮挡等原因引起的噪声,使用本文方法对点云进行去噪,结果如图4所示。
由图4(a)可以看出位于大楼表面的两个大字噪声被完整去除,未出现误删、多删、漏删等情况,痕迹平滑,较好地保持了数据的完整性;从图4(b)及图4(c)可以更为明显地看出去噪后的大楼表面为光滑的曲面,与图3相比更能直观地反映去噪结果。由于字体与大楼表面贴附紧密,如果简单使用人机交互的方式去除,将不可避免地出现误删、多删、漏删等情况,特别是字体轮廓与大楼表面相连接的区域,几乎无法分辨哪些是楼体表面的点,哪些是噪声点。
以上两个实例,通过去噪前后的图像对比定性地验证了本文方法的有效性。
如前所述,目前点云数据去噪结果一般是通过去噪前后的图像对比定性地进行评价,缺少定量评价指标。文献[14]通过信息熵理论描述点云的特征信息,证明某点熵值越大,所在区域的无序程度越高,则该点的信息量越大,细节表现越精确[14]。因此,可以通过计算点云数据的熵值,定量评价去噪结果。计算某点熵值的公式为[14-15]
(5)
式中,Ci表示点i的曲率;Cj表示点i邻近点j的曲率;pi与pj分别为点i及j的曲率概率分布。则点云数据整体熵值等于各点所含的熵值之和
(6)
点云数据整体熵值越大,其所含的特征信息量越大,对物体的细节表现越精确。
同时,采用本文方法与应用比较广泛的拉普拉斯算法、双边滤波算法对2.2节中的点云数据进行去噪。去噪后的点云数据经基于信息熵理论的定量评价方法计算后,结果见表1。
表1 不同方法去噪结果对比
从表1可以看出,本文方法去噪结果的熵值较传统的拉普拉斯算法、双边滤波算法大,说明数据去噪的同时原始细节特征保留较多,整体包含的特征信息量较大,具有更高的可行性。
针对海量点云数据,为了在高效去除噪声的同时又能保持原始点云数据的特征信息,本文提出了基于移动最小二乘曲率特征的点云去噪算法。通过两个实例分析表明,本文提出的基于移动最小二乘曲率特征的点云去噪方法,能够在大量保留原始点云特征信息的同时,高效地去除噪声,且算法原理简单、运行效率快,具有一定的可行性和适用性。