多分辨率多层逐次加点LiDAR滤波算法
【专利摘要】本发明公开了多分辨率多层逐次加点LiDAR滤波算法。该算法包括三层滤波：在初始层，首先借助局部最小值法选择部分地面点作为种子点，然后借助多面函数插值算法构建一定分辨率DEM，再次计算每个采样点与对应DEM网格点高程差，最后将高程差小于初始设定阈值的点为标定为地面点；重复上述过程直至没有采样点被划分为地面点；在第二层及第三层滤波中，首先在前层基础上将DEM分辨率提高一倍，并提高对应残差阈值，然后重复前层计算过程直至没有采样点被标定为地面点，停止迭代，输出地面点数据。本发明的平均计算精度稍高于自适应TIN算法，远高于2004以前发展的7种滤波算法，且高于近年来发展的非参数逐次滤波等算法。
【专利说明】多分辨率多层逐次加点LiDAR滤波算法
【技术领域】
[0001]本发明涉及一种多分辨率多层逐次加点LiDAR滤波算法。
【背景技术】
[0002]LiDAR以其高速度、高精度及高密度等特点已成为DEM数据获取技术中不可或缺的地面补充系统。LiDAR原始数据主要包括地面点和非地面点，故将LiDAR采样数据用于DEM构建前必须进行滤波处理，即剔除掉原始数据中的非地面点。其中，滤波算法的计算效率显著影响DEM构建精度。
[0003]根据LiDAR数据中非地面点滤波基本原理，传统滤波算法可分为插值滤波、形态学滤波、坡度滤波及分割滤波等。大量实例分析表明，每种算法在地形平坦及景观简单区域均可以获得较好的滤波效果，但复杂地形及大密度分布植被显著影响滤波算法精度。为此，众多学者对传统滤波算法进行了大量改进，以期提高滤波效果。其中，插值滤波算法以对复杂地形高效适应能力广泛应用于LiDAR滤波。例如，Kraus and Pfeifer提出了基于线性推估法的滤波算法。该算法首先假设所有采样点具有相同的权，并基于线性推估构建初始DEM ;然后计算每个点残差,并将残差为负值的点赋予较大权重，正值的点赋予较小权重。重复上述过程直至迭代收敛。Axelsson发展了基于自适应TIN插值的滤波算法。该算法首先获取部分地面点作为种子点，构建TIN ;然后比较采样点与对应TIN距离和角度，当两者小于设置阈值时，对应点即为地面点。重复上述过程直至迭代收敛。Haugerud andHarding提出了基于TIN插值的逐次减点算法(Virtual Deforestation, VDF)。VDF首先以所有采样点为已知点，基于TIN构建DEM ;然后比较每个采样点与其所在的3X3移动窗口平均值差值，将差值大于初始设定阈值点剔除；最后将剩余点作为已知点重复上述过程，直至收敛。为了进一步提高VDF滤波精度，vans and Hudak发展了一种多尺度曲率分类算法(MCC)。与VDF相比，MCC主要有三点改进:1)用高精度样条插值代替TIN插值；2)插值分辨率从最低层到最高层由高到低变化；3)同时，剔除非地面点的曲率阈值从低层到高层逐次变大。Mongus and 2alik提出了一种非参数VDF算法。在初始循环中，该算法首先选择部分地面点作为控制点，并基于薄板样条(Thin Plate Spline,TPS)构建DEM;然后比较每个采样点与DEM差值，将差值大于设置阈值的点剔除；在下一循环中，该算法降低DEM网格分辨率，再次检核采样点残差。重复上述过程直至收敛。
[0004]为了客观比较各种算法滤波精度，ISPRS发布了 15组LiDAR实验数据。其中，每组数据均准确表明每个采样点属性(地面点或非地面点)。Sithole and Vosselman比较了 8种滤波算法精度，认为Axlesson提出的自适应TIN算法平均精度远高于其它7种算法。后续尽管相关研究人员对传统算法进行了改进或者提出了一些其它新滤波算法，但以ISPRS提供的15组实验数据计算表明，它们平均计算精度均低于自适应TIN精度。
【发明内容】
[0005]本发明要解决的技术问题是提供一种能够提高滤波算法精度的多分辨率多层逐次加点LiDAR滤波算法。
[0006]为了解决上述技术问题，本发明采用的技术方案是:多分辨率多层逐次加点LiDAR滤波算法，设初始层生成的DEM网格分辨率为h，残差阈值为t，层计数N=I ;具体计算步骤如下:
[0007](I)首先在初始层中，以所有采样点为原始数据，以局部最小值法选择地面点作为种子点；其中移动窗口尺寸要稍大于测区最大地面建筑物尺寸，确保选择的种子点全部为地面点；
[0008](2)基于多面函数构建网格分辨率为h的DEM ；
[0009](3)逐次计算采样点与对应DEM网格点所在的3X3移动窗口中每点高程差，即残差，因此每个采样点有9个残差；如果9个残差中有4个以上残差均小于设定残差阈值，则标定该点为地面点，并将其添加到种子点中；此过程可以有效减弱种子点中非地面点对计算结果影响；
[0010](4)以步骤(3)更新后的种子点作为已知数据，基于局部最小值法重新选择种子点；其中局部最小值值法窗口尺寸为h ;执行步骤(4)的原因为:当LiDAR种子点分布不均匀(任意两点非常接近，而与其它点距离较远)，且以它们为已知数据构建DEM时，多面函数系数矩阵极可能为奇异矩阵，进而影响DEM插值精度；因此，该步骤目的是避免多面函数系数矩阵奇异；
[0011](5)重复步骤(2)至(4)，直至没有采样点被划分为地面点；
[0012](6) N=N+1 ;当N〈4时，提高DEM网格分辨率，增大残差阈值，即h=h / 2，t=t+0.1 ；重复步骤(2)至(5)直至收敛。
[0013]本发明的有益效果是:
[0014]基于传统的插值滤波技术发展了一种多分辨率多层滤波算法(MPioD)。以ISPRS提供的15组数据为实例分析表明，MProD平均计算精度稍高于自适应TIN算法，远高于2004以前发展的7种滤波算法，且高于近年来发展的非参数逐次滤波等算法。
【专利附图】
【附图说明】
[0015]下面结合附图和【具体实施方式】对本发明作进一步详细的说明。
[0016]图1是本发明多分辨率多层逐次加点LiDAR滤波算法实施例的流程示意图。
【具体实施方式】
[0017]I多分辨率多层逐次加点滤波算法原理
[0018]该方法包括三层滤波:在初始层中，首先以局部最小值法选择部分地面点作为种子点，然后以多面函数构建一定分辨率DEM，最后将残差小于设定阈值的采样点标定为地面点。重复上述过程直至收敛。在后续层中，DEM分辨率和残差阈值逐次提高，并重复初始层中的计算过程直至没有采样点被标定为地面点，停止迭代，输出地面点数据。设该方法初始层生成的DEM网格分辨率为h，残差阈值为t，层计数N=I。具体计算步骤如下:
[0019](I)首先在初始层中，以所有采样点为原始数据，以局部最小值法选择地面点作为种子点。其中移动窗口尺寸要稍大于测区最大地面建筑物尺寸，确保选择的种子点全部为地面点。[0020](2)基于多面函数构建网格分辨率为h的DEM。
[0021](3)逐次计算采样点与对应DEM网格点所在的3X3移动窗口中每点高程差，即残差，因此每个采样点有9个残差。如果9个残差中有4个以上残差均小于设定残差阈值，则标定该点为地面点，并将其添加到种子点中。此过程可以有效减弱种子点中非地面点对计
算结果影响。
[0022](4)以步骤(3)更新后的种子点作为已知数据，基于局部最小值法重新选择种子点。其中局部最小值值法窗口尺寸为h。执行步骤(4)的原因为:当LiDAR种子点分布不均匀(任意两点非常接近，而与其它点距离较远)，且以它们为已知数据构建DEM时，多面函数系数矩阵极可能为奇异矩阵，进而影响DEM插值精度。因此，该步骤目的是避免多面函数系数矩阵奇异。
[0023](5)重复(2)-(4)直至没有采样点被划分为地面点。
[0024](6) N=N+1 ;当％4时，提高DEM网格分辨率，增大残差阈值，即h=h / 2,t = t+0.1 ；重复(2)-(5)直至收敛。
[0025]试验流程如图1所示。
[0026]2试验分析
[0027]以ISPRS发布的15组数据进行试验分析，试验精度评价指标为基于分类误差矩阵计算的第一类误差(Tl)、第二类误差(T2)、综合误差(TO)以及kappa系数。第一类误差表示将地面点分类为非地面点比例；第二类误差表示将非地面点分类为地面点比例。设误差矩阵如表 I 表示，则 Tl = b / (a+b) ；T2=c / (c+d) ；T0=(b+c) / (a+b+c+d)。MProD 计算结果如表2所示。由表2可见，MProD在samp31区域计算精度最高，其综合误差为1.11%，Kappa系数为97.76% ;在sampll区域综合误差最大，为13.01 在samp53Kappa系数最小，为46.69%。计算结果表明,MProD在地形平坦,场景简单的samp31计算精度最高，而在地形复杂、植被覆盖密度较高sampll区域计算精度也较高。
[0028]表1误差矩阵
[0029]
【权利要求】
1.多分辨率多层逐次加点LiDAR滤波算法，其特征在于，设初始层生成的DEM网格分辨率为h，残差阈值为t，层计数N=I ;具体计算步骤如下: (1)首先在初始层中，以所有采样点为原始数据，以局部最小值法选择地面点作为种子点；其中移动窗口尺寸要稍大于测区最大地面建筑物尺寸，确保选择的种子点为地面点； (2)基于多面函数构建网格分辨率为h的DEM； (3)逐次计算采样点与对应DEM网格点所在的3X3移动窗口中每点高程差，即残差，因此每个采样点有9个残差；如果9个残差中有4个以上残差均小于设定残差阈值，则标定该点为地面点，并将其添加到种子点中；此过程可以有效减弱种子点中非地面点对计算结果影响； (4)以步骤(3)更新后的种子点作为已知数据，基于局部最小值法重新选择种子点；其中局部最小值值法窗口尺寸为h ； (5)重复步骤(2)至(4)，直至没有采样点被划分为地面点； (6)N=N+1 ;当N<4时，提高DEM网格分辨率，增大残差阈值，即h=h / 2，t=t+0.1 ;重复步骤⑵至(5)直至收敛。
【文档编号】G06T17/05GK103700142SQ201310653679
【公开日】2014年4月2日 申请日期:2013年12月3日 优先权日:2013年12月3日
【发明者】陈传法, 阳凡林, 戴洪磊, 田茂义, 曹学伟 申请人:山东科技大学