0.   引言

利用卫星、无人机等拍摄的矿区遥感图像广泛用于矿区生产管理与监测，如道路规划、无人驾驶、塌陷区扰动识别、矿山植被修复监测、土壤侵蚀估算、地质滑坡监测与识别、污染物识别等^[1-5]。受设备稳定性、电子元件安装、成像及计算误差、信号传输等因素影响，矿区遥感图像会出现噪声污染、模糊化等情况，导致后续难以获得较好的视觉处理效果。因此，去噪成为遥感图像得以有效应用的重要预处理步骤。

现有的遥感图像去噪方法大致可分为基于统计方法、基于域变换方法、基于学习方法3类。基于统计方法根据含噪图像本身的统计特征完成去噪，主要包括MF (Mean Filtering，均值滤波)、NLM (Non−Local Mean，非局部均值)滤波、BF (Bilateral Filter，双边滤波)等。文献[6]采用BF进行结构特征保持，实现了矿区遥感图像初步去噪。基于域变换方法是将图像数据转换到频域^[7-9]，通过阈值函数得到合适的图像分割及处理频率系数，然后进行逆变换重构，得到保留了大部分能量的低频数据对应的去噪图像，主要包括 WT (Wavelet Transform，小波变换)、FT (Fourier Transform，傅里叶变换)等。文献[10]采用提升小波与双参数阈值函数完成矿区遥感图像去噪，结合边缘算子重构了清晰的图像。基于学习方法主要包括2类：① 生成式去噪模型，如字典学习及稀疏表示等，通过学习含噪图像的过完备字典，将图像表示为字典原子稀疏的线性组合，达到去噪目的；② 判别式去噪模型，如深度RNN (Recursive Neural Network，递归神经网络)、CNN (Convolutional Neural Network，卷积神经网络)、GAN (Generative Adversarial Network，生成对抗网络)等，通过学习大量的噪声图像，形成数据映射式的去噪模型^[11-13]。文献[14]提出一种基于聚类的组稀疏字典学习方法，并将其用于多光谱遥感图像去噪。文献[15]采用平移不变性的K−SVD (K−Singular Value Decomposition,K奇异值分解)完成字典学习，结合稀疏表示理论，将去噪过程转换为一个L1范数非凸函数优化问题，实现遥感图像去噪。文献[16]针对残差网络提出空洞卷积方法，并扩大神经网络的感知野，提高了遥感图像去噪性能。文献[17]创新性地结合传统WT和GAN进行遥感图像去噪和超分辨率重建，克服了在非GAN下细节趋于平滑的问题。

上述方法均具有较好的遥感图像去噪效果，但普遍存在细节过度平滑、纹理保持不足等缺点，而纹理细节是遥感图像应用中的重要信息。本文基于引导滤波良好的边缘保持特性，提出了迭代引导滤波方法，可增强遥感图像边缘特征提取效果；将迭代引导滤波与传统去噪方法结合，可提高图像去噪效果；采用迭代引导滤波与BM3D(Block Matching 3D，三维块匹配)融合方法对矿区遥感图像进行去噪处理，结果表明融合方法可在提升去噪性能的同时，很好地保持细节特征。

1.   基于引导滤波的图像边缘保持方法

1.1   引导滤波

引导滤波用1张引导图像对输出图像中像素进行线性表达^[18]，即

式中：$ {w}_{k} $为以引导图像中任意像素点k为中心的窗口；$ {q}_{i} $为输出图像中窗口$ {w}_{k} $内像素点i的值；$ {a}_{k} $，$ {b}_{k} $为线性表达参数；$ {I}_{i} $为引导图像中窗口$ {w}_{k} $内像素点i的值。

当$ {I}_{i} $为待优化输入图像时，输出图像与输入图像之间的梯度存在倍数关系，因此输出图像能够很好地保留输入图像的边缘部分。设待优化输入图像$ {p}_{i} $含有加性噪声$ {n}_{i} $，则其与输出图像的关系为

在窗口$ {w}_{k} $中关于参数$ {a}_{k} $，$ {b}_{k} $的引导滤波最小二乘损失函数为

$式中 \epsilon$为输入图像一阶边缘参数约束项的权重。

式(3)中，ϵa_k²为输入图像的一阶边缘参数约束项。根据式(3)得

式中：$ m $为窗口$ {w}_{k} $内的像素个数；$ {\mu }_{k} $，$ {{\sigma }_{k}}^{2} $分别为输入图像在窗口$ {w}_{k} $中像素值的均值及方差；$ {\bar p}_{k} $为输入图像在窗口$ {w}_{k} $内的像素均值。

对输出图像中包含任意像素点i的窗口w_i的线性表达参数取均值，得到式(1)的输出结果：

1.2   迭代引导滤波

在引导滤波基础上，引入迭代引导滤波，其主要思想是将引导滤波结果与输入图像作差，得到残差图像(主要为输入图像与引导滤波输出图像在边缘上的内容差值及噪声)，再对残差图像进行引导滤波，所得结果与最初的引导滤波结果相加，得到最终输出图像，并将该图像作为下一次引导滤波的输入，直至满足迭代条件。迭代引导滤波过程如图1所示。

图  1  迭代引导滤波过程

Figure  1.  Iterative guided filtering process

下载: 全尺寸图片 幻灯片

迭代过程中$ \epsilon $大小会影响$ {a}_{k} $，进而影响$ {\bar a}_{i} $。图像质量随着迭代次数增加不断提高，为了得到更好的边缘保持效果，应增大$ {a}_{k}{I}_{i} $在损失函数中的权重，因此对$ \epsilon $进行更新：

式中：$\epsilon_{j} $为$\epsilon $在第$j( {j = 1,2}, $…$, n)$ 次迭代中的值；$ n $为迭代次数。

迭代引导滤波具体步骤如下。

(1) 初始化约束参数$ \epsilon $、迭代次数n及其他结束条件(如图像块梯度变化量、PSNR(Peak Sgnal to Noise. Ratio，峰值信噪比)变化值、残差图像像素二范数等)。

(2) 输入残差图像。

(3) 根据式(3)—式(8)进行引导滤波，并将滤波结果与输入残差图像作差，得到残差结果。

(4) 对残差结果进行引导滤波并输出结果。

(5) 将步骤(4)中输出结果与输入残差图像相加。

(6) 判断是否满足迭代次数或其他迭代要求，若是则停止迭代，否则根据式(9)更新参数$ \epsilon $，并转至步骤(2)继续迭代。

选取1张典型的含高斯噪声图像进行迭代引导滤波，设置窗口半径为1，$ \epsilon $=0.64，n=7，滤波结果如图2所示。可看出经7次迭代引导滤波后，取得了较好的去噪效果。对比首次引导滤波后的残差图像及末次迭代引导滤波后的残差图像(为便于观察，对图像进行了像素值翻转处理)，可看出经首次引导滤波后，残差图像中仍存在大量边缘数据，导致滤波结果较为平滑，而末次迭代引导滤波后的残差图像中边缘数据较少，多为离散的噪声数据。以最终输出图像的PSNR为衡量标准，不同迭代次数下PSNR如图3所示。与单次引导滤波相比，迭代引导滤波处理后的图像PSNR提升了2.3 dB。

图  2  噪声图像迭代引导滤波结果

Figure  2.  Iterative guided filtering results of noise image

下载: 全尺寸图片 幻灯片

图  3  迭代引导滤波处理后图像PSNR

Figure  3.  PSNR of the image processed by iterative guided filtering

下载: 全尺寸图片 幻灯片

上述研究表明，迭代引导滤波对于具有稀疏特性的残差数据具有很好的边缘保持效果。本文将迭代引导滤波分别与K−SVD字典学习、小波软阈值、NLM滤波、BM3D滤波相结合，用于提升上述方法的去噪性能。在计算机参数为Intel I5−7200U，4 GB RAM，实验平台为Matlab R2016条件下，以图2(a)所示噪声图像为例，采用上述方法进行去噪处理，将去噪结果与输入图像作差，得到残差数据，对其进行迭代引导滤波。以PSNR，SSIM(Structural Similarity，结构相似性)为指标，迭代引导滤波对不同方法去噪性能的提升结果见表1。

表  1  迭代引导滤波对于典型去噪方法的提升结果

Table  1.  Improvement results of iterative guided filtering to typical denosing methods

指标	K−SVD字典学习	小波软阈值	NLM滤波	BM3D滤波
PSNR增大值/dB	0.6	1.9	3.1	3.0
SSIM	0.1	0.2	0.4	0.4
运算时间增加值/s	0.5	0.7	1.0	1.2

下载: 导出CSV 

| 显示表格

从表1可看出：迭代引导滤波对4种方法的去噪性能均有所提升，对BM3D滤波、NLM滤波的提升效果明显，对K−SVD字典学习提升效果稍差，主要原因是K−SVD字典学习根据系数稀疏性实现去噪，而残差图像中的边缘数据不够突出；融合迭代引导滤波后，4种方法的运算时间增加值均较小。

2.   迭代引导滤波与BM3D滤波融合方法

鉴于迭代引导滤波对BM3D滤波的去噪性能有较好的提升效果，本文采用迭代引导滤波与BM3D滤波融合，对矿区遥感图像进行去噪。BM3D滤波在NLM滤波基础上发展而来，结合了基于统计和基于域变换的图像去噪方法，具有信噪比高、视觉效果好、计算复杂度低等优点。

2.1   BM3D 滤波

BM3D滤波主要包括基础估计和最终估计^[19]。基础估计如下。

(1) 相似块聚类。在输入图像中进行$ {W}_{1} $尺寸网格划分，并在$ {W}_{2} $大小的空间内寻找相似块，将所有相似块组成一个三维块数据。

(2) 对三维块中第二维块数据进行小波变换，并对第三维数据进行Hadamard一维变换。针对变换结果，采用硬阈值函数进行第三维数据的一维反变换及第二维块数据的二维反变换。

(3) 将步骤(2)中反变换结果聚合为原位置的同尺寸图像，每个像素点值等于其所在相似块值的加权和。

最终估计步骤如下。

(1) 对含噪图像及基础估计处理后得到的图像分别进行相似块聚类，得到2组三维块数据。

(2) 对得到的2组三维块数据进行与基础估计中步骤(2)相同的变换处理，并对含噪图像三维块数据的变换系数进行维纳滤波缩放。

(3) 对变换结果加权聚合到原位置，得到输出图像。

2.2   融合方法去噪性能

某矿区遥感图像大小为1 206×765，对其进行方差为$ 20 $的高斯处理，模拟成像及信号传输噪声。分别采用小波软阈值、K−SVD字典学习、非局部相似性K−SVD字典学习、BM3D滤波、迭代引导滤波与BM3D滤波融合方法对含噪图像进行处理。其中，小波软阈值中小波基为经验小波；K−SVD字典学习区域划分大小为$ \sqrt{M}\times \sqrt{N} $($ M $,$ N $分别为图像长度和宽度)，重建方法为OMP(Orthogonal Matching Pursuit，正交匹配追踪)，重建原子个数为20；迭代引导滤波与BM3D滤波融合方法中窗口半径取2，$ \epsilon $=0.64，n=30。图像去噪结果如图4所示，性能指标见表2。可看出迭代引导滤波与BM3D滤波融合方法较其他方法恢复了更多的图像边缘特征，纹理更加丰富，且PSNR更大、SSIM更高。

图  4  遥感图像去噪结果

Figure  4.  Denoising results of remote sensing images

下载: 全尺寸图片 幻灯片

表  2  不同方法的去噪性能对比

Table  2.  Comparison of denosing performance of different methods

指标	小波软阈值	K−SVD 字典学习	非局部相似性 K−SVD字典学习	BM3D 滤波	融合方法
PSNR/dB	23.6	25.2	27.3	28.2	30.8
SSIM	0.76	0.72	0.75	0.86	0.92

下载: 导出CSV 

| 显示表格

3.   融合方法在矿区遥感图像去噪中的应用

以煤矸石场识别及矿区滑坡区域分割为例，验证迭代引导滤波与BM3D滤波融合方法的应用效果。某矿区煤矸石场遥感图像如图5所示。针对原始图像及采用融合方法去噪后的图像，采用文献[20]中方法进行煤矸石场识别，将识别结果与人工标注结果进行对比，可看出采用融合方法去噪后区域识别准确率更高。贵州省毕节地区某矿区附近滑坡区域遥感图像如图6所示，采用Canny算子对经融合方法去噪前后的图像边缘进行识别，进而完成区域分割，结果表明去噪后的图像区域划分更接近实际滑坡区域。

图  5  煤矸石场识别应用效果

Figure  5.  Application effect of coal gangue yard identification

下载: 全尺寸图片 幻灯片

图  6  滑坡区域边缘识别应用效果

Figure  6.  Application effect of landslide area edge recognition

下载: 全尺寸图片 幻灯片

4.   结论

(1) 在引导滤波基础上，提出了迭代引导滤波。通过对残差信息进行引导映射，并迭代进行引导滤波及超参数收缩，在有效去除图像噪声的同时，更好地保持边缘数据，提升去噪效果。

(2) 将迭代引导滤波分别与传统的K−SVD字典学习、小波软阈值、NLM滤波、BM3D滤波相结合对图像进行去噪，结果表明输出图像的PSNR、SSIM均得到提升，NLM、BM3D滤波的去噪性能提升效果更明显。

(3) 将迭代引导滤波与BM3D滤波融合方法应用于矿区遥感图像去噪，结果表明融合方法较其他方法提升了图像PSNR及SSIM，并可获得更好的局部纹理视觉，用于矿区煤矸石场识别及滑坡区域边缘识别场景中取得了较好的效果。