矿区沉陷DEM 多重滤波方法研究

姚顽强 ，蒙延斌 ，郑俊良 ，薛志强

（1.西安科技大学 测绘科学与技术学院，陕西 西安 710054；2.陕西彬长孟村矿业有限公司，陕西 咸阳 713602）

煤炭资源在我国的能源中占据着很重要的位置，有效地推动了国民经济的快速发展[1-2]。但是随着煤炭开采强度的提高导致地表大范围塌陷，使地表产生移动和变形[3-5]，给矿区环境带来严重的影响[6]。在实地监测中，现有的水准测量、GPSRTK 测量及InSAR 等遥感监测方法均具有一定的局限性。前者不仅工作量大[7]，而且仅通过计算各项地表岩移参数与角量信息很难获取到整个移动盆地的形变信息；而后者InSAR 作为1 种全新的对地观测技术，虽然在矿区开采沉陷得到广泛应用[8-10]，但是只能测取地表的微小形变信息，存在大梯度变形区域失相干的缺点，使得InSAR 技术在沉陷监测中的应用受到了限制。

近些年来，无人机遥感技术发展迅速，在矿区监测中取得了一定成效，并在多方面得到应用[11-15]。无人机LiDAR 作为1 种激光遥感技术，具有获取数据速度快，飞行时受地形因素小的特点，可以获取到矿区高精度、高分辨率的三维点云数据，并进行点云滤波和插值生成DEM，并将多期DEM叠加得到沉陷DEM。但是基于点云滤波[16-18]的局限性，仍存在噪声，从而使沉陷DEM 也会包含噪声和误差。有学者通过不同的去噪算法来提高沉陷模型的精度，但是目前还是存在问题：①无人机所获取的数据精度不能确定；②利用现在主流的点云滤波和DEM 插值算法计算结果难以验证，影响后续DEM 的建模精度[19-20]；③没有根据整个移动盆地的误差分布特征进行算法的选择，去噪效果较差，不能满足沉陷监测的需求。

综上，以凉水井煤矿某工作面为研究区域，研究沉陷DEM 的去噪问题；利用无人机采集的2期点云数据进行滤波插值生成DEM，叠加后得到沉陷DEM；根据沉陷盆地的误差分布特征选择多重滤波组合去除噪声，并通过实测数据对比去噪效果，得到最优组合去噪方案，为沉陷监测提供了好的技术支撑。

1 研究区概况

研究区位于陕西省神木市的凉水井煤矿，在陕北黄土高原北部，毛乌素沙漠的南边，属于丘陵区。地理坐标为：东经111°14′22″～111°21′24″，北纬38°47′29″～38°53′24″。工作面内没有河流存在，地貌主要是以风积沙地貌为主，大部分区域为风积沙层所覆盖。工作面地表的其他区域无村庄、耕地、建筑物分布。地势由东北向东南方向逐渐降低，地表整体西部要高于东部，倾角小于1°。431303 工作面位于井田的最南端，是431 盘区的首个回采工作面，工作面为4-3煤层，4-2煤层是上覆岩层，在前期已开采完毕。煤层相对位置情况如图1，上部为盘区的南部采空积水区，西部为煤矿辅运大巷，北部和东部主要是实煤区。矿区工作面的回采区域内煤层属于较稳定煤层，煤层平均采深137.92 m。

图1 凉水井矿双煤层相对位置示意图Fig.1 Diagram of the relative position of the double coal seam

实测观测线整体布设图如图2。

图2 实测观测线布设图Fig.2 Layout of actual observation lines

研究区工作面走向为超充分采动，走向观测线布设在工作面的中心线上。综合考虑工作面地质结构情况，工作面两端走向线覆盖整个工作面，并向两端外侧各延伸150 m，布设了走向观测线、倾向观测线以及道路观测线，其中红色表示地面观测点。

实验所采用的机载LiDAR 系统型号为SZT-250，扫描的测程可以达到250 m，具有高精度的导航定位系统和激光扫描系统。数据采集时通过移动测量系统操控软件来获取点云的姿态和导航定位信息。激光扫描仪的测距精度为15 mm，视场角为360°，扫描速度在10～100 线/s，在进行多目标探测时，每束激光可以达到5 次回波，每期无人机航飞都采集出2 组数据，这样可以保证所采集的点云具有足够的密度，点云的质量也可以得到很好的保证。所扫描出的每2 期点云数据的平均密度为60～700 脉冲/m2，经过检验后，点云的质量达到标准，可进行点云数据的处理工作。

2 数据获取与误差分析

2.1 沉陷DEM 获取

在本次数据采集中,无人机LiDAR 对凉水井矿区工作面分别在2020 年8 月15 日（回采位置1 024.6 m）、2020 年11 月7 日（回采位置1 384.6 m）和2021 年5 月19 日（回采完毕）一共进行了3 次观测,无人机飞行时工作面开采位置如图3。

图3 无人机飞行时工作面开采位置Fig.3 Working face mining position at the time of the drone flight

点云处理及成果分析流程如图4。通过对期不同的点云数据进行滤波、插值,然后利用第1 期和第3 期的DEM 数据进行相减，得到累积的沉陷DEM。并视为这个DEM 为原始的沉陷DEM，由于利用现有的滤波算法不能准确的分离出地面点和非地面点，并且对点云进行插值的时候,还会存在插值误差，所以使得初始沉陷模型也会包含误差，对后期的实验处理结果会有很大的影响。将沉陷DEM 利用ArcScene 进行场景属性设置，基于范围进行计算，将沉陷DEM 高程方向拉伸40 倍，以突出沉陷模型的误差特征。沉陷DEM三维模型如图5。

图5 去噪前沉陷DEM 三维视图Fig.5 3D view of DEM before noise removal

2.2 误差分析

通过对原始点云的预处理，利用点云插值和滤波生成的DEM 仍然存在误差，主要来源有3 个方面：

1）滤波算法时未除去的非地面点误差。对点云进行滤波时，不同的滤波算法都会有一定的限制性，不能完全剔除一些误差对于模型的影响。在不同的季节，植被的生长，从而使一些植被点不重叠，导致沉陷DEM 会出现噪声。

2）不同点云插值生成DEM 引起的模型误差。每一期所采集的点云数据都会有所差异，例如点云的分布位置，点云的厚度和宽度，在内插生成DEM 时,会对精度造成一定的损失。采集点云数据的区域情况，当地可能会出现地形起伏,坡度和坡向变化的情况。这样会使得沉陷DEM 有很大的模型误差,要想减小这种误差,可以利用最优的点云插值算法或者去噪方法,同时也可以适当增加点云的宽度和厚度，让点云的分布更均匀,方便数据的处理。

3）无人机LiDAR 外业扫描引起的误差。无人机在外业扫描时，会受到外界环境和人为因素的影响。例如无人机的飞行姿态时刻都在变化，以及无人机在扫描地面时会和地面有1 个倾角范围，角度过大或过小会陷入盲区，导致有些地物或者植被扫描不到，点云也会有缺失和分布稀疏的情况。不同的航线之间采集数据具有重叠度，在没有误差的情况下,2 个航带所采集的数据是相互交融到一起的。这样在进行两期DEM 影像相减的时候也会出现误差。

然而减小这种误差的途径有很多，操作人员可以掌控无人机和地面的飞行高度，通过增加扫描范围获取更多点云的数量，利用TerraSolid 软件可以将范围线导入点云中，看范围线是否被包含在内，检查点云的质量。可以利用PCL 点云库对点云进行配准，根据源点云点数，目标点云点数进行重采样，最后输出变换矩阵，得到最终的配准点云数量，减小系统性误差和产生的一些随机性误差。

根据所分析的几种误差可知,不同的误差类型都可以影响沉陷DEM 的精度,所以为了后续不再对沉陷模型产生误差，直接对沉陷DEM 进行去噪。

3 沉陷DEM 滤波及实验分析

3.1 经典滤波和多重滤波原理

在图像中，由于成像系统、传输介质的不完善，图像在形成过程中往往会受到多种噪声的污染，扰乱图像的可观测的信息。沉陷DEM 与图像类似，采用图像滤波方法可实现沉陷DEM 中噪声的去除，细节特征更加明显。沉陷DEM 去噪流程图如图6。

图6 沉陷DEM 去噪流程Fig.6 Subsidence DEM denoising process

经典滤波目前常用的方法有高斯滤波、维纳滤波、中值滤波、均值滤波。4 种滤波都有不同的特性。均值滤波对高斯噪声有良好的去噪能力，但是在消除噪声的同时也会对图像的高频细节成分造成破坏和损失。高斯滤波根据图像的边缘特征，既含有低频分量又含有高频分量，平滑时不会被不需要的高频信号污染，保留了大部分信号。维纳滤波对高斯噪声有明显抑制作用，但是容易失去边缘信息。中值滤波让图像周围的像素值接近真实值，消耗孤立的噪声点，能保持良好的边缘特性，但是会失真过多。

多重滤波是根据误差的分布情况和算法特点选择不同的组合算法，以及根据沉陷DEM 的误差分布特性，并结合地表沉陷盆地本身的特征分布进行算法的选择叠加，本研究主要是选择了中值滤波+维纳滤波，维纳滤波+中值滤波+高斯滤波，均值滤波+维纳滤波，中值滤波+均值滤波这4 种组合方式，并对算法进行叠加，对地表沉陷未除去的噪声误差继续实验验证。

3.2 滤波方法对比

在MATLAB 的工作环境下,为了验证几种滤波去除噪声效果的好坏，从主观和客观的角度综合进行分析评价。分别使用中值滤波、维纳滤波、均值滤波、高斯滤波进行滤波处理。对比分析4种滤波的去噪效果；得到的中值滤波算法、维纳滤波算法、均值滤波算法、高斯滤波算法的PSNR 值分别为27.526 8、28.439 1、28.614 2、26.992 7。

PSNR 值是1 种图像客观的评价指标，基于误差敏感的图像质量评价，其PSNR 的计算公式为：

式中：PSNR 为图像的峰值信噪比；MSE 为原图像与处理图像之间均方误差。

式中：I、K分别为原始影像和处理后影像；i、j为像素值的行和列；m、n为表示大小的尺寸值。

在几种滤波方法中，对于含有噪声的沉陷DEM，PSNR 值有所不同，对于突变噪声，高斯和维纳滤波可以取得更好地去除效果，但是各种滤波算法对噪声伤害都很大，在尽量保留图像细节特征的条件下对噪声进行去除。而每种滤波特性有所不同，处理不好会使图像质量下降，失真，PSNR 值越大，表示图像的失真越少。

通过4 组沉陷DEM 去噪实验，得到的去噪后沉陷DEM 三维视图如图7。由图7 可以看出：去噪后的沉陷DEM 相对于初始的沉陷DEM 去噪效果更好,并且在沉陷DEM 中的非沉陷区域中服从高斯分布的突变噪声信号被去除了,保留了下沉盆地的一些真实的细节特征。也在一定程度上反映出了矿区沉陷模型的真实情况。

图7 去噪后沉陷DEM 三维视图Fig.7 3D view of DEM subsidence after noise removal

3.3 与实测下沉量对比

在开采沉陷期间,研究区可采煤层包括4-2煤层和4-3煤层，工作面实际煤层采高与观测下沉结果之间不符合开采沉陷的客观规律，是由于研究区多煤层重复采动异常造成的。4-2煤层虽然在开采完毕后经历长时间塌陷，但仍可能存在一定的裂隙和空洞。上覆岩层在短期内很难达到稳定状态。4-3煤层的开采对4-2煤层上覆岩层造成进一步扰动，破坏进一步加剧，下沉量更大。为了更直观地反映本研究不同方法进行沉陷DEM 去噪的效果，基于处理前、后的沉陷DEM 沿下沉盆地走向线提取下沉数据，以横坐标表示控制点点号,纵坐标表示下沉量，绘制下沉量曲线图。走向观测线不同滤波去噪下沉情况如图8。

图8 走向观测线不同滤波去噪下沉情况Fig.8 De-noising subsidence of strike observation line with different filters

在研究区选取54 个观测点，通过与实测数据对比可以看出：从控制点A173 到A182 的曲线走势情况来看,在这段区间地表有明显的下沉,形成1 个下沉盆地。综上，沉陷模型在去噪前的误差基础上，地表下沉盆地在边缘的部分相对于来说较平缓，在矿区工作面的中心区域下沉具有明显的突变特征。总体来预判，高斯滤波和维纳滤波的去噪效果相对较好。

3.4 去噪后误差分析和精度验证

通过上述下沉量曲线可以看出，高斯滤波和维纳滤波有着很好地去噪效果。这说明沉陷DEM中的一些误差服从高斯分布，去除了沉陷边缘的突变噪声。而均值滤波和中值滤波的效果却一般。分析得出，根据沉陷DEM 的误差特性，后两者的去噪范围并没有将此误差包含在内，不符合地表沉陷的分布特征，从而使得去除噪声效果较差。所以需要对4 种滤波算法组合去除噪声，将沉陷模型精度进一步提高。滤波后的直方图如图9。

图9 多重滤波直方图Fig.9 Multi-filter histograms

因为高斯噪声的范围会遍布到所有的灰度级，以及其他噪声分布也遵循一定的规律，所以为了达到更好地去除噪声效果，针对不同噪声采用不用的滤波器。相比而言，维纳滤波+中值滤波的对于高斯噪声的效果最好，维纳滤波对于服从高斯分布的突变误差更好，而其他几种方法不但消除了突变误差，还消除了边缘误差，非地面点误差，以及非开采沉陷因素所引起的误差。使得沉陷DEM的边缘信息更好，沉陷模型的精度相对于4 种独立去噪方法的精度提高了很多，去噪效果的好坏的顺序为维纳滤波+中值滤波＞维纳滤波+中值滤波+高斯滤波＞均值滤波+维纳滤波＞中值滤波+均值滤波。

对多重滤波的去噪结果进行验证与分析，走向观测线多重滤波下沉情况如图10，经典滤波去噪误差指标见表1，多重滤波去噪误差指标见表2，不同滤波算法误差分布如图11。

表1 经典滤波去噪误差指标Table 1 Classical filter denoising error indexes

表2 多重滤波去噪误差指标Table 2 Multi-filter denoising error indicators

图10 走向观测线多重滤波下沉情况Fig.10 Multiple filter sinking of strike observation line

图11 不同滤波算法误差分布图Fig.11 Error distribution of different filtering algorithms

实验表明，根据沉陷DEM 的误差分布特征，所组合起来的几种滤波去噪方法有很大的成效，相比于上述的其他4 种单独滤波方法不仅去除了沉陷DEM 服从高斯分布的突变误差，还有非地面点误差，植被覆盖度所引起的重叠误差，以及边缘误差，使边缘轮廓更加清晰。沉陷区域的下沉量也明显减少了很多。沉陷DEM 的精度得到了改善，表1 和表2 中的误差指标也符合要求，均方根误差和平均绝对误差都有所减小，也体现出了去噪前后误差的分布情况，与预想是一致的，并且在上述的去噪实验中，维纳+中值滤波去噪效果是最好的，也证实了本研究对于沉陷DEM 去噪的效果具有较高的可靠性。

4 结语

利用无人机LiDAR 技术进行矿区地表沉陷监测时，可以获取到高精度的点云数据，但是目前主流的点云插值与滤波算法生成沉陷DEM 时会有精度不足的现象。其中，非地面点噪声、点云插值误差以及外业扫描所引起的点云密度不足等误差都会使沉陷DEM 有噪声出现。基于矿区沉陷DEM 的建模精度和误差分布情况，根据几种单独滤波和多重滤波结果表明，多重滤波不仅去除了服从高斯分布的突变误差，还去除了非沉陷区域中的边缘误差，非地面点误差等误差。其中维纳+中值滤波去噪效果最佳，下沉盆地的细节突出更加明显，说明基于沉陷DEM 误差所采取的去噪方法具有一定的合理性和可靠性。但是在本研究中，所得到的实验结果是在榆神矿区得到的，地理环境相对较复杂，还要其他未知误差的来源，所以利用无人机LiDAR 对开采沉陷还需要进一步研究。