刘哲延,刘 璐
(贵州省第一测绘院 贵州 贵阳 550025)
机载激光雷达(light detection and ranging,LiDAR)是一种集激光测距系统、全球定位系统(global navigation satellite system,GNSS)和惯性导航系统(inertial navigation system,INS)于一身的新型测量系统,用于获取被测对象表面三维坐标并生成精确的数字三维模型。机载LiDAR技术具有较大程度穿透植被获取真实地表坐标的能力,被广泛用于大范围地面数字高程模型(digital elevation model,DEM)制作。由机载LiDAR获取的点云数据主要由地面点、非地面点和噪声组成,复杂的地形、形态各异的非地面点以及大量的噪声极大增加了通过滤波获取正确地面点的难度。在实际生产中,通过专业点云处理软件滤波得到的地面点云通常无法直接用于制作DEM,而是需经过人工编辑后才能做进一步使用。因此,从机载LiDAR点云中通过滤波尽可能多地获取正确的面点以减少后期人工编辑具有非常重要的现实意义。
在实际航空摄影过程中,由于存在悬空云雾或仪器自身问题等情况,机载LiDAR点云可能有大量噪点。针对去噪和地面点滤波算法,已经有学者做过大量研究[1-3]。但该类研究主要侧重于理论,多探讨算法本身,与工程实践结合不够紧密,在工程实践中遇到的情况极其复杂,单一算法难以满足实践需求,通常是通过相关商业软件设定滤波算法组合以达到期望效果。Terrasolid软件是目前被广泛用于处理机载LiDAR点云实际工程数据,在国内得到了较多应用[4-5]。采用Terrasolid进行常规组合滤波的流程通常为去除飞点、孤立点和低点,通过软件内置的ground算法(分类地面点)得到地面点,再在此地面点基础上进行建筑、植被等点云的分类。但Terrasolid的地面点分类采用的是渐进三角网加密算法,经典的渐进三角网加密算法基本原理为:首先对数据区域划分规则格网,认为格网中最低点为地面点,并作为初始化的种子点[6-7];然后构建不规则三角网,再遍历非地面点,若三角面的坡度小于该阈值,则计算待判断点到三角面的距离,以及待判断点与距离最近的三角面的顶点的连线与三角面之间的夹角,示意图如图1所示,如果这两者均小于给定的阈值,则将待判断点视为地面点,计入地面点集[8];重复以上步骤,直到不再有新的地面点生成。
图1 经典的渐进三角网加密算法中迭代加密步骤对待判断点计算距离d和夹角α
经典的渐进三角网加密算法虽然普适性较好,但是该算法执行效率较低且对孤立点、飞点、低点、低矮植被、低矮构筑物等较为敏感,导致滤波结果出现错误。本文基于Terrasolid软件,结合点云回波特性、地形数据和噪声的相对关系,提出了一种组合滤波的方法,提高了地面点滤波的准确度,为实际工程中点云组合滤波方法提供参考。
就机载LiDAR点云而言,通过噪声点与地表点的位置关系可以将噪声点分为三类:独立于地表点的噪声点、与地表点相交的噪声点、与地表点混杂的噪声点。少量的孤立独立噪声点可以通过点的离群距离去除,但大量聚集的噪声点却经常被误判为地物,相交噪声点和混杂噪声点仅难以通过噪声点的离群距离剔除,该类噪声点通常由人工编辑去除,这往往耗时费力,且难以明确界定交织点为地面点还是噪声点。如图2所示[4]。
图2 噪声点类别
利用噪声点和地表点的位置关系可以在Terrasolid中快速实现去噪功能,但由于该算法中地表点的提取是与去噪交替进行的,故该算法的去噪效果与地表点的选择有较大关系。机载LiDAR点云数据通常都有多次回波,首次回波一般为树冠、裸露地表、楼顶、高空噪声等区域,末次回波则大多数为地面点,非首、末次回波多为冠层下的树叶、树干、镂空的地面构筑物等。利用回波特性,可以去除大量相交噪声点和混杂噪声点。基于机载LiDAR点云的回波特性,算法如图3所示。
图3 算法流程图
(1)初始化机载LiDAR点云。滤波前初始化准备,把所有类别的点云划分到Default类中。
(2)孤立点、飞点、低点分类。利用孤立点、飞点、低点与点云主体的位置关系,剔除离群较远的噪声点。
(3)利用回波特性构造初步地表点集合。将首次回波的点云和末次回波的点云归为一类,再对该点集提取地表点,这些表面点能够包含绝大部分的地面点和地物点,同时去掉了中间回波的干扰,可以用作去噪的基准面。
(4)以地表点为基准区分高于地表的点、地表附近的点和低于地表的点。借助粗略划分地表点点集合为基准,利用噪声点和地表点的位置关系去噪。
(5)把除噪声以外的点归入Default。后续进行地面点滤波需要去噪后的点,该类点包含了首次提取的地表点和原本Default图层中的点。
(6)提取地面点。使用Terrasolid内置的地面点滤波算法提取地面点。
试验区位于贵州省六盘水市境内,数据由Riegl 580Ⅱ机载激光雷达设备于2021年采集所得。试验区最大点密度为13.5 个/m2,平均点密度5.09 点/m2。取1 km×1 km为试验区,该区域高差达500 m,属于典型的高地形。该区域密林区,植被覆盖高,由于测区高差大,在航摄的时候遇到的云雾,导致有较多的连片悬空噪点和连片低点。如图4、图5所示。
图4 试验区范围
本文所述的组合滤波方法在Terrasolid中实现如表1所示。
运行本文组合滤波方法后,生成对应的DEM山体阴影图,并将本文所得结果和Terrasolid常规组合滤波算法结果对比如图6~图8所示。
图6 常规滤波后DEM山体阴影图(45°俯视图)
图7 常规滤波后DEM山体阴影图
图8 本文组合滤后DEM山体阴影图
从图6~图8中可以看出,本文提出的组合滤波方法相对Terrasolid中传统滤波方法在抵抗噪声和建筑物滤除方面有着更好的效果,但也将部分山体被误分为非地面点。但总体来说,本文组合算法与Terrasolid中传统滤波方法相比,更能兼顾抗噪和建筑滤除的效果。
综上所述,本文基于Terrasolid软件,结合点云回波特性、地形数据和噪声的相对关系,提出了一种组合滤波的方法,并将其与常规采用Terrasolid进行组合滤波的流程方法进行了对比,验证了本文组合滤波方法在地面点分类中更能兼顾抗噪和建筑滤除效果。同时本文给出了Terrasolid中组合滤波方法的具体实现,具有较强的可操作性,可对广大专业技术人员的工程实践提供直接参考。但在工程实践中可能遇到不同的场景,专业技术人员应对不同滤波方案进行比选,选取最合适的点云滤波方法。