Airborne LiDAR 3D Filtering based on Intensity Voxel Primitive

  • WANG Liying , * ,
  • WANG Sheng ,
  • LI Yu
Expand
  • School of Geomatics, Liaoning Technical University, Fuxin 123000, China
WANG Liying, E-mail:

Received date: 2019-02-15

  Request revised date: 2019-07-22

  Online published: 2019-12-25

Supported by

Natural Science Foundation of Liaoning Province of China(20170540419)

Scientific Research Fund of Liaoning Provincial Education Department(LJ2019JL015)

Copyright

Copyright reserved © 2019

Abstract

The existing binary voxel primitive based 3-Dimensional (3D) filtering algorithms for airborne Light Detection And Ranging (LiDAR) data, which use only elevation features, cannot distinguish between connected ground and non-ground objects. As a result, an airborne LiDAR 3D filtering algorithm based on intensity voxel primitive was proposed in the present study. First, airborne LiDAR data were regularized into intensity voxel structure based on computational geometry theory, in which intensity value of the voxel corresponds to the quantized intensity of the LiDAR point(s) within the voxel. Second, based on the theory of 3D connected region construction, the non-zero voxels with the lowest local elevation was selected as ground seeds, and then ground seeds and their connected regions where the voxels are 3D connected and have similar grayscales and slope with seeds were labelled as ground voxels. The proposed algorithm makes comprehensive use of the features of elevation, reflection intensity, and slope, supports 3D filtering in areas where the ground are adjacent to non-ground objects but with different intensities, and provides more effective information for the accurate distinction between the connected ground and non-ground objects. The proposed algorithm is helpful to improve the filtering accuracy and extend voxel primitive based 3D filtering algorithm for more complex scenes. The International Society for Photogrammetry and Remote Sensing (ISPRS) benchmark dataset, which contains a variety of features that is expected to be difficult for automatic filtering, were used to analyze the sensitivity of “spatial adjacency size” parameter in the proposed algorithm and to assess the accuracy of the proposed algorithm quantitatively. Results show: (1) The 51-adjacency was the optimal spatial adjacent size. (2) The average Kappa coefficient of the proposed algorithm was 0.9380, 0.7749, and 0.6866 in relatively flat, steep slope, and discontinuous terrain areas, respectively. (3) In terms of total error, the proposed algorithm improved the accuracy of 7 out of 15 samples and had a higher accuracy than all other binary voxel primitive based 3D filtering algorith-ms on average.

Cite this article

WANG Liying , WANG Sheng , LI Yu . Airborne LiDAR 3D Filtering based on Intensity Voxel Primitive[J]. Journal of Geo-information Science, 2019 , 21(12) : 1945 -1954 . DOI: 10.12082/dqxxkx.2019.190152

1 引言

随着数字城市建设的快速发展,如城市规划及管理、交通管理以及抢险救援等方面都要求快速提供城市地形的3D信息。如何快速、大批量地获取高精度的地形数据是数字城市建设的中心任务之一,也是目前研究的热点。机载LiDAR可直接获取地表地形的3D点云数据及其多次回波、反射强度等辅助信息,该项技术的出现使得地形3D信息的快速提取成为可能。
传统的旨在分离地形信息的机载LiDAR数据地面滤波算法按照采用的基本数据单元(简称基元)类型可分为5类[1]
(1)面向点基元的滤波。算法以点为基元,根据点的高程[2,3,4,5,6,7,8,9,10]、反射强度[6]或坡度[2,4,7-8]等特征,利用分类器完成地面目标点的分离。点基元为数据的原始表达、具有真3D特性,但离散的空间LiDAR点无法明晰表达各点间的邻接和拓扑关系[4],导致算法设计难度增加、执行效率降低。另外,点基元仅具备高程、多次回波、反射强度及其统计特征,导致该类算法的可利用特征固定、不全面。
(2)面向像元基元的滤波。算法首先将LiDAR数据插值为高程[11,12,13,14,15,16,17,18]或强度影像[11,12],然后以像元为基元,依据其高程或反射强度特征,利用分类器完成面向像元的分类。其缺点在于数据插值丢失了大量点云各个LiDAR点之间相互位置关系信息,由此造成的信息损失会影响后继滤波的准确度且其滤波结果为2D 形式,难以直接构建地面目标的3D 几何结构。另外,同点基元类似,像元基元的可利用特征同样不完整。
(3)面向对象基元的滤波。算法首先分割LiDAR数据,然后以各分割对象为基元,依据其强度[19]、高 程[20,21]、形态[22,23,24]、回波[24]、拓扑[23]及其统计特征,利用分类器区分地面目标对应的分割对象。对象基元的特征相对全面,但其特征的可靠性严重依赖于点云分割的质量。
(4)面向剖面基元的滤波。算法首先将LiDAR数据表示为正交剖面的数据结构,然后进行分割及后继的地面目标分离[25,26,27]。其缺点在于剖面能顾及的邻域有限、可利用的特征不足。
(5)面向体元的滤波。算法首先将LiDAR数据规则化为以二值体元为基元的规则空间结构,然后利用地面体元间的3D连通性、基于3D连通集合构建理论完成滤波[28,29]。体元基元具备真3D特性且各体元间隐含有几何拓扑关系,基于体元的滤波方法设计更容易、执行效率更高。另外,虽然体元基元可用的特征类型也不全面,但其可以通过空间邻域构建提取全面且可靠的可利用特征。因而基于体元基元的滤波算法对比其他基元的滤波算法具有明显优势。但是,已有的基于体元基元的滤波算法[28,29]仍存在缺陷。① 未充分利用LiDAR数据的所有信息。二值(1和0,分别依据体元内是否包含LiDAR点赋值)的赋值方案导致算法仅可利用LiDAR数据的高程信息、仅支持高度不同的地面目标的分离。因此必须在滤波算法中引入更有效的目标区分信息,以便使其适用于更复杂的应用场景。② 特征完整性差。现有二值体元滤波算法仅利用了地面体元间空间连通的几何特征,由此导致算法无法有效区分连通的非地面和地面体元,因此有必要引入其他更全面的特征以支持更精确的地面目标区分。
有鉴于此,本文提出将反射强度信息融入体元基元构建LiDAR数据的强度体元模型,进而综合利用高程、反射强度及坡度等特征完成地面滤波。该算法将几何、辐射及几何形态信息用于滤波的辅助决策,在提高滤波精度的同时扩展了基于体元的机载LiDAR 3D滤波方法的应用场景。

2 研究方法与实验数据

本文算法的技术流程如图1所示,包括数据的强度体元模型构建、3D连通区域构建下的地面体元检测等步骤。
图1 提出的算法的技术流程

Fig. 1 Flowchart of the proposed algorithm

2.1 数据的强度体元模型构建

将离散、非规则分布的LiDAR点集表示为以体元为基元的规则空间结构。首先,将场景空间依据体元大小划分成3D体元网格;然后,将点云映射到3D体元网格;最后,根据体元中LiDAR点的平均反射强度的离散化程度为各体元赋值为{0, 1, …, 255}不同强度等级,如图2所示。
图2 强度体元模型构建原理示意

Fig. 2 Schematic diagram of intensity voxel model construction

其中,场景空间可由包含LiDAR数据的轴向平行包围盒确定。体元大小依据LiDAR数据在各个平面的平均点间距[29]由式(1)确定。式(1)中取最小值是为了使构建的体元模型与LiDAR数据间的精度损失更少。
x = y = A xy n , z = min A xz n , A yz n
式中:(Δx, Δy, Δz)为体元沿xyz方向的大小;AxyAxz, Ayz)是LiDAR数据沿XYXZ,YZ)平面投影所得2D点集占据的面积;n为点数。
另外,原始LiDAR数据中通常包含高程、强度异常数据。高程异常数据为非真实目标信息,其存在会影响体元模型构建的准确性(见3.2节),必须在体元模型构建前予以剔除。剔除方案:统计LiDAR数据中各LiDAR点高程值的频次,并以直方图的形式可视化显示;确定与真实目标对应的最高(低)高程阈值TheTle);针对LiDAR数据中各LiDAR点,若其高程值高于The或低于Tle,则判定为高程异常点,予以剔除,否则保留。强度异常数据不能反映目标的真实激光反射能量信息,必须在体元赋值前予以剔除。剔除方案类似于高程异常数据剔除方案,不同的是此处需统计LiDAR数据中各个LiDAR点强度值的频次。

2.2 3D连通区域构建下的地面体元检测

局部地面目标特性包括:连续、光滑的曲面;反射强度值一致性强;坡度变化不大。经强度体元模型离散化,局部地面体元表现为空间连通、强度值近似且坡度变化小的3D连通区域。因此,提出基于3D连通区域构建的思想分离地面体元,即首先利用地面的局部高程最低特性选取地面种子,然后标记地面种子及其3D连通区域为地面体元。
2.2.1 地面种子体元选取
地面体元具有局部高程最低的特征,因此可以从强度体元模型中搜寻高程最低的非0值体元作为地面种子。但是,鉴于实际地面的分段连续、光滑性,本文采用分块选取的方案:在水平方向上依据设定的分块尺寸对强度体元模型进行分块,并取各块内高程最低的非0值体元为种子。其中,分块尺寸可依据LiDAR数据中最大目标结构的尺寸确定,其原因为:尺寸过小,则最大目标结构内部的体元被错判为地面种子;尺寸过大,则地面种子数量减少、无法准确地表达复杂地形。
2.2.2 地面种子的3D连通区域构建
对任一种子体元Vss=1, 2, …, Ns,Ns代表种子体元的个数),搜寻并标记强度体元模型中与其空间连通、强度值接近且坡度变化较小的所有未标记体元,直至完成所有种子的3D连通区域标记,所得为地面体元集合,伪码见算法1。
算法1: LABEL(Vs, u, Lg, V)程序的伪码
Input:Vs=当前地面种子,u=种子体元的强度值,Lg=地面体元标签,V=强度体元模型
Output:被标记的Vs的3D连通区域
1 标记Vs的标签为Lg
2 初始化一个空栈,并把Vs放入栈中
3 if 栈 = Φ, then
4 终止程序
5 else
6 从栈中弹出栈顶元素te
7 标记te的空间邻域内未标记的、强度值接近u、坡度变化小于给定阈值Ts的体元的标签为Lg,并将上述体元存入栈中
8 跳入3
算法1实现的关键在于确定强度值相似准则、坡度阈值及最优空间邻域尺寸。其中,前者可统计获取,方案如下。统计强度体元模型中非0值体元的强度值频率,如图3所示(以Samp23为例)。
图3 Samp23非0值体元的强度值频率直方图

Fig. 3 Non-zero greyscales frequency histogram for the Samp23

图3可知,强度分布呈现多峰性;各目标的分布呈正态多峰分布。因此,可用高斯混合模型[30]图1中的直方图进行拟合及高斯分解,得到各个高斯分布的分布特征参数(均值μ及标准差σ)。鉴于本文的目的是滤波,所以与地面对应的高斯分布的特征参数(μg, σg)及分布范围[g1, g2](可利用局部极小值法确定)被用于确定地面目标的强度范围,具体如下:若令μg-ml×σg=g1,μg+mr×σg=g2,则可唯一确定mlmr。若基于高斯分布的对称性,令乘数m= max{ml, mr},则可由此确定LiDAR点云中地面目标的强度范围[μg-m×σg,μg+m×σg]。因此,算法1中的“强度值接近u”条件实际计算中转换为“位于地面目标的统计强度范围内”。
坡度阈值Ts的自适应确定方案为:对当前体元Vs,检测其空间邻域内是否含有已标记的地面体元,若有,则将已有地面体元之间的最大坡度值确定为该邻域内的Ts;否则,将地Ts设置为90°。其中,基于Vs和其空间邻域内已标记的地面体元计算坡度E
E = r e - r p c e - c p 2 + l e - l p 2 e = 1 , L , t
式中:(rp, cp, lp)、(re, ce, le)分别为Vs和已标记的地面体元坐标;t为空间邻域内已标记地面体元数。
另外,在3D连通区域的标记过程中应用不同的空间邻域尺寸会得到不同的滤波结果。最佳邻域尺寸将在实验中确定。

2.3 实验数据

采用国际摄影测量测量与遥感协会第三工作组提供的专门用于滤波算法测试的数据(http://www.itc.nl/isprswgIII-3/filtertest/)检验提出算法的有效性及可行性。该数据包含了滤波中可能遇到的各种典型困难(表1)。采用网站提供的地面真实数据作为标准数据定量评价提出的算法的精度。
表1 测试数据集特征

Tab. 1 Features of the test dataset

样本 特征 点数/个 密度/点/m-2 地形特征
Samp11 城区 38 010 1.066 陡坡上混合有植被和建筑物
Samp12 城区 52 119 1.036 小物体(车辆等)、混合植被和建筑物
Samp21 城区 12 960 1.096 狭窄的桥梁、植被
Samp22 城区 32 706 1.039 桥梁(西南方向)、通道(东北方向)
Samp23 城区 25 095 1.198 大的不规则形状的建筑物
Samp24 城区 7492 1.172 斜坡
Samp31 城区 28 862 0.976 复杂屋顶、高低混合地物
Samp41 城区 11 231 1.553 大的数据空洞、不规则建筑物
Samp42 城区 42 470 1.084 铁路、高频率的地形起伏
Samp51 山区 17 845 5.509 斜坡、密集植被
Samp52 山区 22 474 5.984 断层、矮小植被
Samp53 山区 34 378 5.689 间断的陡峭地形
Samp54 山区 8608 5.758 地形表面有密集覆盖物
Samp61 山区 35 060 6.382 不连续的陡坡、山脊、沟渠
Samp71 山区 15 645 5.403 间断地形、地下通道

3 实验结果及分析

3.1 空间邻域尺寸的敏感性测试

考察空间邻域尺寸参数对滤波结果的影响,并由此得出最佳空间邻域尺寸。在分块尺寸、地面目标强度范围等参数设置均相同条件下,对各样本采用不同空间邻域尺寸的空间搜索及坡度自适应确定方案,其滤波结果的Kappa系数如表2所示。
表2 不同空间邻域尺寸的滤波算法的Kappa系数

Tab. 2 Kappa coefficients of the proposed algorithm with different spatial adjacency sizes

样本 Kappa系数
6邻域 18邻域 26邻域 56邻域 64邻域 51邻域
城区 Samp11 0.4271 0.5110 0.5303 0.6495 0.6551 0.6729
Samp12 0.8779 0.9042 0.9059 0.9232 0.9228 0.9232
Samp21 0.9303 0.9494 0.9485 0.9473 0.9462 0.9538
Samp22 0.7958 0.8375 0.8624 0.9062 0.9026 0.9064
Samp23 0.8499 0.8855 0.8880 0.9144 0.9148 0.9188
Samp24 0.7274 0.8006 0.8143 0.8643 0.8612 0.8892
Samp31 0.9544 0.9624 0.9646 0.9701 0.9694 0.9718
Samp41 0.9040 0.9386 0.9398 0.9389 0.9389 0.9409
Samp42 0.9571 0.9771 0.9775 0.9790 0.9789 0.9806
山区 Samp51 0.7763 0.7937 0.7911 0.7670 0.7647 0.7751
Samp52 0.2107 0.4741 0.6038 0.7606 0.7607 0.7623
Samp53 0.3138 0.4430 0.4977 0.5049 0.5000 0.5162
Samp54 0.8626 0.9020 0.9062 0.9065 0.9058 0.9086
Samp61 0.5342 0.6877 0.7337 0.7564 0.7485 0.7648
Samp71 0.6652 0.7263 0.7417 0.6917 0.6942 0.7787
平均 0.7191 0.7862 0.8070 0.8320 0.8309 0.8442

注:加粗数值表示不同空间邻域尺寸的滤波算法的Kappa系数最大值。

需要说明的是,本文算法所得为地面目标体元,而标准数据中则为LiDAR点。为了与标准数据对比,以定量评价提出算法的精度,首先统计滤波所得地面体元内的LiDAR点,然后与标准数据对比,进而利用Kappa系数定量评价提出算法的精度。由表2可知:
(1)Kappa系数来看,56邻域为最佳空间邻域尺寸。空间邻域尺寸的增加并不意味着滤波精度的必然提高。原因在于:空间邻域尺寸过大时,非地面体元被误判为地面体元的概率随之增加,从而在滤波结果中引入更多的误差,这可以解释为何64邻域的滤波结果的平均Kappa系数反而比56邻域的有所降低。考虑到滤波算法的精度降低主要由空间邻域过大、空间邻域的顶层体元造成,为了进一步提高滤波精度,实验继续考察了51邻域(56邻域的顶层体元剔除后所剩空间邻域尺寸)的滤波算法的Kappa系数(表2最后一列),其平均值为0.8442,比56邻域精度提高了约1%。可见,51邻域为最佳空间邻域尺寸。
(2) 针对点云密度约为1点/m-2的城区数据,56、51邻域的滤波算法的平均Kappa系数分别为0.8992、0.9064;针对点云密度约为6点/m-2的山区数据,56、51邻域的滤波算法的平均Kappa系数分别为0.7312、0.7510。由此可见,针对上述2种密度的点云数据,51邻域均为最佳空间邻域尺寸。

3.2 滤波结果及讨论

考察异常数据的存在对强度体元模型构建及滤波结果的影响,并给出最优空间邻域尺寸下的滤波结果。以Samp23为例,LiDAR点数为25 094(图4(a)),经高程异常剔除后,25 077个剔除异常数据(图4(b))被映射到大小为135×189×174的强度体元模型中(图4(c)-(d))。其中体元尺寸为1 m ×1 m×0.4 m,图中仅示出了19 344个非0值体元。利用本文提出的滤波算法从强度体元模型中分离地面体元,图4(e)即为分离所得10 351个地面体元。
图4 Samp23及其剔除异常数据、灰度体元模型、滤波结果

Fig. 4 Samp23 and its outliers removal, the grayscale voxel model, and the filtering result

为了考察高程异常数据对滤波结果的影响,实验分别统计了基于原始LiDAR数据和剔除高程异常数据进行体元模型构建及滤波的结果(表3)。
表3 Samp23中高程异常数据对体元模型构建及滤波结果的影响

Tab. 3 Influence of elevation outliers in Samp23 on the construction of the voxel model and filtering results

点云数据 点数/个 场景空间范围/m 体元尺寸/m 体元模型尺寸/个 非0值体元数/个 地面体元数/个 总误差/%
原始 25 094 146×205.9×86.0 1.1×1.1×0.5 135×189×174 18 973 10 319 0.0505
剔除异常数据 25 077 146×105.9×43.5 1.1×1.1×0.4 135×189×110 19 344 10 351 0.0404
表3可知,高程异常数据的存在对强度体元模型的影响主要反映在z方向,如场景空间的高度、z方向的体元尺寸及由此导致的强度体元模型层数的差异(174对比110)。因此,基于原始LiDAR数据对比基于剔除高程异常数据所得强度体元模型中虽然非0值体元数差别不大(18 973对比19 344),但前者包含更多的空白体元(4 439 610对比2 806 650),因而基于该数据的处理效率更低。基于上述强度体元模型做滤波,分别分离出10 319个和10 351个地面体元。分别统计上述体元中的激光点,并和标准数据做对比进而评价滤波算法总误差(表3),基于原始点云的滤波总误差(0.0505)明显大于异常数据剔除后的滤波总误差(0.0404)。这也证明了异常数据剔除的必要性。
将本文算法的滤波结果和标准数据对比,进而用I类误差、II类误差、总误差及Kappa系数等测度定量评价提出算法的精度,如表4所示。
表4 算法误差和Kappa系数

Tab. 4 Errors and Kappa coefficients of the proposed algorithm

样本 Ie IIe Te Kappa系数
Samp11 0.1546 0.1693 0.1609 0.6729
Samp12 0.0410 0.0356 0.0384 0.9232
Samp21 0.0107 0.0344 0.0160 0.9538
Samp22 0.0363 0.0450 0.0405 0.9064
Samp23 0.0246 0.0580 0.0404 0.9188
Samp24 0.0250 0.0933 0.0438 0.8892
Samp31 0.0094 0.0194 0.0140 0.9718
Samp41 0.0211 0.0380 0.0296 0.9409
Samp42 0.0123 0.0063 0.0080 0.9806
Samp51 0.0398 0.2016 0.0751 0.7751
Samp52 0.0113 0.2921 0.0409 0.7623
Samp53 0.0133 0.5212 0.0338 0.5162
Samp54 0.0158 0.0714 0.0457 0.9086
Samp61 0.0026 0.3253 0.0137 0.7648
Samp71 0.0221 0.2124 0.0437 0.7787
平均 0.0293 0.1416 0.0430 0.8442

注:Ie表示将地面点错分为非地面点比例;IIe表示将非地面点错分为地面点比例;Te表示错分的激光点的比例。

表4可知:
(1) IeIIeTe范围分别为0.01~0.16、0.01~0.53、0.01~0.17,这说明提出的算法表现为最小化Ie
(2)从Te来看,城区内Samp11总误差最大(0.1609)。为了探究原因,考察了其分布情况,如图5(a)所示。由图5(a)可知,Samp11为阶梯式的陡坡地形,存在部分非地面(如植被)与地面体元高程且强度接近的情况,从而造成非地面目标的误判(即IIe的增大)。另一方面,由于目标遮挡、地形较陡等原因导致Samp11中的地面点分布不均匀,由此导致地面目标形成过多的3D连通区域。要想完整分离地面对应的3D连通区域则需要较小的分块尺寸以便获得更多的地面种子,但是如前所述,分块尺寸不可小于场景内最大结构单元的尺寸。因此很容易造成地面3D连通区域的错分,即Ie的增大。山区内Samp51 Te最大,为0.0751。Samp51和Samp11的地形特征类似,均为地表覆盖有浓密植被的陡坡地形,如图5(b)所示。
图5 样本Samp11、Samp51的I类和II类误差分布顶视图

Fig. 5 Top view of type I and II error distribution of the proposed algorithm for Samp11 and Samp51

(3)从Kappa系数看,所有样本的平均Kappa系数高达0.8442。但该测度并不能示出提出的算法针对不同地形的精度,因此将15个样本依据地形特征分为3类进行精度评价(图6)。由图6可知,本文提出算法在相对平坦、陡坡和不连续地形中的平均Kappa系数分别为0.9380、0.7749、0.6866,说明算法在相对平坦地形表现最佳、陡坡地形中表现较令人满意、不连续地形中表现最差。原因在于:激光点分布相对均匀,与提出的算法的理论假设切合度较高,因而滤波精度很高。不连续地形中地面点分布也不均匀,造成了较大的分块尺寸(由场景内最大目标结构的尺寸决定)无法适应较多地面种子需求间的矛盾,因此该类地形的精度稍差,如Samp53的Kappa系数最低(0.5162)。陡坡地形精度较差的原因见上文Samp11的误差分析。
图6 3类地形的精度对比

Fig. 6 Accuracy comparison between the three terrain types

通过本文算法与其他经典滤波算法误差对比(表5)可知,相比于经典的Axelsson算法[18],本文算法降低了7个样本的总误差。对比Wang等[30]、王丽英等[31]算法,本文算法在各种地形中均表现更优。
表5 本文算法与经典滤波算法总误差比较

Tab. 5 Total error comparison between the proposed algorithm and three classical filtering algorithms

样本 Axelsson[18] Wang等[30] 王丽英等[31] 本文算法
Samp11 0.1076 0.1949 0.1862 0.1609
Samp12 0.0325 0.0402 0.0487 0.0384
Samp21 0.0425 0.0205 0.0174 0.016
Samp22 0.0363 0.0497 0.0460 0.0405
Samp23 0.0400 0.0591 0.0504 0.0404
Samp24 0.0442 0.0634 0.0623 0.0438
Samp31 0.0478 0.0158 0.0203 0.0140
Samp41 0.1391 0.0217 0.0209 0.0296
Samp42 0.0162 0.0107 0.0098 0.0080
Samp51 0.0272 0.0809 0.0818 0.0751
Samp52 0.0307 0.0490 0.0424 0.0409
Samp53 0.0891 0.0346 0.0453 0.0338
Samp54 0.0332 0.0562 0.0508 0.0457
Samp61 0.0208 0.0193 0.0231 0.0137
Samp71 0.0163 0.0542 0.0614 0.0437
平均 0.0482 0.0513 0.0511 0.0430

注:Axelsson[2]为商业软件TerraScan采用的算法;Wang等[30]、王丽英等[31]为基于二值体元模型的3D滤波算法。

本文算法的滤波结果由输入参数(如坡度阈值Ts、分块尺寸、地面的统计强度范围和空间邻域尺寸)共同决定。其中,分块尺寸、Ts和地面的统计强度范围是根据实际数据源设置的。使用本文给出的设置方案,可以很容易地确定“数据源”型参数,并不能限制本文算法应用到其他源数据的适用性。空间邻域尺寸针对密度分别为1/m-2和6点/m-2的点云数据可以直接使用51邻域,因为实验已证明51邻域为上述2种点云密度数据的滤波最佳邻域尺寸。

4 结论与讨论

4.1 结论

本文针对机载LiDAR数据提出了一种基于强度体元模型的3D滤波算法。该算法首先将机载LiDAR点云数据规则化为强度体元结构。然后,基于地面点的局部高程最低特性选取地面种子体元。最后,标记与地面种子空间连通、强度值接近且坡度变化不大的非0值体元为地面体元。其中,坡度阈值通过统计空间邻域内已有地面种子间的坡度确定;强度值的一致性通过拟合地面强度信息的高斯混合模型确定。
基于ISPRS提供的标准测试数据定量评价了算法精度,并和已有经典滤波算法进行对比。实验结果表明:① 平均Kappa系数在相对平坦、陡坡及不连续地形分别为0.9380、0.7749和0.6866,在激光点分布相对均匀的平坦地区表现最优、激光点分布不均匀的陡坡及不连续地形中表现则较差;② 算法对比经典的Axelsson算法[2]改进了15个样本中的7个样本精度;③ 对比已有的基于二值体元模型的经典滤波算法,算法在复杂城区环境、不连续及坡度急剧改变地形中均表现更优。

4.2 讨论

该算法以3D连通区域构建理论为基础,将反射强度及坡度信息作为辅助信息用于点云数据滤波,其对比基于二值体元模型的3D滤波算法[30]可为地面和非地面目标的区分提供更有效的信息,因而显著提高了滤波精度并将基于体元模型的3D滤波方法扩展应用至更复杂的场景。但是,该算法仅支持高程或强度不同的地面目标的有效分离。在后续的研究中可以尝试和其他光学影像等数据融合生成多值3D图像(即体元赋值体元内激光点的RGB等色彩信息)以支持更复杂场景的分类。另外,针对算法在陡坡和不连续地形中表现较差,后续研究可考虑根据点云局部密度情况自适应设置体元分辨率方案以提高上述2种地形的点云数据形成3D连通区域的准确性,从而有效提高滤波精度。
[1]
张继贤, 林祥国, 梁欣廉 . 点云信息提取研究进展和展望. 测绘学报, 2017,46(10):1460-1469.

[ Zhang J X, Lin X G, Liang X L . Advances and prospects of information extraction from point clouds[J]. Acta Geodaetica et Cartographica Sinica, 2017,46(10):1460-1469. ]

[2]
Axelsson P . DEM generation from laser scanner data using adaptive TIN models[J]. International Archives of Photogrammetry and Remote Sensing, 2000,33(B4):111-118.

[3]
Zhang W, Qi J, Wan P , et al. An easy-to-use airborne LiDAR data filtering method based on cloth simulation[J]. Remote Sensing, 2016,8(6):501.

DOI

[4]
郭波, 黄先锋, 张帆 , 等. 顾及空间上下文关系的Joint Boost点云分类及特征降维[J]. 测绘学报, 2013,42(5):715-721.

[ Guo B, Huang X F, Zhang F , et al. Points cloud classification using JointBoost combined with contextual information for feature reduction[J]. Acta Geodaetica et Cartographica Sinica, 2013,42(5):715-721. ]

[5]
Liu X Q, Chen Y M, Cheng L , et al. Airborne laser scanning point clouds filtering method based on the construction of virtual ground seed points[J]. Journal of Applied Remote Sensing, 2017,11(1):016032.

DOI

[6]
Morsy S, Shaker A, El-rabbany A . Multispectral LiDAR data for land cover classification of urban areas[J]. Sensors, 2017,17(5):958.

DOI PMID

[7]
Nie S, Wang C, Dong P L , et al. A revised progressive TIN densification for filtering airborne LiDAR data[J]. Measurement, 2017,104:70-77.

DOI

[8]
黄作维, 刘峰, 胡光伟 . 基于多尺度虚拟格网的LiDAR点云数据滤波改进方法[J]. 光学学报, 2017,37(8):346-355.

[ Huang Z W, Liu F Z, Hu G W . Improved method for LiDAR point cloud data filtering based on hierarchical pseudo-grid[J]. Acta Optica Sinica, 2017,37(8):346-355. ]

[9]
Chen C, Li Y, Zhao N , et al. A fast and robust interpolation filter for airborne lidar point clouds[J]. Plos One, 2017,12(5):1-20.

DOI PMID

[10]
朱笑笑, 王成, 习晓环 , 等. 多级移动曲面拟合的自适应阈值点云滤波方法[J]. 测绘学报, 2018,47(2):153-160.

[ Zhu X X, Wang C, Xi X H , et al. Hierarchical threshold adaptive for point cloud filter algorithm of moving surface fitting[J]. Acta Geodaetica et Cartographica Sinica, 2018,47(2):153-160. ]

[11]
Wang C K, Seng Y H, Chu H J . Airborne dual-wavelength LiDAR data for classifying landcover[J]. Remote Sensing, 2014,6(1):700-715.

DOI

[12]
Axelsson A, Linberg E, Olsson H . Exploring multispectral ALS data for tree species classification[J]. Remote Sensing, 2018,10(2):183.

DOI

[13]
Hui Z, Hu Y, Yevenyo Y , et al. An improved morphological algorithm for filtering airborne LiDAR point cloud based on multi-level kriging interpolation[J]. Remote Sensing, 2016,8(5):1-16.

DOI

[14]
Pingel T J, Clarke K C, Mcbride W A . An improved simple morphological filter for the terrain classification of airborne LIDAR data[J]. ISPRS Journal of Photogrammetry and Remote Sensing, 2013,77(1):21-30.

DOI

[15]
Chen Q, Gong P, Baldocchi D , et al. Filtering airborne laser scanning data with morphological methods[J]. Photogrammetric Engineering and Remote Sensing, 2007,73(2):175-185.

DOI

[16]
隋立春, 张熠斌, 柳艳 , 等. 基于改进的数学形态学算法的LiDAR点云数据滤波[J]. 测绘学报, 2012,39(4):390-396.

[ Sui L C, Zhang Y B, Liu Y , et al. Filtering of airborne LiDAR point cloud data based on the adaptive mathematical morphology[J]. Acta Geodaetica et Cartographica Sinica, 2012,39(4):390-396. ]

[17]
潘锁艳, 管海燕 . 机载多光谱LiDAR 数据的地物分类力法[J]. 测绘学报, 2018,47(2):198-207.

[ Pan S Y, Guan H Y . Object classification using airborne multispectral LiDAR data[J]. Acta Geodaetica et Cartographica Sinica, 2018,47(2):198-207. ]

[18]
Huo L, Silva C A, Klauberg C , et al. Supervised spatial classification of multispectral LiDAR data in urban areas[J]. Plos One, 2018,13(10):e0206185.

DOI PMID

[19]
Zou X, Zhao G, Li J, , et al. 3D land cover classification based on multispectral lidar point clouds[J]. International Archives of Photogrammetry Remote Sensing and Spatial Information Sciences 2016, XLI-B1:741-747.

[20]
Matikainen L, Karila K, Hyyppa J , et al. Object-based analysis of multispectral airborne laser scanner data for land cover classification and map updating[J]. ISPRS Journal of Photogrammetry and Remote Sensing, 2017,128:298-313.

DOI

[21]
Hu J, Yang L, Shen J , et al. Filtering of LIDAR based on segmentation[J]. Geomatics and Information Science of Wuhan University, 2012,37(3):318-321.

[22]
Filin S . Surface clustering from airborne laser scanning data[J]. International Archives of Photogrammetry Remote Sensing and Spatial Information Sciences, 2002,34(3/A):119-124.

[23]
Yan M, Blaschke T, Liu Y , et al. An object-based analysis filtering algorithm for airborne laser scanning[J]. International Journal of Remote Sensing, 2012,33(22):7099-7116.

DOI

[24]
Lin X, Zhang J . Segmentation-based filtering of airborne LIDAR point clouds by progressive densification of terrain segments[J]. Remote Sensing, 2014,6(2):1294-1326.

DOI

[25]
Shan J, Sampath A . Urban DEM generation from raw lidar data: a labeling algorithm and its performance[J]. Photogrammetric Engineering and Remote Sensing, 2005,71(2):217-226.

DOI

[26]
Hu X, Li X, Zhang Y . Fast filtering of LiDAR point cloud in urban areas based on scan line segmentation and GPU acceleration[J]. IEEE Geoscience and Remote Sensing Letters, 2013,10(2):308-312.

DOI

[27]
Sithole G, Vosselman G. Automatic structure detection in a point-cloud of an urban landscape[C]//Proceedings of the 2nd GRSS/ISPRS Joint Workshop on Remote Sensing and Data Fusion over Urban Areas, Berlin: IEEE, 2003: 67-71.

[28]
Wang L, Xu Y, Li Y . Aerial LIDAR point cloud voxelization with its 3D ground filtering application[J]. Photogrammetric Engineering and Remote Sensing, 2017,83(2):95-107.

DOI

[29]
王丽英, 徐艳, 李玉 . 机载LiDAR点云体元化及其在3D滤波中的应用[J]. 仪器仪表学报, 2018,39(7):173-182.

[ Wang L Y, Xu Y, Li Y . Aerial LIDAR point cloud voxelization and its 3D filtering application[J]. Chinese Journal of Scientific Instrument, 2017,83(2):95-107. ]

[30]
赵泉华, 石雪, 王玉 , 等. 可变类空间约束高斯混合模型遥感图像分割[J]. 通信学报, 2017,38(2):34-43.

[ Zhao Q H, Shi X, Wang Y , et al. Remote sensing image segmentation based on spatially constrained gaussian mixture model with unknown class number[J]. Journal on Communications, 2017,38(2):34-43. ]

Outlines

/