数学物理学报  2018, Vol. 38 Issue (4): 823-832   PDF    
扩展功能
加入收藏夹
复制引文信息
加入引用管理器
Email Alert
RSS
本文作者相关文章
汪晓艳
丁义明
一种去除图像中Cauchy噪声的滤波算法
汪晓艳1, 丁义明2     
1. 武汉大学数学与统计学院 武汉 430072;
2. 武汉理工大学数学科学研究中心 武汉 430070
摘要:数字图像在采集、传输等过程中会产生各种噪声,噪声特征不同,处理方法也不同,如何尽可能恢复被强噪声干扰的图像是一个有意义的研究课题,因为传统的滤波算法在强噪声情况下,难以得到理想的结果.该文在中值滤波的基础上,结合局部区域内像素聚类的思想,提出一种对图像边缘进行修正的滤波算法,对三幅被Cauchy噪声干扰的图像处理的结果表明,该算法弥补了中值滤波在细节处理方面的不足,和其它方法相比,在滤除强噪声并保护图像细节和边缘方面有明显提高,是一种有效的去除强噪声的滤波算法.
关键词图像去噪    强噪声    中值滤波    K-均值聚类    
A Filtering Algorithm for Removing Cauchy Noise in Images
Wang Xiaoyan1, Ding Yiming2     
1. Department of Mathematics and Statistics, Wuhan University, Wuhan 430072;
2. Center for Mathematical Sciences, Wuhan University of Technology, Wuhan 430070
Abstract: Digital images may be disturbed by different noises during acquisition and transfer processes. Treatment methods are varied with the characteristics of noises. It is interesting to develop algorithms to recover the images from strong noises, because the traditional filtering algorithms do not work well in strong noises. A new filtering algorithm combined median filtering with local area pixel clustering is proposed to filter strong Cauchy noises from digital images. Based on the processing of three images with Cauchy noise, it is demonstrated that the algorithm can obtain more details of images than the median filtering. In comparison with other methods, our filtering algorithm is an effective strong noise removal filtering algorithm because it protects image details and the edges significantly.
Key words: Image denoising     Strong noise     Median filtering     K-means clustering    
1 引言

图像是人们从客观世界获取信息的重要来源之一, 是人类进行信息传播和交流的主要媒介.由于图像在生成、传输过程中具有的某些未知因素和不确定因素, 导致图像在包含丰富信息的同时, 还可能存在多种噪声如高斯噪声(主要由阻性元器件内部产生)、椒盐噪声(主要是图像切割引起的黑图像上的白点噪声或光电转换过程中产生的泊松噪声), 以及乘性噪声等, 具体的噪声种类取决于噪声的分布情况, 噪声的特征与噪声分布直接相关.一般情况下, 在对图像进行分析之前, 图像去噪是数字图像处理中的重要环节和步骤, 作为图像预处理步骤, 需要在尽可能去除噪声的同时, 尽量保留更多原图细节信息, 因为去噪结果直接影响到图像的后续分析等操作.

针对不同性质的噪声污染的图像应该采用不同的滤波方法, 一般说来, 图像滤波方法的选择取决于噪声与图像的关系.描述噪声的方法一般借用随机过程的描述, 即用其概率分布函数和概率密度分布函数, 如高斯噪声服从正态分布, 椒盐噪声则是由图像传感器、传输信道、解码处理等产生的黑白相间的亮暗点噪声, 是一种非连续分布.柯西分布是数学期望和方差不存在的连续型分布, 由于柯西分布的特殊性, 因此该文研究图像在加上服从标准柯西分布的噪声后, 如何较好地去除图像噪声.

20世纪70年代, Tukey提出的基于排序统计的中值滤波是当前使用最为广泛的非线性抑制噪声的方法之一.从中值滤波的去噪原理我们知道, 在图像噪声较弱时, 中值滤波去噪效果较好, 而当噪声较强时, 中值滤波受噪声干扰较大, 损失大量的细节信息, 造成图像边缘模糊及细节缺失.为了改进这些问题, 出现了多种基于中值滤波的改进算法研究, 如自适应中值滤波(Adaptive Median Filter, AMF)[1], 开关中值滤波(Switch Median Filter, SWF)[2], 极值中值滤波(Extremum Median Filter, EMF)[3].这些算法在实际图像去噪过程中对中值滤波处理部分噪声有一定的改进, 但仍有局限性, 噪声强度不同, 这些改进滤波性能表现效果不同[4].基于物理学中热传导方式提出的图像增强算法各向异性扩散模型[14], 是一种有选择性的平滑过程, 这种平滑在均匀区域不受限制, 只在图像边缘区域受到抑制, 从而有选择性的降低了噪声的干扰, 并且对边缘细节有一定的保护作用.而对于高斯噪声去噪方法的研究有很多, 文献[5-9]都是较为经典的高斯噪声去噪方法.本文基于柯西噪声提出一种在中值滤波基础上的改进算法, 该算法由中值滤波及边缘修复两部分组成, 以达到较好的去除柯西噪声效果.

2 噪声模型

$Z=\{Z_{ij} | i \in N, j\in N \}$为原图像, $Y=\{Y_{ij} | i \in N, j\in N \}$为被噪声污染后的图像, $\varepsilon$为所加噪声.定义图像模型

$ Y=Z+\varepsilon, $ (2.1)

其中$\varepsilon$服从标准柯西分布, 为数学期望和方差不存在的连续型分布, 具有自己的分布密度, 满足

(1) 柯西分布分布函数

$ F(x)=\frac{1}{2}+\frac{1}{\arctan(x)}, \quad x\in \mathbb{R} . $ (2.2)

(2) 柯西分布密度函数

$ f(x)=\frac{1}{\pi\cdot(1+x^2)}, \quad x\in \mathbb{R}. $ (2.3)

与正态分布相比, 柯西分布尾部概率更大, 下降至0的速度较慢, 产生强干扰的概率较大.

标准柯西分布也是自由度为$1$$t$ -分布, 因此用Matlab模拟生成标准柯西分布随机数直接生成一组服从$ t(1)$分布的随机数即可. 图 1为标准正态分布与标准柯西分布密度函数图.柯西分布较正态分布尾部厚, 柯西分布尾部概率相对较大, 下降至$0$的速度慢些, 但又与椒盐噪声分布不同, 保持一定的拖尾性, 图 2为图像加了柯西噪声后的效果图, 将Cameraman、Lena两幅图像分别加上2倍、4倍标准柯西分布噪声, 加上噪声后的图像效果与椒盐噪声有些类似, 均产生黑白相间的亮暗点噪声, 但柯西分布所加的部分噪声点与局部区域像素值接近, 因此在用以往中值滤波等线性去噪方法时, 会产生了一定的干扰, 造成一个噪声点"污染"一片的情况.

图 1 标准正态分布与柯西分布密度函数图

图 2 图像加噪效果图

加上柯西噪声后, 图像产生黑白相间的亮暗点噪声, 所加噪声强弱不同, 噪声点的密集程度也不同.

3 处理方法

图像去噪方法较多, 针对不同噪声选取不同的去除噪声方法.中值滤波是基于排序统计理论的一种能有效抑制噪声的非线性信号处理技术.其实现原理为:将某个像素邻域中的像素按灰度值进行排序, 然后选择该序列的中间值作为输出的像素值, 让周围像素灰度值的差比较大的像素改取与周围的像素值接近的值, 从而可以消除孤立的噪声点.利用中值滤波算法可以很好地对图像进行平滑处理.这种算法简单, 时间复杂度低, 但其对点、线和尖顶多的图像不宜采用中值滤波.自适应中值滤波(AMF)是在中值基础上, 根据不同情况, 选择合适的网格大小进行去噪, 是中值滤波的一种改进算法.

各向异性扩散模型(Anisotropic Diffusion)[14-15]模型用下面的方程描述

$ \begin{equation} \frac{\partial{u}}{\partial{t}}={\rm div}(g(|\nabla{u}|)|\nabla{u}|). \end{equation} $ (3.1)

其基本原理是利用梯度算子来辨别由噪声引起的图像梯度变化和由边缘引起的图像梯度变化, 然后用邻域加权平均去除由噪声引起的小梯度变化, 同时保留由边缘引起的大梯度变化, 通过迭代, 直至图像中的噪声被去除.

非局部平均(Non-Local Means, NLM)[16-17]去噪算法, 通过对这些自相似结构块做加权平均来估计参考块的中心点, 从而降低噪声(零均值的高斯白噪声), 尽管NLM取得了很好的去噪效果, 但对原图像的结构信息保护仍不够, 根据灰度分布的相似性决定权值对离散噪声图像$v=\{v_{i} | {i}\in I \}$, 用NLM方法去噪结果计算${NL[v](i)}$, 如下

$ \begin{equation} {NL[v](i)}=\sum \limits_{j=I}{w^{t}(i, j)v(j)}, \end{equation} $ (3.2)

其中, $w^{t}(i, j)$为归一化后权重, 其以来像素与像素之间的相似性, 且$0\leq |w^{t}(i, j)|\leq 1, $ $ \Sigma{w^{t}(i, j)}=1$.像素$i$与像素$j$之间的相似性由它们各自的领域$N(i)$$N(j)$的灰度值相似性决定

$ \begin{equation} w(i, j)=e^{\frac{\|N(i)-N(j)\|_{2}^{2}}{h^2}}, \end{equation} $ (3.3)

其中$\parallel{\cdot}\parallel_2$表示欧式距离, $h $与噪声方差相关, 权重和为${Z(i)=\sum \limits_{j}{w(i, j)}}$, 归一化后的权重为: $w^t(i, j)=\frac{w(i, j)}{Z(i)}$.

4 改进方法

由于图像被加上柯西噪声后, 有一些类似椒盐噪声的黑白点, 首先运用$3*3$中值滤波去除较为极端的噪声点, 图 3将标准柯西分布密度函数分为三个区域, 区域1范围内噪声点易被中值滤波去除, 区域2、区域3中噪声点在被中值滤波处理过程中, 噪声值易对去噪结果产生影响, 接下来我们着重改善这部分被污染的点.

图 3 噪声分布

如果区域2、区域3范围内强噪声点不在图像边缘区域, 可以被中值滤波去除, 当强噪声点处于图像边缘区域时, 中值滤波会产生较大的偏差; 区域1中噪声值较小, 加到图像像素中依然与像素接近, 在中值滤波处理过程中, 噪声值影响不大.

通过中值滤波去噪后的残差分析, 我们发现$3*3$中值滤波后, 残差较大的问题点主要集中在边界以及图像细节较多的部分, 需要针对边界处理做一些修正.然而, 如果对所有存在问题的点都进行处理, 工作量会很大.从理论上来看, 若残差仅为噪声时, 可以计算残差大于30的概率为0.01, 即如果能够较好处理了这些点, 我们就处理了约${1\%}$的点.当然, 根据不同的需要, 也可以选择不同的阈值, 调节需处理点的数目.

$ \begin{equation} {P(residual>30)}=\int^{\infty} _{30}{f(x){\rm d}x}+\int^{-30} _{\infty}{f(x){\rm d}x}=2\times{\int^{\infty} _{30}{f(x){\rm d}x}}=0.0106. \end{equation} $ (4.1)

图 4标注了三幅图像加上2倍、4倍标准柯西噪声, 可以看出中值滤波去噪后残差大于30的点基本处于图像像素突变较大的临界区域以及细节较多的位置, 其中红色像素点代表残差大于30, 绿色像素点代表残差小于$-30$.

图 4 (Color Online)残差标记图像

三幅图像加上2倍、4倍标准柯西噪声中值滤波去噪后, 残差大于30的点集中处于图像像素突变较大的临界区域以及细节较多的位置.

图像边缘总以图像中像素的突变形式出现的, 所以图像边缘中包含着大量有效的信息.由于景物的边缘具有十分复杂的形态, 涉及像素突变, 因此, 最常用的边缘检测方法是"梯度检测法".

4.1 检测边缘, 找到图像边缘点集

边缘检测的基本算法有很多, 有梯度算子、方向算子、拉普拉斯算子以及坎尼(Canny)算子等.几种常用的边缘检测方法有属于梯度算子的Roberts算子、Sobel算子和Prewitt算子.

在算法实现过程中, 通过$2*2$矩阵(Roberts算子)或者$3*3$模板作为核与图像中的每个像素点做卷积和运算, 然后选取合适的阈值以提取边缘. Roberts算子定位比较精确, 但由于不包括平滑, 所以对于噪声比较敏感. Prewitt算子和Sobel算子都是一阶的微分算子, 而前者是平均滤波, 后者是加权平均滤波且检测的图象边缘可能大于2个像素.这两者对灰度渐变低噪声的图象有较好的检测效果, 但是对于混合多复杂噪声的图象, 处理效果就不理想了.梯度算子计算简单, 但精度不高, 只能检测出图象大致的轮廓, 而对于比较细的边缘可能会忽略. Prewitt和Sobel算子比Roberts效果要好一些.文献[10]尽管如此, 考虑到该文提出的后续计算方法, 找到点后需要往四周扩充, 所以此处采用Roberts算子.

Roberts边缘算子是一个$2\times2$的模板, 采用的是对角方向相邻的两个像素之差.从图像处理的实际效果来看, 边缘定位较准, 对噪声敏感.设$f=\sqrt{G}*x^2+G*y^2$, 其中$Gx $为水平梯度, $Gy$为垂直梯度.

表 1 Robert算子

通过设置每一个指标的阈值, 我们可以得到边界细节情况, 图 5显示三幅图像运用Roberts边缘算子找到的图像边缘.将找到的边缘点集向左右上下各扩充一个单位, 即找到这些点周围$9$个点, 将这些点的集合记为$A$. 图 6展示了三幅图像在加上2倍、4倍、8倍标准柯西噪声并用中值滤波做初步去噪处理后, 我们找到的边界点集合$A$元素个数及中值滤波后残差大于30的目标点个数的柱状图对比.

图 5 利用Roberts算子找到的图像边缘图

图 6 边界点个数与残差大于$30$目标点个数对比

三幅图像在加上2倍、4倍、8倍标准柯西噪声并用中值滤波做初步去噪处理后, 我们找到的边界点集合$A$元素个数及中值滤波后残差大于$30$的目标点个数的柱状图对比, 可以看出我们找到的边界点远远多于目标点个数.

表 2 图像边缘覆盖率

覆盖率是指找到的图像边缘点对残差大于30的点的覆盖比率.除图像Hat外, 其他两幅图像边缘点覆盖目标点覆盖率均接近60%, 而Hat图像覆盖率较低, 但并不是所有选出来的边缘点的像素值都需要修正, 修正方法中在4.3节中有详细的描述.

4.2 聚类处理

剔除A中的重复点, 剩余点各取$3*3$区域, 即若$(x, y)\in{A}$, 则取以$(x, y)$为中心的$3*3$矩阵, 即$(x-1, y-1)$$(x-1, y)$$(x, y-1)$$(x, y)$$(x, y+1)$$(x+1, y-1)$$(x+1, y)$$(x+1, y+1)$为一组, 共有$a=\#\{{(x, y)\mid (x, y)\in A}\}$组数据, 放入矩阵$B$中, $B$$a*9$矩阵.

然后运用聚类的方法, 首先将矩阵$B$中为负的去掉(不参与聚类, 但保留位置), 分为两类$C_{big}$$C_{small}$, 记$C_{big}$组中值为$med_{big}$, 记$C_{small}$组中值为$med_{small}$.

4.3 像素取值的修正

计算各个点与8个方向的梯度差, 记正梯度个数$N_+$, 负梯度个数$N_-$, 根据8个梯度的正负个数确定该点取值到底为:中值、原值、大组中值、小组中值.

(ⅰ)若$|N_{+}-N_{-}|<3$, 则取加噪后图像对应像素值;

(ⅱ)若$N_{+}\geq {C_{big}}$, 则取偏小一类的中值${med_{small}}$;

(ⅲ)若$N_{+}\geq {C_{small}}$, 则取偏大一类的中值${med_{big}}$;

(ⅳ)若同时满足(ⅱ)、(ⅲ), 则保留$3*3$中值滤波得到的值;

(ⅴ)若$N_{+}=0$$N_{-}=0$, 则保留$3*3$中值滤波值.

4.4 运行结果

运用本文提出的滤波算法, 先采用中值滤波去除图像中较为明显的噪声点, 而后采用修复边界的思想对部分像素点进行修正, 去噪结果见不同去噪算法信噪比结果图中所示.中值滤波对于柯西噪声已经有较好的处理, 在对此方法改进后信噪比的提升空间并非是无穷大, 如何描述去噪方法对中值滤波的改进效果是一个有趣的问题.

表 3 不同去噪算法信噪比结果

峰值信噪比(Peak Signal-to-Noise Ratio, PSNR)是衡量去噪图像质量的常用定量指标.峰值信噪比定义为

$ MSE=\frac{1}{mn}*\sum \limits_{i=0}^{m-1}\sum\limits_{j=0}^{n-1} {\|{Y(i, j)-Z(i, j)}\|}^2, $ (4.2)
$ PSNR=10*\log_{10}(\frac{255^2}{MSE})=20*\log_{10}(\frac{255}{MSE}). $ (4.3)

不同去噪方法去噪结果来看, 在强噪声情况下, NLM算法边缘模糊较为严重, 本文方法去噪结果最好.

图 7 加2倍柯西噪声去噪算法信噪比结果

去噪结果显示各去噪方法都导致边缘不同程度的模糊, 其中NLM算法边缘模糊较为严重, 本文方法去噪效果较好.

图 8 加4倍柯西噪声去噪算法信噪比结果

在强噪声下, 去噪算法造成图像细节丢失较多, NLM算法图像细节损失严重, 本文方法细节保护较好.

峰值信噪比是一个绝对的概念, 不同的图像, 改进后的峰值信噪比增加$1$, 对不同的图像的改进效果并不相同.为了更加直观描述不同去噪方法的去噪效果, 我们引进改进比例的直观概念来衡量去噪效果.如果中值滤波去噪后, 将残差大于$30$的像素点均恢复为原始图像真值, 此时的图像信噪比是任何调整方法所能达到的峰值信噪比的最大值$c$, 它是我们去噪的一个目标参考值, 可以用来衡量我们改进中值滤波的程度.改进后的图像峰值信噪比越接近$c$, 改进效果越好.令我们的去噪算法得到的图像峰值信噪比为$a$, $3*3$中值滤波峰值信噪比为$b$, 定义改进比率$r$

$ \begin{equation} r=\frac{a-b}{c-b}\times {100 \%}. \end{equation} $ (4.4)

改进比率$r$是一个相对值, 通常$b\le a\le c$, $r$位于$0\sim 100\%$之间, 其大小反应了$a$在最小值$b$和最大值$c$之间的相对位置.

表 4 极端值表格及方法改进百分比

图像Cameraman、Hat改进比率均超过30%, Lena图像细节分布较为集中, 改进效果不大.

5 总结

标准柯西分布与正态分布密度函数形状类似, 但尾部概率不同, 并且柯西分布期望方差不存在的特征决定不能直接采用高斯噪声的去噪算法, 对三幅被Cauchy噪声干扰的图像进行处理和分析, 加噪后效果类似椒盐噪声, 采用中值滤波在较弱噪声时去噪效果较好, 但当噪声较强时, 中值滤波使得图像边缘及细节丢失较为严重, 如何修补细节是研究的侧重点.我们主要通过引进聚类思想, 对中值滤波模糊边缘的问题进行修正.通过边界点局部区域相似像素值的聚类, 采用八个方向的梯度来确定边界点像素取值.我们借助极端值, 引入改进比率的概念来衡量去噪效果.该算法弥补了中值滤波在细节处理方面的不足, 和其它方法相比, 在滤除强噪声并保护图像细节和边缘方面有明显提高, 是一种有效的去除强噪声的滤波算法.从处理结果来看, 该法对中值滤波在强噪声情况下有较好的改进作用.

参考文献
[1] Hwang H, Haddad R A. A Adaptive median filters:new algorithms and results. IEEE Transaction on Image Processing, 1995, 4(4): 499–502. DOI:10.1109/83.370679
[2] Wang Z, Zhang D. Progressive switching median filter for the removal of impulse noise from highly corrupted images. IEEE Trans on Circuits and Systems-Ⅱ: Analog and Digital Signal Processing, 1999, CAS-Ⅱ, 46(1): 78-80
[3] Xing Z J, Wang S J, Dengh J, et al. A new filtering algorithm based on extremum and median value. Journal of Image and Graphics, 2001, 6(6): 533–536.
[4] 武英. 一种有效去除椒盐噪声的滤波算法. 南京晓庄学院学报, 2008, 6: 62–65.
Wu Y. An effective salt and pepper noise filter. Journal of Nanjing Xiaozhuang University, 2008, 6: 62–65. DOI:10.3969/j.issn.1009-7902.2008.06.017
[5] Buades A, Coll B, Morel J M. A review of image denoising methods, with a new one. Multiscale Model Simul, 2005, 4(2): 490–530. DOI:10.1137/040616024
[6] Kervrann C, Boulanger J. Optimal spatial adaptation for patchbased image denoising. IEEE Trans Image Process, 2006, 15(10): 2866–2878. DOI:10.1109/TIP.2006.877529
[7] Takeda H, Farsiu S, Milanfar P. Kernel regression for image processing and reconstruction. IEEE Trans Image Process, 2007, 16(2): 349–366. DOI:10.1109/TIP.2006.888330
[8] Elad M, Aharon M. Image denoising via sparse and redundant representations over learned dictionaries. IEEE Trans Image Process, 2006, 15(12): 3736–3745. DOI:10.1109/TIP.2006.881969
[9] Chatterjee P, Milanfar P. Clustering-based denoising with locally learned dictionaries. IEEE Trans Image Process, 2009, 18(7): 1438–1451. DOI:10.1109/TIP.2009.2018575
[10] 马艳, 张治辉. 几种边缘检测算子的比较. 工矿自动化, 2004, 2004(1): 54–56.
Ma Y, Zhang Z H. The comparison of several kinds of edge detection operator. Industry and Mine Automation, 2004, 2004(1): 54–56. DOI:10.3969/j.issn.1671-251X.2004.01.023
[11] Aharon M, Elad M, Bruckstrein A M. The K-SVD:An algorithm for designing of over-complete dictionaries for sparse representation. IEEE Transactions on Image Processing, 2006, 54(11): 4311–4322. DOI:10.1109/TSP.2006.881199
[12] Elad M, Aharon M. Image denoising via sparse and redundant representations over learned dictionaries. IEEE Transaction on Image Processing, 2006, 15(12): 3736–3745. DOI:10.1109/TIP.2006.881969
[13] Li H B, Liu F. Image Denoising Via Sparse and Redundant Representations Over Learned Dictionaries in Wavelet Domain. IEEE Trans Image Process, 2009, 2009: 754–758.
[14] Perona P, Malik J. Scale space and edge detection using anisotropic diffusion. IEEE Transactions on Pat-tern Analysis and Machine Intelligence, 1990, 12(7): 629–639. DOI:10.1109/34.56205
[15] 余锦华, 汪源源. 基于各向异性扩散的图像降噪算法综述. 电子测量与仪器学报, 2011, 25(2): 105–116.
Yu J H, Wang Y Y. Image noise reduction based on anisotropic diffusion:A survey. J Electronic Meas Instr, 2011, 25(2): 105–116.
[16] 王敏, 王洪剑, 孙光英, 等. 基于NLM的图像三维去噪算法. 红外技术, 2013, 35(4): 238–241.
Wang M, Wang H J, Sun G Y, et al. Three-dimensional image denoising algorithm based on non-local means. Infrared Technology, 2013, 35(4): 238–241.
[17] Charles K, Boulanger J. Optimal spatial adaptation for patch-based image denoising. IEEE Transactions on Image Processing, 2006, 15(10): 2866–2878. DOI:10.1109/TIP.2006.877529