基于区域预测的LiDAR点云数据形态学滤波方法
【专利摘要】本发明涉及基于区域预测的LiDAR点云数据形态学滤波方法。方法包括以下步骤：步骤一，构建滤波初始表面模型；步骤二，去除初始表面模型的粗差点，得到滤波表面模型；步骤三，对滤波表面模型进行分块：步骤四，滤波处理。本发明的方法克服了现有技术中形态学滤波不能很好对复杂场景滤波的缺点,以及使用恒定的人工设置的地形参数进行滤波处理而导致的自适应不强的问题，利用了分块原理和使用各个分块区域预测地形参数s能够根据地形的起伏情况调整高差阈值等优点，从而能够自适应地对点云数据进行滤波，最终获得较好的滤波结果。
【专利说明】基于区域预测的LiDAR点云数据形态学滤波方法
【技术领域】
[0001]本发明属于机载LiDAR数据处理【技术领域】，更进一步涉及机载LiDAR数据滤波【技术领域】中基于区域预测的形态学点云数据滤波算法。可应用于机载LiDAR扫描地面得到的点云数据，可以有效的将点云数据滤波并提取DEM。
【背景技术】
[0002]机载LiDAR是一种主动式遥感技术，能够扫描大面积区域直接获取地面三维信息，其数据具有高精度、高密度、离散不规则等特点。
[0003]LiDAR点云数据包含有真实地面、建筑物和植被等多种信息，对点云数据进行滤波得到真实地面点，形成数字高程模型(DEM)具有重要意义。目前许多点云滤波算法采用基于形态学的滤波方法。
[0004]赵明波等人在文献“赵明波，何峻，田军生，付强，‘基于改进的渐进多尺度数学形态学的激光雷达数据滤波方法’ [J]，光学学报，33 (3)，292-301 (2013) ”中提出了改进的基于渐进形态学的激光雷达数据滤波算法。该方法的通过改进形态学开运算来处理大面积空白区域的点云数据，改进的开运算具体是腐蚀运算时把当前网格邻域内的空白区域高程值用最大高程值替代，同理，膨胀运算则是用最小高程值替代。该方法首先将原始激光雷达点云数据进行格网化并剔除粗差，然后使用窗口尺寸混合增长的方式和改进的开运算对格网进行迭代的形态学滤波，最终得到原始点云数据的分类。此算法中需人工设置大量参数，自适应性不强。
[0005]隋立春等人在文献“隋立春，张熠斌，柳艳等，‘基于改进的数学形态学算法的LiDAR点云数据滤波’ [J]，测绘学报，39 (4)，390-396 (2010) ”中提出了改进的形态学LiDAR数据滤波方法。此方法对形态学算法进行扩展和改进，提出了针对不同地形特点的自适应滤波算法，具体的是在数学形态学“开”算子的基础上，提出了增加一个“带宽”参数用于点云数据滤波的方法，最后对再山包地形和山坡地形进行了判别和处理。该方法提出的动态滤波结合“带宽”参数的改进的数学形态学算法可以有效地用于LiDAR点云数据的滤波，并在一定程度上避免和降低了地面存在山包时错误滤波的后果。但是此算法没有对ISPRS提供的测试数据进行实验，无法与其他算法比较，并且算法在不同地形、不同高度和不同形状山包情况下还要进行进一步试验和完善。
[0006]综上所述，基于形态学的点云滤波方法在点云滤波方面虽然获得了较好的效果，但是一般来说，基于形态学的点云滤波方法中开运算后的高差阈值的确定需要许多参数，这些参数需根据不同地形设置，并且假设地形坡度恒定，针对复杂场景这是不合理的，导致算法自适应性不强，给对不同地形的区域滤波带来困难。
【发明内容】
[0007]本发明的目的是针对现有技术基于形态学的点云数据滤波方法中需人工设置滤波阈值导致自适应不强的缺点，提出了基于区域预测的形态学LiDAR点云数据滤波算法。[0008]在本发明中，对点云数据的基于区域预测的滤波原理包括了四部分:点云数据格网化、去除粗差点、滤波准则和地面点的确定。
[0009]本发明充分利用形态学滤波算法自下而上、从局部出发扩展到全局的特点和速度快、效率较高的优点，同时利用标准差代表地形的起伏情况的特点，使用分块区域的标准差来预测地形坡度参数的取值以进行滤波。最终得到较好的滤波结果，又能最大限度地保留实验区域的地形特征。
[0010]为此，本发明提供的基于区域预测的LiDAR点云数据形态学滤波方法包括以下步骤:
[0011 ] 步骤一，构建滤波初始表面模型；
[0012]步骤二，去除初始表面模型的粗差点，得到滤波表面模型；
[0013]步骤三，对滤波表面模型进行分块:
[0014]按适当尺寸对滤波表面模型进行分块，其中分块η记为Rn，η = 1,2, 3，…，N，N为分块的总数；该分块Rn中包含M个点云数据，其中点云数据m记为rm，m = I, 2, 3，…，M ;
[0015]步骤四，滤波处理:
[0016]对各分块分别进行滤波处理，其中对分块Rn，η = 1,2, 3，…，N的处理过程如下:
[0017]从最小结构元素窗口尺寸开始执行以下操作，之后结构元素窗口的半径每增加一米，执行一次以下操作直至最大结构元素窗口尺寸:
[0018]分别对分块Rn中所有点云数据进行开运算，得到各点云数据开运算后的高程值，
[0019]分别计算各点云数据的高程值与其自身的开运算后的高程值的差值，差值大于6?的点云数据标记为非地面点，且将各点云数据的高程值更新为开运算后的高程值，其中:en = SnXwXc, w为当前结构元素窗口的半径，c为格网边长，且1米≤c≤5米，sn = aX σ n+b, σ η为分块Rn中所有点云数据高程值的标准差,0.005 ≤a ≤ 0.02,
0.005 ≤ b ≤ 0.02 ；
[0020]滤波处理后未标记的点云数据为地面点。
[0021]进一步，本发明的基于区域预测的LiDAR点云数据形态学滤波方法还包括:
[0022]将步骤四确定的地面点插值得到数字高程模型，
[0023]接着分别确定所有点云数据的初始高程值与数字高程模型之间的垂直距离，其中垂直距离小于0.5米的点云数据为最终地面点。
[0024]优选的，步骤二中去除初始表面模型的粗差点是将所有点云数据的高程值取负值，然后使用高差阈值为5米、结构元素窗口半径为I米的形态学滤波来标记低位粗差点为非地面点。
[0025]优选的，步骤三中按尺寸大小为(20mX 20m)~(50mX50m)对滤波表面模型进行分块。
[0026]优选的，步骤四中最小结构元素窗口半径为I米，最大结构元素窗口半径为5~30米。
[0027]优选的，对于分块时不能被整分的滤波表面模型，被分块后剩余的格网，采用如下方法进行滤波处理:
[0028]从最小结构元素窗口尺寸开始执行以下操作，之后结构元素窗口的半径每增加一米，执行一次以下操作直至最大结构元素窗口尺寸:[0029]分别对剩余格网中所有点云数据进行开运算，得到各点云数据开运算后的高程值，
[0030]分别计算各点云数据的高程值与其自身的开运算后的高程值的差值，差值大于e的点云数据标记为非地面点，且将各点云数据的高程值更新为开运算后的高程值，其中:e=0.3XwXc，w为当前结构元素窗口的半径，c为格网边长，且I米≤c≤5米，滤波处理后未标记的点云数据为地面点。
[0031]与现有技术相比，本发明具有以下优点:
[0032]第一，本发明在对点云数据进行形态 学滤波之前，结合先验知识判断实验区域的地貌特征，对规则化DSM按一定的尺寸进行了分块。由于实验区域一般地形复杂多变，将整个区域分块能够针对分块中的相似地形进行滤波，以便有效的滤去地物点保留地面点。克服了现有技术中形态学滤波不能很好对复杂场景滤波的缺点，使得使用本发明中的方法可以获得较好的滤波结果。
[0033]第二，本发明在对点云数据进行形态学滤波处理时，使用了各个分块区域预测的地形参数来确定高差阈值。由于标准差代表了地形的起伏情况，所以以标准差对地形参数进行预测从而得到高差阈值可以自适应地对数据滤波。克服了现有技术中的地形参数需要人工设定并为恒定，自适应不强的问题，使得使用本发明中的方法可以根据地形自适应地对点云数据进行滤波。 [0034]本发明应用于点云数据滤波处理中，能够有效去除地物点并保留地面点，尤其能有效减小第II类误差，同时维持较低的第I类误差和总误差。
【专利附图】
【附图说明】
[0035]图1为实施例1中对sample 12滤波的仿真图；
[0036]图2为实施例2中对samp I e31滤波的仿真图；
[0037]图3为实施例3中对sample53滤波的仿真图；
[0038]图4为实施例4中对sample54滤波的仿真图。
【具体实施方式】
[0039]本发明的构建初始表面模型可采用本领域常规的方法，具体包括对点云数据格网化:根据点云数据横坐标X和纵坐标y的最大值和最小值确定点云数据的Xy平面范围，并选择合适的格网尺寸构成规则格网，使得每个格网中至少包含一个点云数据，合适的格网尺寸一般为(ImXlm)~(5mX5m)，且格网为正方形；如果在一个格网中没有点云数据，使用插值来得到缺失的数据点；如果一个格网中有多个点云数据，则选取格网中最低的点作为格网的取值。
[0040]本发明的去除初始表面模型的粗差点的方法可采用本领域传统的方法，优选以下方法:粗差点分为高位粗差点和低位粗差点，高位粗差点在滤波处理过程中可以被滤掉，而低位粗差点不能，需对低位粗差点进行单独处理，处理方法是将所有点云数据的高程值取负值，然后使用高差阈值为5米、结构元素窗口半径为I米的形态学滤波来标记低位粗差点为非地面点。
[0041]本发明中对滤波表面模型进行分块时，通过观察实验区域的地貌特征，对滤波表面模型按适合的尺寸进行分块，尽量保证块内的地形变化有相似性，并且分块面积不能小于最大建筑物面积以避免某个分块落入最大建筑物内，分块时的合适尺寸(20mX 20m)~(50mX 50m)。
[0042]本发明的滤波使用结构元素窗口尺寸线性增长的方式来迭代地进行滤波，初始结构元素窗口半径为lm，最大结构元素窗口尺寸大于实验区域中最大建筑物的尺寸,一般最大结构元素窗口半径为大于5~30m。
[0043]本发明所述的点云数据的高程值是指该点云数据的Z坐标值，对于点云数据来说，其高程值为Z值。各点云数据的初始高程值为其自身的Z坐标值。在滤波处理过程中各点云数据的最初高程值为其自身的Z坐标值，之后，每进行一次开运算，其高程值便更
新一次。
[0044]本发明所述的插值为文献参见文献:Thomas J.Pingel, Keith C.Clarke, WilliamA.McBride.An improved simple morphological filter for the terrainclassification of airborne LIDAR data[J].1SPRS Journal of Photogrammetry andRemote Sensing, 2013(77):21-30.中公开的 spring-metaphor inpainting technique。
[0045]本发明所述点云数据开运算后的高程值计算方法是首先对点云数据开运算指先对当前分块中所有点云数据进行腐蚀运算，再对分块中所有点云数据进行膨胀运算。具体包括:
[0046]首先对当前分块中所有点云进行腐蚀运算:
[0047]遍历每个点云数据，
[0048]设当前点云数据为rmn，
[0049]
【权利要求】
1.一种基于区域预测的LiDAR点云数据形态学滤波方法，其特征在于，方法包括以下步骤: 步骤一，构建滤波初始表面模型； 步骤二，去除初始表面模型的粗差点，得到滤波表面模型； 步骤三，对滤波表面模型进行分块: 按适当尺寸对滤波表面模型进行分块，其中分块n记为Rn，n = 1,2, 3，…，N，N为分块的总数；该分块Rn中包含M个点云数据，其中点云数据m记为rm，m =1, 2, 3，…，M ;步骤四，滤波处理: 对各分块分别进行滤波处理，其中对分块Rn，n = 1,2, 3，…，N的处理过程如下: 从最小结构元素窗口半径开始执行以下操作，之后结构元素窗口的半径每增加一米，执行一次以下操作直至最大结构元素窗口尺寸: 分别对分块Rn中所有点云数据进行开运算，得到各点云数据开运算后的高程值， 分别计算各点云数据的高程值与其自身的开运算后的高程值的差值，差值大于％的点云数据标记为非地面点，且将各点云数据的高程值更新为开运算后的高程值，其中:en =SnXwXc,w为当前结构元素窗口的半径，c为格网边长,且1米≤c≤5米，sn = aX σ n+b,σ n为分块Rn中所有点云数据高程值的标准差，0.005 ≤ a ≤ 0.02，0.005 ≤ b ≤0.02 ；滤波处理后未标记的点云数据为地面点。
2.如权利要求1所述的基于区域预测的LiDAR点云数据形态学滤波方法，其特征在于，方法还包括: 将步骤四确定的地面点插值得到数字高程模型， 接着分别确定所有点云数据的初始高程值与数字高程模型之间的垂直距离，其中垂直距离小于0.5米的点云数据为最终地面点。
3.如权利要求1所述的基于区域预测的LiDAR点云数据形态学滤波方法，其特征在于，步骤二中去除初始表面模型的粗差点是将所有点云数据的高程值取负值，然后使用高差阈值为5米、结构元素窗口半径为I米的形态学滤波来标记低位粗差点为非地面点。
4.如权利要求1所述的基于区域预测的LiDAR点云数据形态学滤波方法，其特征在于，步骤三中按尺寸大小为(20mX 20m)~(50mX50m)对滤波表面模型进行分块。
5.如权利要求1所述的基于区域预测的LiDAR点云数据形态学滤波方法，其特征在于，步骤四中最小结构元素窗口半径为I米，最大结构元素窗口半径为5~30米。
6.如权利要求1所述的基于区域预测的LiDAR点云数据形态学滤波方法，其特征在于，对于分块时不能被整分的滤波表面模型，被分块后剩余的格网，采用如下方法进行滤波处理: 从最小结构元素窗口尺寸开始执行以下操作，之后结构元素窗口的半径每增加一米，执行一次以下操作直至最大结构元素窗口尺寸: 分别对剩余格网中所有点云数据进行开运算，得到各点云数据开运算后的高程值， 分别计算各点云数据的高程值与其自身的开运算后的高程值的差值，差值大于e的点云数据标记为非地面点，且将各点云数据的高程值更新为开运算后的高程值，其中:e =.0.3 XwX c，w为当前结构元素窗口的半径，c为格网边长，且1米≤c≤5米，滤波处理后未标记的点云数据为地面点。
【文档编号】G06T5/00GK103745436SQ201310732382
【公开日】2014年4月23日 申请日期:2013年12月23日 优先权日:2013年12月23日
【发明者】苗启广, 宋建锋, 郭雪, 许鹏飞, 陈为胜, 宣贺君, 刘如意, 张萌 申请人:西安电子科技大学