CN101493933A - 一种局部结构自适应的图像扩散去噪方法 - Google Patents

Abstract

本发明提供一种局部结构自适应的图像扩散去噪方法。特别涉及一种图像局部结构的分类方法、以及一种根据局部结构确定扩散传导系数的方法。图像扩散是指模拟物理学的热扩散形式，对图像像素值进行扩散平滑处理，达到去除噪声的目的。各向异性扩散能够在平滑图像的同时保持边缘、线条等细节信息。本发明所提出的是一种能够根据图像局部细小结构调整不同方向扩散强度的各向异性扩散方法。首先根据局部结构将像素分类，对不同类型定义不同的传导系数函数，再由这些函数根据梯度计算出扩散流量，进而实现保持边缘信息的图像扩散去噪。

Description

一种局部结构自适应的图像扩散去噪方法

技术领域

本发明属于信息技术中的数字图像处理技术领域，涉及图像像素的分类特征值的提取、图像像素的分类方法，以及图像扩散流量的计算方法。

背景技术

图像扩散是指的是模拟物理学的热扩散形式，对图像像素值进行扩散平滑处理，达到去除噪声的目的。

传导系数函数是指用来根据图像梯度计算扩散流量的具体函数形式，一般采用负指数函数等非单调负递减函数。

梯度是指图像中相邻像素之间的灰度值差异，是一个具有一定方向的矢量，包括水平方向和垂直方向两个分量。

扩散流量是指根据图像中相邻像素之间的梯度等信息计算出某个像素受其邻域像素值影响而改变的量值大小。

基于偏微分方程的图像扩散去噪方法的研究始于1990年，经过近20年的发展，已涌现出许多图像扩散算法。然而过去的研究主要集中于如何利用不同传导系数函数根据梯度值计算扩散流量，扩散流量只依赖于梯度。但相同的梯度在图像中不同的位置代表着不同的信息，这些算法不能很好地区分图像中的有用细节信息和噪声，造成去噪能力的下降和图像的模糊。

文献1(IEEE Trans.Pattern Analysis and Machine Intelligence，1990，12(7)：629)报道的非线性P-M扩散算法，利用图像中的梯度信息和一个单调递减传导系数函数进行非线性扩散平滑，可在平滑噪声的同时保留图像边缘和线条等细节信息。但它的传导系数函数只依赖于最近邻的梯度值，它无法滤除边缘强度较大的脉冲噪声(如椒盐噪声)，而且会模糊图像中的弱边缘。

文献2(SIAM J.Numerical Analysis，1992，29(1)：182)报道的正则化P-M扩散算法，在计算图像梯度时使用了高斯滤波器，使得梯度的计算更为可靠，利用图像中的梯度信息和一个单调递减传导系数函数进行非线性扩散平滑，可以去除孤立点噪声，但它同样会模糊图像中的边缘和线条。

文献3(International Journal of Computer Vision，1999，31(2)：111)报道了一种采用矩阵张量D(对2-D图像是2*2的矩阵)作为传导系数函数的一致性扩散算法，根据局部梯度信息，在不同的方向采用不同的传导系数函数。算法可较好保留图像中的纹理信息，包括有规则的线条和边缘，但有时会在没有纹理的地方生成奇怪的纹理。

文献4(IEEE Trans.Image Processing，2002，11(7)：689)报道了一种提出一种前后向扩散(FAB：Forward-And-Backward)算法，可以在平滑扩散的同时增强图像中的边缘和线条细节，但也往往会增强一些不需要的噪声。另外，和P-M扩散相似，对一些弱对比度的边缘和线条，也会被模糊。

文献5(IEEE Trans.Pattern Analysis and Machine Intelligence，2005，27(10)：1552)报道了一种基于图像局部不一致性的扩散算法。根据上下文环境，在不同的位置采用不同的尺度。扩散流量由不一致性和梯度共同决定。但算法中定义的不一致性并不能很好地区分图像结构信息与噪声，它对较大幅度的噪声不能有效地去除。

现有的算法在对传导系数函数自适应中，只考虑局部梯度方向、方差、均值等统计量，而没有对局部结构进行更深入和具体的分析。

发明内容

本发明的目的在于克服上述几种算法没有利用图像局部结构信息的缺点。基于偏微分方程的图像扩散去噪方法有许多，然而过去的研究主要集中于如何利用不同传导系数函数根据梯度值计算扩散流量，扩散流量只依赖于梯度值。但相同的梯度在图像中不同的位置代表着不同的信息，本发明首次将基于局部结构分类与图像扩散相结合，首先对像素进行分类，再对不同的类型采用不同的传导系数函数，即可以有效地保留了图像中的各种细节信息，如边缘、细线条等，也可以同时有效去除不同类型的混合噪声，如高斯噪声和脉冲噪声。

本发明的具体技术方案包括以下步骤：

(1)确定用于图像像素分类的局部结构块。

用于进行像素分类的局部结构图像块如附图1所示，为当前分类像素及其8-邻域像素组成的3*3图像区域。图中’O’表示当前像素，’P₁～P₈’表示它的8-邻域像素。

(2)提取用于图像像素分类的局部结构的分类特征值。局部结构的分类特征值包括像素的局部结构中其他像素与它的灰度最大距离、最小距离、最小绝对距离，以及绝对距离值差；所述的像素分类方法根据局部结构的分类特征值将像素分为3类：孤立点、图像结构、平坦或缓变区；

用于进行像素分类的特征量包括以下几个变量：

a)邻域灰度距离：D_i＝f(P_i)-f(O)(i＝1，2，...，8)

b)最大距离：D_max＝max(D_i)(i＝1，2，...，8)

c)最小距离：D_min＝min(D_i)(i＝1，2，...，8)

d)最小绝对距离：AD_min＝min(abs(D_i))(i＝1，2，...，8)

e)绝对距离值差：AD_min max＝max(abs(D_i))-min(abs(D_i))(i＝1，2，...，8)

(3)根据图像某个像素的分类特征值对像素进行分类。

根据上述分类特征将每个像素分为孤立点用Type I表示，图像结构用Type II表示、平坦或缓变区用TypeIII表示三种类别之一，像素分类规则如下：

if D_min＞T or D_max＜-T

TypeI(孤立点)

else if AD_min＜T and AD_min max＞T

TypeII(图像结构)

else

Type III(平坦或缓变区)

T为分类阈值，物理意义上表示我们识别图像细节信息的最小灰度差异。对不同类型的图有一定差异，可根据经验设定(比如取灰度值范围的10％)，也可根据预先对图像的某种统计信息(如方差)来估计。

像素分类的示例图如附图2～图7所示。

(4)对不同类别的像素设置不同类型的传导系数函数。

传导系数函数，包括3种类型，对应于3种像素类别，分别是最小最大函数、指数函数的幂函数、以及指数函数；

根据上述分类结果，对三种不同类别的像素，分别定义其传导系数函数f(s)如下：

TypeI

$f\left(s\right)=\left\{\begin{array}{c}\min \left(s\right),D_{\min }>T\\ \max \left(s\right),D_{\max }<-T\end{array}\right.$

TypeII

$\left\{\begin{array}{c}f\left(s\right)={\left(c\left(s\right)\right)}^2,D_i\&GreaterEqual;T\\ f\left(s\right)=\sqrt{c\left(s\right)},D_i<T\end{array}\right.$

TypeIII

f(s)＝c(s)

其中

$c\left(s\right)=\frac{1}{1+{(s/k)}^2},$ 或 $c\left(s\right)=e^{-{(s/k)}^2}$

对于孤立点(第I类)，属于噪声，而且其与周围像素的最小绝对差值越大(D_min或D_max的绝对值越大)，是噪声的概率就越大，因此将传导系数定义为该点与周围像素的最小绝对差值；对于包含图像细节信息的第II类，局部结构中即有较强的边缘(D_i≥T)，又有较弱的边缘(D_i＜T)。此时，对较强的边缘采用较小的扩散系数(f(s)＝(c(s))²)，而对较弱的边缘采用较大的扩散系数) $(f\left(s\right)=\sqrt{c\left(s\right)});$ 对第III类，当前像素处于平坦或缓变区域，选择与P-M算法相同的扩散系数f(s)＝c(s)。

(5)根据传导系数函数计算扩散流量并对图像进行图像扩散平滑。扩散流量，根据确定的传导系数函数计算，利用图像梯度值计算任意一个像素与其8-邻域像素的扩散流量对图像进行图像扩散平滑，包括3种类型，对应于3种像素类别，分别是最小最大函数、指数函数的幂函数、以及指数函数；图像扩散平滑，根据确定的邻域各个方向的扩散流量，修改当前像素的灰度值，以达到在保持边缘的情况下平滑去噪的目的；

扩散方程的一般形式定义为：

$I_t=\mathrm{div}(f(\&dtri;I)\&dtri;I)$

其中f(s)为步骤(4)中定义的传导系数函数。具体的数字实现过程如下：

表示图像8个方向的梯度值。若记(x，y)点的8方向边缘为E(x，y，i)，i＝1，...，8，在离散坐标的数字图像中分别如下：

E(x，y，1)＝I(x，y-1)-I(x，y)

E(x，y，2)＝I(x，y+1)-I(x，y)

E(x，y，3)＝I(x-1，y)-I(x，y)

E(x，y，4)＝I(x+1，y)-I(x，y)

E(x，y，5)＝I(x-1，y-1)-I(x，y)

E(x，y，6)＝I(x-1，y+1)-I(x，y)

E(x，y，7)＝I(x+1，y-1)-I(x，y)

E(x，y，8)＝I(x+1，y+1)-I(x，y)

div表示散度，在离散坐标的数字图像中，平滑扩散流量计算公式为：

$I_t=\underset{i=1}{\overset{8}{\mathrm{\&Sigma;}}}f(E(x,y,i)\&CenterDot;E(x,y,i))$

其中c(s)函数中的系数k为各方向边缘绝对值均值乘以某个[1～5]之间的常数。

根据流量对图像进行平滑的随迭代时间演化公式为：

Iⁿ⁺¹＝Iⁿ+λI_i

λ(≤0.125)是一个控制扩散速度的常数。

本发明的整个实现流程图如附图8所示。对输入的数字图像，首先根据每个像素及其邻域像素构成的局部结构块对它进行分类，同时计算其8个方向的梯度值；接着，根据像素的分类结果选择不同的传导系数函数；第三步，根据传导系数函数和梯度值计算扩散流量；第四步，利用扩散流量对图像进行扩散处理，并以处理结果更新替代原来的图像。重复以上四步操作，直到达到指定的迭代次数。

本发明的所提出方法的处理效果如附图11所示。

本发明所提出的算法正是为了克服这一局限性，在对图像局部结构进行初步识别的基础上进行传导系数函数调整。

附图说明

图1用于进行像素分类的局部结构块示意图。

图2孤立点示例1。

图3孤立点示例2。

图4图像结构示例1。

图5图像结构示例2。

图6平坦或缓变区示例1。

图7平坦或缓变区示例2。

图8局部自适应的图像扩散方法流程图。

图9带有各种不同强度结构信息的测试图像。

图10加入高斯噪声和椒盐噪声后的测试图像。

图11扩散去噪处理结果。

具体实施方式

(1)读入一幅带有噪声的灰度数字图像，如附图5所示，记为J0；重复(2)~(5)步迭代10次；

(2)对图像J0中的每个像素计算其8-邻域的梯度值，超出边界部分以边界像素值代替；

(3)计算所有梯度值的绝对值的平均值k；

(4)对每个像素：

a)根据其8-邻域梯度值统计分类特征值；

b)根据其分类特征对此像素进行分类；

c)根据不同类别得到不同方向的传导系数函数；

d)利用传导系数函数、梯度值和k值计算每个方向的扩散流量；

e)根据扩散流量得到更新后的像素值，保存于另一幅图像J1中；

(5)令J0＝J1，返回(2)；

(6)输出结果图像J1。

扩散去噪处理结果如图11所示。

Claims (1)

1、一种局部结构自适应的图像扩散去噪方法，其特征在于，方法包括以下计算机可实现的步骤：

步骤一、确定用于图像像素分类的局部结构：局部结构是选取图像中当前处理的像素及其周围8邻域，其计9个像素所组成的像素集中像素灰度值的分布形式；

步骤二、提取用于图像像素分类的局部结构的分类特征值，局部结构的分类特征值包括像素的局部结构中其他像素与它的灰度最大距离、最小距离、最小绝对距离，以及绝对距离值差；所述的像素分类方法根据局部结构的分类特征值将像素分为3类：孤立点、图像结构、平坦或缓变区；用于进行像素分类的特征量包括以下几个变量：

a)邻域灰度距离：D_i＝f(P_i)-f(O)(i＝1，2，...，8)

b)最大距离：D_max＝max(D_i)(i＝1，2，...，8)

c)最小距离：D_min＝min(D_i)(i＝1，2，...，8)

d)最小绝对距离：AD_min＝min(abs(D_i))(i＝1，2，...，8)

e)绝对距离值差：AD_min max＝max(abs(D_i))-min(abs(D_i))(i＝1，2，...，8)；

步骤三、根据图像某个像素的分类特征值对像素进行分类；

根据上述分类特征将每个像素分为孤立点用Type I表示，图像结构用Type II表示、平坦或缓变区用Type III表示三种类别之一，像素分类规则如下：

if D_min＞T or D_max＜-T

TypeI

else if AD_min＜T and AD_min max＞T

TypeII

else

TypeIII

T为分类阈值，物理意义上表示我们识别图像细节信息的最小灰度差异；

步骤四、对不同类别的像素设置不同类型的传导系数函数；传导系数函数，包括3种类型，对应于3种像素类别，分别是最小最大函数、指数函数的幂函数、以及指数函数；

根据上述分类结果，对三种不同类别的像素，分别定义其传导系数函数f(s)如下：

TypeI

$f\left(s\right)=\left\{\begin{array}{c}\min \left(s\right),D_{\min }>T\\ \max \left(s\right),D_{\max }<-T\end{array}\right.$

TypeII

$\left\{\begin{array}{c}f\left(s\right)={\left(c\left(s\right)\right)}^2,D_i\&GreaterEqual;T\\ f\left(s\right)=\sqrt{c\left(s\right)},D_i<T\end{array}\right.$

TypeIII

f(s)＝c(s)

其中 $c\left(s\right)=\frac{1}{1+{(s/k)}^2},$ 或 $c\left(s\right)=e^{-{(s/k)}^2}$

对于孤立点，将传导系数定义为该点与周围像素的最小绝对差值；对于包含图像细节信息的图像结构，当D_i≥T为较强的边缘采用较小的扩散系数，当D_i＜T为较弱的边缘采用较大的扩散系数；对平坦或缓变区，当前像素处于平坦或缓变区域，选择与P-M算法相同的扩散系数；

步骤五、根据传导系数函数计算扩散流量并对图像进行图像扩散平滑，扩散流量，根据确定的传导系数函数计算，利用图像梯度值计算任意一个像素与其8-邻域像素的扩散流量对图像进行图像扩散平滑，包括3种类型，对应于3种像素类别，分别是最小最大函数、指数函数的幂函数、以及指数函数；图像扩散平滑，根据确定的邻域各个方向的扩散流量，修改当前像素的灰度值，以达到在保持边缘的情况下平滑去噪的目的；

扩散方程的形式定义为：

$I_t=\mathrm{div}(f(\&dtri;I)\&dtri;I)$

其中f(s)为步骤(4)中定义的传导系数函数，具体的数字实现过程如下：

表示图像8个方向的梯度值，记(x，y)点的8方向边缘为E(x，y，i)，i＝1，...，8，在离散坐标的数字图像中分别如下：

E(x，y，1)＝I(x，y-1)-I(x，y)

E(x，y，2)＝I(x，y+1)-I(x，y)

E(x，y，3)＝I(x-1，y)-I(x，y)

E(x，y，4)＝I(x+1，y)-I(x，y)

E(x，y，5)＝I(x-1，y-1)-I(x，y)

E(x，y，6)＝I(x-1，y+1)-I(x，y)

E(x，y，7)＝I(x+1，y-1)-I(x，y)

E(x，y，8)＝I(x+1，y+1)-I(x，y)

div表示散度，在离散坐标的数字图像中，平滑扩散流量计算公式为：

$I_t=\underset{i=1}{\overset{8}{\mathrm{\&Sigma;}}}f\left(E(x,y,i)\right)\&CenterDot;E(x,y,i)$

其中c(s)函数中的系数k为各方向边缘绝对值均值乘以某个[1～5]之间的常数；

根据流量对图像进行平滑的随迭代时间演化公式为：

Iⁿ⁺¹＝Iⁿ+λI_t

λ(≤0.125)是一个控制扩散速度的常数。

Images