自适应调整窗口维度的动态中值滤波算法

董林鹭，杨平先，陈明举，林国军

（四川轻化工大学自动化与信息工程学院，四川自贡643000）

由于互联网的快速发展，获取信息时，在传播过程中会受到各种干扰，使信号质量变差. 造成图像随机产生暗点和亮点，称为椒盐噪声[1-3]，该噪声严重影响了图像信息的表达，除去图像椒盐噪声的研究从未停止过[4-5]，学者们提出了许多滤波方法[6-9]，如 自 适 应 中 值 滤 波（AMF）[10]、加 权 中 值 滤 波（WMF）[11]等，利用图像内部数据，对噪声点进行滤波处理，其缺点是滤波处理时间较长，无法满足实际工程需要. 王志军[12]提出一种自适应中值梯度倒数加权的图像滤波算法.该算法自适应调节窗口大小检测出椒盐噪声并进行中值滤波处理，同时通过设定的阈值对高斯噪声进行梯度倒数加权平滑，以达到滤除混合噪声的目的. 兰霞[13]等人提出了一种迭代的中值滤波算法，该方法首先依据像素的强度值判断噪声点的位置，然后在循环迭代的处理框架内，对噪声像元进行逐步恢复，但是对复杂模型的噪声效果还不够理想. 这些滤波算法都有一个特点，即计算复杂度高且滤波时间长，不能满足实际工程应用中的大批量实时滤波要求. 自适应中值滤波（AMF）[10]算法首先选定大小为3×3 的滤波窗口，利用该窗口内像素中值分别与窗口内像素最大值、像素最小值做差，若前者之差＜0且后者之差＞0，则用窗口中间位置处像素值与窗口内像素最大值、像素最小值分别做差，若前者之差＜0且后者之差＞0，就判定该窗口中心处像素点为原图像素，直接出该像素值.否则增加窗口维度，继续重复以上步骤，直到窗口维度达到设定值，完成自适应中值滤波过程.

王拓[14]等人提出一种消除椒盐噪声的迭代自适应中值滤波算法，该算法主要是利用椒盐噪声的特性即噪声像素值为（0,255），先判断3×3窗口内噪声是否为百分之百，若不是，则选取该窗口内非噪声像素的中值替换窗口中心处噪声像素，若是，扩大窗口为5×5. 重复以上过程，若是百分之百被污染，则标记该窗口中值处像素坐标，进行迭代处理，完成滤波过程. 周华[15]提出基于动态窗口的自适应中值滤波算法，该算法先进行噪声检测，即利用噪声点与原像素点的差值关系，设立阈值，区分出椒盐噪声像素点与原图像像素点，判断出噪声点后，利用比较AM值与像素极大值Gmax与极小值Gmin的关系，确定是否调整窗口大小，达到动态窗口自适应中值滤波的目的.

本文利用传统中值滤波（MF）[16]能快速滤波的优点，对其不足之处，即需要手动调整滤波窗口和滤波过程造成图像细节模糊化的缺点进行改进，提出一种基于椒盐噪声自适应调整滤波窗口的算法，其思路是估计出图像的噪声密度，并在手动调整中值滤波窗口情况下进行大量实验，得出不同噪声密度下的最优滤波窗口，建立两者的函数关系，达到自适应选取在该密度下滤波效果最好的滤波窗口维度，再从滤波后的图像中分离出噪声修复像素和滤波前未被噪声污染的细节像素，将细节像素与噪声修复像素镶嵌融合，达到保护图像细节的作用.

1 经典中值滤波原理

中值滤波就是用一个含有奇数点的移动窗口，将窗口中心点的值用窗口内各点的中值代替，即：

式中：Sxy是以像素(x,y)为中心的含奇数个像素点的邻域窗口，{ f(s,t)}是窗口Sxy内所有像素灰度值按大小排序后的一维数据序列. 图1 是一幅行与列为5×5 的图像，对其采用滤波窗口分别为3×3、5×5中值滤波处理，达到经典中值滤波算法原理的目的.

图1 中值滤波模板

例如，对于像素(x,y)的一个3×3 的邻域内所有像素值(10,20,20,20,15,20,20,25,100)，对这些值排序后为(10,15,20,20,20,20,20,25,100)，那么，其中值是第5个值，而在一个5×5的邻域中，中值就是第13个值，等等. 可见该算法简单，滤波处理时间快，但容易造成图像细节模糊[17-18]、细节丢失和在不同噪声密度下需要手动调整窗口维度来获得最优输出结果等缺点.

2 本文算法

本文算法分为两部分：①估计噪椒盐噪声密度P，并根据噪声密度P%自适应调节最优中值滤波窗口对噪声图像进行滤波处理；②优化滤波处理后图像的细节.

2.1 估计椒盐噪声密度

分析椒盐噪声特点，将图像的像素归为两类，第一类是图像细节像素点记为fi,j，第二类是椒盐噪声污染像素点记为Zi,j，椒盐噪声污染图像数学表达式为：

利用经典中值滤波对噪声图像进行滤波预处理，目的是为了快速得到一幅近似原图的滤波图，其表达式为：

式中，f"i,j表示滤波处理后原细节像素，由于包含了误处理像素点，所以区别于fi,j且fi,j≈f"i,j.Z"i,j表示噪声修复区域像素.为了减少滤波误处理带来的误差影响，本文在提取噪声时加入阀值，经过大量实验后，阈值设为40 能得到最优结果. 利用噪声图像减去预处理图像. 当各像素点差值大于阈值时，判定该点为盐（像素值255）噪声，并令其值为1. 当各像素点差的绝对值小于0时，判定该点为椒（像素值0）噪声，为了在估计噪声密度时方便统计噪声像素点个数，令所有小于0 的像素点其值为1. 提取噪声像素Zi,j的过程满足表达式（4）：

利用MATLAB 统计出图像中像素Zi,j的个数，即寻找差值图像中像素值为1 的个数（记为S1）. 再利用图像像素维度的高M、宽N，求出图像的总像素点个数S=M×N，则椒盐噪声估计密度P 表达式满足：

2.2 本文算法对滤波后图像细节的保护

由经典中值滤波算法（MF）滤波原理可知，在滤波过程中是遍历整幅图像，对每个像素点取中值，图像细节部分同时被模糊化. 针对这一问题，在估计出噪声密度后利用本文算法对噪声图再次进行滤波，中值滤波窗口为最优输出效果时窗口大小，然后提取滤波后图像中噪声修复区域像素Z"i,j，表达式如（6）所示：

式（6）中设立了阀值15，目的是减少fi,j≈f"i,j的误差对Z"i,j的影响.fi,j由式（2）、（4）联合求得：

利用式（7）得到未被噪声污染的像素fi,j与式（6）得到的滤波后图像中噪声修复像素点Z"i,j镶嵌融合，保证了滤波后未被噪声污染区域像素的细节不会受到滤波算法的影响，表达式为：

Gi,j即为本文算法滤波后的最终输出结果.

为了更加直观地反映出本文算法对滤波后图像细节的保护，可以应用图例来具体介绍该过程. 首先对大小为256×256 的8 bit 标准MATLAB 实验图Baboon（图2a），加入20%密度的椒盐噪声（图2b），用窗口大小为5×5 的中值滤波（MF）对其滤波处理（图2c），整个过程如图2所示.

图2 中值滤波过程

由公式（4）提取噪声Zi,j过程如图3所示.

由公式（7）提取未被污染区域像素fi,j的过程如图4所示.

图3 噪声提取过程

图4 未被污染区提取过程

由公式（6）提取噪声修复区域像素Z"i,j过程如图5所示.

图5 噪声修复区域提取

由公式（8），将噪声修复区域Z"i,j与未被噪声污染区域fi,j镶嵌融合得到本文算法的输出结果，如6所示.

图6 算法输出过程

3 本文算法流程

首先对噪声图进行中值滤波预滤波处理，然后利用噪声图Ui,j与滤波处理后的图像ui,j做差与阀值比较提取出噪声Zi,j，由噪声图像Ui,j与噪声Zi,j做差取出未被污染区域fi,j，再利用未被污染区域像素点个数与图像总像素点个数估计出噪声密度P，根据噪声密度与中值滤波窗口的函数关系，自适应调整中值滤波窗口对噪声图像进行滤波处理，最后结合提取出来的未被污染区域fi,j与阀值提取出本文算法滤波后的噪声修复区域Z"i,j，将fi,j与Z"i,j镶嵌融合，得到一张完整的滤波处理图，图7为本文算法的输出结果.

图7 本文算法流程图

4 实验结果与分析

4.1 噪声密度与对应最优滤波窗口关系确定

本文选用了30 幅大小为256×256 的8 bit 标准MATLAB实验图片，加10%～90%密度的椒盐噪声，分别用窗口大小为W(3×3)、W(5×5)、W(7×7)、W(9×9)、W(11×11)、W(13×13)、W(15×15)的中值滤波（MF）对椒盐噪声图像进行滤波处理，并采用峰值信噪比对各中值滤波窗口在不同椒盐噪声下的效果做出客观评价，取30 幅图像的平均评价值（表1）来确定不同噪声密度下中值滤波算法的最优窗口大小，并建立函数关系.其中峰值信噪比(PSNR)原理如下：

式中，M 和N 分别为图像横向和纵向的像素个数；X(i,j)为原始图像的灰度值；Y(i,j)为经过滤波器处理后的图像灰度值；F 为图像中可能出现的最大灰度值（本文选用8 bit 灰度图像，则S=255）. 由式（9）可知，PSNR越大，滤波图像与原图越接近.

表1 各种中值滤波窗口对30幅图像滤波后的平均PSNR值的比较

分析表1 可见，噪声密度在10%、20%时窗口大小为W(3×3)的PSNR 值最高、当椒盐噪声密度在30%、40%时窗口大小为W(5×5)的PSNR 值最高，当椒盐噪声密度在50%、60%时窗口大小为W(7×7)的PSNR 值最高，当椒盐噪声密度在70%时窗口大小为W(11×11)的PSNR值最高，当椒盐噪声密度在80%时窗口大小为W(13×13)的PSNR值最高，当椒盐噪声密度在90%时窗口大小为W(15×15)的PSNR值最高，为了直观地观察不同椒盐噪声下中值滤波窗口对滤波效果的影响，将表1绘制成折线图如图8所示.

根据表1得到了不同密度下椒盐噪声滤波效果最好时所对应中值滤波窗口大小的关系，从而建立自适应调取窗口大小的函数关系如下：

图8 各种中值滤波窗口对30幅图像滤波后的平均PSNR值的比较

W(n×n)即为同椒盐噪声密度下自适应选取中值滤波窗口的大小.

4.2 本文算法的主观评价

为了证明本文对中值滤波算法细节的保护优势，选用4 幅标准MATLAB 实验图加入30%密度的椒 盐 噪 声，分 别 用 窗 口 大 小 为W(3×3)、W(5×5)、W(7×7)、W(9×9)、W(11×11)、W(13×13)、W(15×15)和本文的算法对其滤波处理，得到的滤波效果如图9所示.

图9 30%噪声密度下不同窗口中值滤波效果与本文算法对比

结合公式（10）与图9 可以看出，在噪声密度为30%时，中值滤波窗口滤波效果最好的是W(5×5)，同时本文算法是基于中值滤波在该噪声密度下滤波效果最好的窗口大小，但本文算法的滤波效果在细节保护上比中值滤波在该密度下的最优窗口保护细节的能力更强.

4.3 本文算法的客观评价

对一种滤波算法的评价常常以峰值信噪比、结构相似度的大小以及滤波处理时间t（单位为s）来对其进行评价，其中时间的，如果只看滤波效果而忽略滤波处理时间的长短，将无法满足实际工程的要求.

因此首先验证本文算法对噪声图像估计出的椒盐噪声密度与实际椒盐密度的误差，再采用客观的评价函数图像的峰值信噪比与图像的结构相似度来对本文算法自适应调用滤波库中算法进行客观验证.其中峰值信噪比的评价原理在前文有详细的介绍，对于结构相似度的图像评价原理，其表达式为：

式中，μx是x的平均值，μy是y的平均值，是x的方差，是y 的 方 差，σxy是x 和y 的 协 方 差. c1=(k1L)2、c2=(k2L)2是用来维持稳定的常数，其L 是像素值的动态范围，k1=0.01，k2=0.03. 由式（11）可得，在0 到1 之间，SSIM 值越接近1，两张图片结构越相似，当值为1时表示两图完全一致.

本文提出的方法运用MATLAB 2016 编程，以2.60 GHz 的英特尔CPU、运行内存为8 GB 的笔记本平台进行仿真，分别对265×265 大小8 bit 的灰度图lena、einstein、Baboon、Cameraman加10%～90%的椒盐噪声，再用窗口大小为W(3×3)、W(5×5)、W(7×7)、W(9×9)、W(11×11)、W(13×13)、W(15×15)的中值滤波算法和本文算法对图像进行滤波处理，得到客观评价值后进行比较，结果如表2 所示，其中Wz表示噪声图与原图的峰值信噪比，客观评价值是多次试验后取得的平均值.

从表2 可以看到，本文对应噪声的估计误差在3%内波动，为了满足函数的调用公式，本文采用了四舍五入的方法，使滤波窗口与噪声密度相匹配，达到快速调整滤波窗口的目的.为了直观地了解本文自适应选择滤波窗口与滤波处理中对图像细节的保护优势，将表2 中峰值信噪比和图像相似度绘制成折线图，结果如图10、11所示.从图中可以看到在不同噪声密度下本文算法客观评价值都是最优，同时随着噪声密度的增加，本文算法对于同一大小窗口的中值滤波的客观评价的提升是逐渐变弱的. 因为，噪声密度越大图像未被污染的像素就越少，所以对细节保护的比例就越小.在滤波时间上本文算法只是在10-3s级别增加变化，完全可以忽略这样的变化. 因此本文算法在实际工程中，对于大批量的椒盐噪滤波处理有着重要的意义.

图10 不同窗口中值滤波与本文算法的PSNR评价比较

图11 不同窗口中值滤波与本文算法的SSIM评价比较

5 结束语

中值滤波因为算法简单、滤波处理时间快的特点被广泛的应用于椒盐噪声的滤出，但该算法在不同的椒盐噪声密度下滤波窗口大小对滤波效果的影响非常大，为了得到最好的滤波效果常常需要手动调整滤波窗口大小，严重限制了该算法在工程上的应用. 针对该问题，本文提出的基于椒盐噪声密度自适应调整中值滤波（MF）窗口的算法，根据噪声的密度自适应调整中值滤波的窗口，使得滤波效果达到最理想的结果，而且对中值滤波本身的滤波时间只是在10-3s级别上略微增加，同时细节保护能力更强，使算法能够满足实际工程应用的要求.