基于改进Directionlet域的广义高斯模型图像去噪方法

摘要

本发明公开了一种基于改进Directionlet域的广义高斯模型图像去噪方法，主要解决现有去噪方法边缘细节丢失严重、均匀区域过平滑的问题。其实现步骤为：(1)选测试图像，加入高斯噪声，得到噪声图像；(2)对噪声图像进行子图分割，确定各子图的变换矩阵；(3)对子图采样，得到陪集；(4)对各陪集进行各向异性小波变换；(5)估计高频子带广义高斯模型的形状参数和局部标准差；(6)由含噪系数估计无噪系数；(7)对无噪系数进行各向异性小波逆变换；(8)根据变换矩阵加权综合，重构各子图；(9)将重构的子图合成，得到去噪结果。本发明具有边缘细节保持好、均匀区域失真少和峰值信噪比高的优点，可用于去除自然图像中的高斯噪声。

说明书

技术领域

本发明属于数字图像处理技术领域，具体地说是一种图像去噪方法，该方法可用于去除被污染的自然图像中的噪声。

背景技术

图像在获取和传输过程中，常会受到噪声的污染。如何对图像进行去噪并尽可能多地保留图像的细节信息以改善图像的质量，成为图像处理中的一个重要任务。人们根据图像的特点、噪声的统计特征和频谱的分布规律，发展了多种去噪方法，其中最为直观的是根据噪声一般集中在高频，而图像频谱则分布在一个有限区间这一特点，采用低通滤波来进行去噪，如滑动平均窗滤波，Wiener滤波等。其他的去噪方法还有基于秩-阶滤波的方法，基于马尔可夫场模型的方法等，这些方法仅具有空域或频域的局部分析性能。

近年来，由于小波具有良好的时频特性和多分辨特性，小波阈值方法被广泛应用到各种去噪处理中，但对经过正交小波变换后的高频系数进行阈值去噪时，会出现“过扼杀”，去噪图像边缘处会出现振荡现象。在小波基础上发展而来的平移不变的平稳小波变换，克服了正交小波变换去噪时的不足，但在阈值去噪时均匀区域出现过平滑现象。为解决二维或更高维奇异性而出现的多尺度几何工具Contourlet，对图像中的奇异曲线的逼近接近最优，Contourlet阈值去噪时边缘信息丢失较少，但在均匀区域会产生虚假成分，即蚊状噪声很明显。

Vladan Velisavljevi′提出的Directionlet变换，是一种新的多尺度几何分析工具。当Directionlet基函数的方向与图像中各向异性目标的方向匹配时，对图像的逼近效果较好，可以稀疏的表示图像的边缘等细节信息；不匹配时Directionlet则退化为小波。但现有的Directionlet变换，变换方向和队列方向都是任意选取的，比如0度，90度，45度和-45度，并没有根据图像的特征进行方向自适应的变换，导致在图像去噪时边缘信息丢失严重，出现模糊失真现象。

发明内容

本发明的目的在于克服上述已有技术的不足，提供一种基于改进Directionlet域的广义高斯模型图像去噪方法，以根据图像的特征进行方向自适应的变换，实现在去除噪声的同时，尽可能保持图像边缘等细节特征的目的，提高图像的质量。

实现本发明目的的技术方案是，首先自适应确定Directionlet变换的变换矩阵，使它根据图像的主要方向进行变换，在去噪处理时，用广义高斯模型对Directionlet系数建模，根据广义高斯模型形状参数的大小，结合局部标准差采取不同的无噪系数估计策略，对图像进行去噪处理，其具体去噪步骤包括如下：

(1)选取测试图像，加入零均值的高斯噪声，得到噪声图像；

(2)对噪声图像进行64×64的子图分割，并用二进小波变换自适应确定各分割子图的Directionlet变换矩阵M_Λ；

(3)利用变换矩阵M_Λ，对各分割子图进行采样，得到分割子图的|det(M_Λ)|个陪集，|det(M_Λ)|是矩阵M_Λ行列式的绝对值；

(4)对各分割子图的每个陪集沿Directionlet变换矩阵M_Λ的变换方向和队列方向分别进行n₁＝2与n₂＝1次的一维小波变换，得到Directionlet变换的高频和低频子带系数；

(5)对各个高频子带，利用该子带所有的变换系数，估计广义高斯模型的形状参数υ和局部标准差σ_x；

(6)对各高频子带的广义高斯模型形状参数υ进行判断：

若0＜υ＜0.5，按下式对噪声图像的无噪系数进行估计，

$\hat{x}\left(y\right)=\left(\begin{array}{cc}0,& \left|y\right|<T_{\upsilon }\\ y-{\sigma }^2{{\sigma }_x}^{-\upsilon }\mathrm{\upsilon \eta }\left(\upsilon \right)y^{\upsilon -1}+o\left(y^{2(\upsilon -1)}\right),& \left|y\right|\ge T_{\upsilon }\end{array}\right)$

其中是无噪系数的估计，y是含噪系数，T_υ＝C_υσ^2/(2-υ)σ_x^-(υ/(2-υ))，C_υ＝(2-υ)(2-2υ)^-(1-υ/2-υ)η(υ)^υ/(2-υ)，$\eta \left(\upsilon \right)=\sqrt{\Gamma (3/\upsilon )/\Gamma (1/\upsilon )},$Γ是Gamma函数，$\Gamma \left(t\right)={\int }_0^{\infty }{e}^{-u}{u}^{t-1}\mathrm{du},$o(y^2(υ-1))是y^2(υ-1)的高阶无穷小量，σ是噪声标准差，Y＝(y₁，y₂，...，y_N)是含噪高频子带系数，N是高频子带系数的数目；

若0.5≤υ＜1，则采用阈值$T=\frac{{\sigma }^2}{{\sigma }_x}$进行软阈值处理；

(7)对低频子带和估计得到的无噪高频子带沿矩阵M_Λ的变换方向和队列方向分别进行n₁＝2与n₂＝1次的一维小波逆变换；

(8)根据变换矩阵M_Λ的变换方向和队列方向加权综合，重构各个分割子图；

(9)将重构的各分割子图按其在原图像中的位置合成，得到去噪后的图像。

本发明与现有技术相比具有如下优点：

(1)本发明采用Directionlet变换的矩阵是根据图像的特征自适应确定的，由于Directionlet变换系数对图像中的边缘各向异性特征逼近较好，对图像去噪时边缘细节信息保持较好，因而解决了现有Directionlet变换方法由于任意选取变换矩阵而退化成小波变换，在图像去噪处理中因对边缘等各向异性特征逼近不够精确，边缘信息丢失严重而引起的边缘模糊问题；

(2)本发明由于在无噪系数估计时充分考虑了广义高斯模型形状参数的影响，采用的局部标准差具有局部适应性，能更精确地估计高频子带的各个系数，解决了现有小波阈值方法对含噪系数的“过扼杀”问题和平稳小波阈值方法均匀区域过平滑问题，且本发明方法去噪后图像的均匀区域更清晰，失真较少。

仿真实验结果表明，本发明在对加有高斯噪声的自然图像去噪中，边缘等细节信息保持较好，均匀区域更清晰，得到了更高的峰值信噪比，提高了图像的质量。

附图说明

图1是本发明的主要操作过程示意图；

图2是用本发明方法与已有去噪方法对测试图像Lena的去噪效果对比图；

图3是用本发明方法与已有去噪方法对测试图像Barbara的去噪效果对比图；

图4是用本发明方法与已有去噪方法对测试图像Peppers的去噪效果对比图。

具体实施方式

参照图1，本发明的具体实现步骤如下：

步骤1：选取测试图像，加入高斯噪声，得到噪声图像，该噪声图像的数学表达式为：

f＝g+n

其中g＝{g(i₁，j₁)|i₁，j₁＝1，2，...N}表示测试图像，n＝{n(i₁，j₁)|i₁，j₁＝1，2，...N}表示均值为零方差为σ²的高斯噪声，噪声图像记为f＝{f(i₁，j₁)|i₁，j₁＝1，2，...N}，N表示图像大小。

步骤2：把噪声图像分割成互不重叠，且大小为64×64的子图，寻找每个分割子图的两个主要方向θ₁和θ₂，根据主要方向θ₁和θ₂确定各分割子图的Directionlet变换矩阵M_Λ，具体步骤如下：

2a)对分割子图进行二进小波变换，得到其水平细节图h(i，j)和垂直细节图v(i，j)，其中(i，j)是二进小波变换系数的位置，i，j＝1，2，...64；

2b)根据h(i，j)和v(i，j)，计算分割子图在(i，j)处的方向θ(i，j)：

若h(i，j)＜＜v(i，j)，即$\left|\right|\frac{\left|v(i,j)\right|-\left|h(i,j)\right|}{\left|v(i,j)\right|}|-1|<0.05,$$\theta (i,j)=\frac{\pi }{2};$

若h(i，j)＞＞v(i，j)，即$\left|\right|\frac{\left|v(i,j)\right|-\left|h(i,j)\right|}{\left|h(i,j)\right|}|-1|<0.05,$θ(i，j)＝0；

若$\left|\right|\frac{\left|v(i,j)\right|-\left|h(i,j)\right|}{\left|v(i,j)\right|}|-1|\ge 0.05$或$\left|\right|\frac{\left|v(i,j)\right|-\left|h(i,j)\right|}{\left|h(i,j)\right|}|-1|\ge 0.05,$$\theta (i,j)=\arctan \frac{v(i,j)}{h(i,j)};$

2c)统计分割子图的方向θ(i，j)的分布，找出出现次数最多的两个方向θ₁和θ₂；

2d)分别求两个主要方向θ₁和θ₂的反正切值，得到两个近似有理斜率r₁和r₂，r₁≈arctanθ₁＝b₁/a₁，r₂≈arctanθ₂＝b₂/a₂，根据有理斜率r₁和r₂，构造变换矩阵$M_{\Lambda }=\left(\begin{array}{cc}{a}_1& b_1\\ a_2& b_2\end{array}\right),$其中沿r₁的方向称为变换矩阵M_Λ的变换方向，沿r₂的方向称为队列方向，a₁，a₂，b₁，b₂都是整数。

步骤3：利用变换矩阵M_Λ，对各分割子图分别沿r₁和r₂方向进行采样，得到分割子图的|det(M_Λ)|个陪集，|det(M_Λ)|是矩阵M_Λ行列式的绝对值。

步骤4：对各分割子图的每个陪集沿r₁方向进行n₁＝2次的一维小波变换，沿r₂方向进行n₂＝1次的一维小波变换，得到Directionlet变换的高频和低频子带系数。

步骤5：由含噪高频子带系数估计广义高斯模型的形状参数υ和局部标准差σ_x按如下步骤进行：

3a)设高频子带的大小为l₁×l₂，由含噪的高频子带系数Y按下式计算出含噪高频子带系数的二阶矩σ_Y²和含噪高频子带系数的峰度k_Y：

${{\sigma }_Y}^2=\frac{1}{l_1\times l_2}\underset{i_2,j_2=1}{\overset{l_1,l_2}{\Sigma }}{Y^2}_{i_2,j_2}$

$k_Y=\frac{1}{{{\sigma }_Y}^4}\frac{1}{l_1\times l_2}\underset{i_2,j_2=1}{\overset{l_1,l_2}{\Sigma }}{Y^4}_{i_2,j_2}$

σ_Y⁴＝(σ_Y²)²

式中是(i₂，j₂)处含噪高频系数的平方，是(i₂，j₂)处含噪高频系数的4次方；

3b)由含噪高频子带系数Y，利用中值估计法求出噪声标准差：$\sigma =\frac{\mathrm{Median}\left(\right|Y\left|\right)}{0.6745},$再根据噪声标准差σ求出无噪系数的标准差：${\sigma }^{\prime }=\sqrt{\max ({{\sigma }_Y}^2-{\sigma }^2,0)};$

3c)由以上得到的k_Y，σ_Y²和σ值，利用数值计算法按照如下公式求出形状参数υ：

$k_Y=\frac{1}{{{\sigma }_Y}^4}[{{6\sigma }_Y}^2{\sigma }^2-3{\sigma }^2+{({{\sigma }_Y}^2-{\sigma }^2)}^2\frac{\Gamma \left(\frac{1}{\upsilon }\right)\Gamma \left(\frac{5}{\upsilon }\right)}{{\Gamma }^2\left(\frac{3}{\upsilon }\right)}];$

3d)由含噪高频子带系数Y，利用最小均方误差估计法得到无噪系数的初始估计值$X_{i_2,j_2}=\frac{{\sigma }^{\prime 2}}{{\sigma }^{\prime 2}+{\sigma }^2}{Y}_{i_2,j_2}$

其中，是(i₂，j₂)处的含噪高频系数，i₂＝1，2，...l₁，j₂＝1，2，...l₂；

3e)由形状参数υ和无噪系数的初始估计值X，计算出局部标准差：σ_x＝SE

其中S是只与形状参数υ有关的因子，$S={\upsilon }^{1/\upsilon }{\left[\frac{\Gamma (3/\upsilon )}{\Gamma (1/\upsilon )}\right]}^{1/2},$

E是大小为5×5的局部邻域内无噪系数初始估计值的绝对均值，

$E=\frac{\underset{p=-2,q=-2}{\overset{\mathrm{2,2}}{\Sigma }}\left|X_{i_3+p,j_3+q}\right|}{5\times 5},$是(i₃+p，j₃+q)处无噪系数的初始估计值，

i₃＝2，...l₁-1，j₃＝2，...l₂-1。

步骤6：将形状参数υ分成两个区间0＜υ＜0.5和0.5≤υ＜1，使用以下估计策略由含噪系数估计无噪系数：

若0＜υ＜0.5，按下式对噪声图像的无噪系数进行估计，

$\hat{x}\left(y\right)=\left(\begin{array}{cc}0,& \left|y\right|<T_{\upsilon }\\ y-{\sigma }^2{{\sigma }_x}^{-\upsilon }\mathrm{\upsilon \eta }\left(\upsilon \right)y^{\upsilon -1}+o\left(y^{2(\upsilon -1)}\right),& \left|y\right|\ge T_{\upsilon }\end{array}\right)$

其中是无噪系数的估计，y是含噪系数，T_υ＝C_υσ^2/(2-υ)σ_x^-(υ/(2-υ))，C_υ＝(2-υ)(2-2υ)^-(1-υ/2-υ)η(υ)^υ/(2-υ)，$\eta \left(\upsilon \right)=\sqrt{\Gamma (3/\upsilon )/\Gamma (1/\upsilon )},$Γ是Gamma函数，$\Gamma \left(t\right)={\int }_0^{\infty }{e}^{-u}{u}^{t-1}\mathrm{du},$o(y^2(υ-1))是y^2(υ-1)的高阶无穷小量，σ是噪声标准差；

若0.5≤υ＜1，则采用阈值$T=\frac{{\sigma }^2}{{\sigma }_x}$对含噪系数进行软阈值处理，得到无噪系数的估计。

步骤7：对低频子带和估计得到的无噪高频子带沿r₁方向进行n₁＝2次的一维小波逆变换，沿r₂方向进行n₂＝1次的一维小波逆变换。

步骤8：对步骤7得到的结果沿r₁方向和r₂方向进行加权综合，重构出各个分割子图。

步骤9：将重构的各分割子图按其在原图像中的位置合成，得到去噪后的图像。

本发明的效果通过以下仿真进一步说明。

一、仿真条件

采用图像去噪中常用的标准图像Lena 512×512、Barbara 512×512和Peppers512×512图像，并对这三幅图像分别加入标准差σ为10，15，20，25，30，35，40，45和50的零均值加性高斯白噪声，分别用本发明去噪方法、小波软阈值去噪方法、平稳小波硬阈值去噪方法、Contourlet硬阈值去噪方法和Contourlet软阈值去噪方法进行仿真去噪处理。

二、仿真结果分析

仿真1，用本发明方法与现有去噪方法对测试图像Lena的去噪效果对比图，结果如图2，其中：

图2(a)是测试图像Lena；图2(b)是σ＝20的Lena含噪图像，峰值信噪比PSNR＝22.13dB；图2(c)是小波软阈值的去噪结果图，PSNR＝26.74dB；图2(d)是平稳小波硬阈值去噪结果图，PSNR＝28.38dB；图2(e)是Contourlet硬阈值去噪结果图，PSNR＝28.93dB；图2(f)Contourlet软阈值去噪结果图，PSNR＝26.76dB；图2(g)是本发明方法去噪结果图，PSNR＝30.19dB；

从图2(a)和图2(c)可以看出，现有小波软阈值方法去噪后Lena的帽沿处出现震荡现象，脸部和其它均匀区域失真严重。

从图2(d)可以看出，现有平稳小波硬阈值方法去噪后虽然Lena的帽沿比较清晰，但帽子上的纹理有些模糊，均匀区域出现过平滑。

从图2(e)和图2(f)可以看出，现有Contourlet阈值方法去噪后，产生了虚假成分，即蚊状噪声很明显。

从图2(g)可以看出本发明方法去噪后Lena的帽沿及帽子上的纹理比较清晰，脸部和其它均匀区域较平滑，失真较少，并且有较高的峰值信噪比，有效去除了图2(b)的噪声。

仿真2，用本发明方法与现有去噪方法对测试图像Barbara的去噪效果对比图，结果如图3，其中：

图3(a)是测试图像Barbara；图3(b)是σ＝20的Barbara含噪图像，PSNR＝22.19dB；图3(c)是小波软阈值方法的去噪结果图，PSNR＝23.67dB；图3(d)是平稳小波硬阈值方法去噪结果图，PSNR＝25.39dB；图3(e)是Contourlet硬阈值方法去噪结果图，PSNR＝25.86dB；图3(f)Contourlet软阈值方法去噪结果图，PSNR＝23.87dB；图3(g)本发明方法去噪结果图，PSNR＝27.87dB。

从图3(a)和图3(c)可以看出，现有小波软阈值方法去噪后Barbara裤子和围巾上的纹理信息丢失严重。

从图5(d)可以看出，现有平稳小波硬阈值方法去噪后Barbara裤子和围巾上的纹理信息有些丢失，桌子腿部分出现震荡。

从图3(e)和图3(f)可以看出，现有Contourlet阈值去噪方法去噪后，Barbara纹理保持较好，图像较清晰，但产生了严重的虚假成分，即蚊状噪声很明显。

从图3(g)可以看出，本发明方法在对图3(b)去除噪声的同时更多的保持了Barbara裤子，围巾等处的纹理细节，可见本发明在细节保持方面有相当大的优势，均匀区域也更清晰。

仿真3，用本发明方法与现有去噪方法对测试图像Peppers的去噪效果对比图，结果如图4，其中：

图4(a)是测试图像Peppers；图4(b)是σ＝20的Peppers含噪图像，PSNR＝22.21dB；图4(c)是小波软阈值方法的去噪结果图，PSNR＝27.21dB；图4(d)是平稳小波硬阈值方法去噪结果图，PSNR＝28.59dB；图4(e)是Contourlet硬阈值方法去噪结果图，PSNR＝28.65dB；图4(f)Contourlet软阈值方法去噪结果图，PSNR＝26.51dB；图4(g)本发明方法去噪结果图，PSNR＝30.82dB。

从图4(a)和图4(c)可以看出，现有小波软阈值方法去噪后Peppers边缘处出现震荡现象。

从图4(d)可以看出，现有平稳小波硬阈值方法去噪后，Peppers均匀区域有过平滑现象。

从图4(e)和图4(f)可以看出，现有Contourlet阈值去噪方法去噪较彻底，但蚊状噪声很明显。

从图4(g)可以看出，本发明方法对图4(b)去噪后各种蔬菜的边缘很清晰，均匀区域较清晰，和测试图像图4(a)相比差别更小，明显比其他方法去噪后的效果要好。

综上，本发明方法不仅视觉效果好，而且峰值信噪比也比较高，在边缘保持方面有较大优势。