一種震動(dòng)耦合擴(kuò)散的去除乘性噪聲模型

,

(中原工學(xué)院，鄭州 450007)

乘性噪聲廣泛存在于現(xiàn)實(shí)生活的很多應(yīng)用當(dāng)中,如核磁共振、合成孔徑雷達(dá)、遙感、超聲波、激光等相干成像領(lǐng)域當(dāng)中。與加性高斯噪聲不同,乘性噪聲服從Gamma分布。目前去除加性噪聲的研究較為成熟,傳統(tǒng)的去除加性高斯噪聲的模型對(duì)于去除乘性噪聲并不適用。而乘性噪聲對(duì)圖像的污染很嚴(yán)重,降低了圖像的畫(huà)面質(zhì)量,對(duì)圖像分割、目標(biāo)檢測(cè)以及邊緣提取等都有不同程度影響[1-2]。

Rudin L I等在經(jīng)典的去除加性噪聲的全變分模型(ROF模型)基礎(chǔ)上首次提出了去除乘性噪聲的RLO模型[3-4]:

(1)

Aubert G和Aujol J A提出了利用最大后驗(yàn)概率(Maximum a posterioro,MAP)估計(jì)方法去除乘性噪聲的AA模型[5-6]:

(2)

Jin Z M等提出了利用對(duì)數(shù)變換z=logu去除乘性噪聲的JY模型[7]:

(3)

Huang Y M等提出了去除乘性噪聲的HNW模型[8]:

(4)

上述這些模型的正則項(xiàng)均為各項(xiàng)同性的全變分,在去除噪聲的同時(shí)存在一些不足之處：邊緣模糊、紋理不清晰、平滑區(qū)域產(chǎn)生“階梯效應(yīng)”。本文將震動(dòng)濾波引入去除乘性噪聲模型中,去噪模型由擴(kuò)散項(xiàng)、震動(dòng)濾波和保真項(xiàng)構(gòu)成。BZN用來(lái)表示圖像中一個(gè)小鄰域內(nèi)由類似成分組成的像素的個(gè)數(shù)。采用BZN可區(qū)分邊緣、噪聲以及平滑區(qū)域[9],進(jìn)而在用此模型去噪的同時(shí)能夠較好地保留圖像的邊緣和細(xì)節(jié)信息。

1 震動(dòng)濾波

各向異性擴(kuò)散方程在進(jìn)行反向擴(kuò)散的時(shí)候可以對(duì)圖像增強(qiáng),但是這個(gè)二階非線性方程是不穩(wěn)定的,不能保證問(wèn)題的解是唯一的。針對(duì)此,Rudin L I等提出了一個(gè)穩(wěn)定的雙曲型偏微分方程,用來(lái)增強(qiáng)圖像,這就是經(jīng)典的震動(dòng)濾波模型[10]:

ut=-sign(uη η)|u|

(5)

其中,sign(uη η)是符號(hào)函數(shù),uη η表示圖像沿梯度方向的二階導(dǎo)數(shù)。

在含有噪聲的圖像中,震動(dòng)濾波模型式(5)不能區(qū)分邊緣和噪聲,從而震動(dòng)濾波在增強(qiáng)邊緣的同時(shí)也增強(qiáng)了噪聲。采用圖像的二階導(dǎo)數(shù)與低通濾波作卷積,能夠克服震動(dòng)濾波這一缺陷。如采用圖像的二階導(dǎo)數(shù)與高斯濾波作卷積,得到如下模型[11]：

ut=-sign(Gσ*uη η)|u|

(6)

2 各向異性擴(kuò)散模型

經(jīng)典的Gauss濾波為各向同性低通濾波,在去噪的同時(shí)會(huì)模糊圖像的邊緣。Perona P和Malik J首先提出了基于PDE的各向異性擴(kuò)散模型(PM模型)[12]：

(7)

其中,擴(kuò)散速度函數(shù)g(|u|)是一個(gè)單調(diào)非增的函數(shù),在梯度模值較大的邊緣,擴(kuò)散速度較小,可以保護(hù)邊緣；在平滑區(qū)域,擴(kuò)散速度較大。常用的擴(kuò)散速度函數(shù)如下：

(8)

其中,k為閾值。

在PM模型中,擴(kuò)散速度函數(shù)g(|u|)對(duì)噪聲很敏感。為了得到更好的去噪效果,近年來(lái)很多學(xué)者都在對(duì)擴(kuò)散速度函數(shù)g(·)進(jìn)行改進(jìn)[13-15]。

3 震動(dòng)耦合擴(kuò)散去噪模型

3.1 BZN的定義

文獻(xiàn)[13-15]中改進(jìn)的擴(kuò)散速度函數(shù)g(·)依賴于圖像的梯度模值|u|,從而改進(jìn)后的擴(kuò)散速度函數(shù)在保護(hù)邊緣的同時(shí)減少了對(duì)噪聲的擴(kuò)散強(qiáng)度,不能有效地去噪。在邊緣、噪聲和平滑區(qū)域的一個(gè)小鄰域內(nèi),由類似成分組成的像素的個(gè)數(shù)有很大的區(qū)別,比如對(duì)于噪聲來(lái)說(shuō),在它的一個(gè)小鄰域內(nèi),由類似成分組成的像素的數(shù)目是最小的。設(shè)BZN表示一個(gè)小鄰域內(nèi)由類似成分組成的像素的個(gè)數(shù),從而采用BZN可以區(qū)分邊緣、噪聲和平滑區(qū)域[9]。由于p點(diǎn)的像素值u(p)的不確定性,文獻(xiàn)[9]中計(jì)算BZN的算法會(huì)直接影響區(qū)域的劃分。因此，對(duì)計(jì)算BZN的算法進(jìn)行改進(jìn),定義如下:

定義 令p=(i,j)表示一個(gè)像素點(diǎn)的坐標(biāo),Mp(w)={(k,l)∶|k-i|≤w,|l-j|≤w}表示以p點(diǎn)為中心的(2w+1)×(2w+1)(w>0)的鄰域內(nèi)的像素點(diǎn),m表示所有對(duì)應(yīng)于Mp(w)的像素值的均值。?q∈Mp(w)的灰度強(qiáng)度Np(q)可分為兩類：

(9)

其中,T為預(yù)定閾值(T>0),u(q)表示q點(diǎn)的像素值。那么像素點(diǎn)p的BZN被定義為：

(10)

通過(guò)分析,可以得到像素點(diǎn)Mp(w)的個(gè)數(shù)N=(2w+1)2,BZNp(u,w,T)∈[0,N]。從BZN的定義中可以得到BZN有以下性質(zhì):

(1)噪聲像素點(diǎn)的BZN最小,接近于0;

(2)邊緣像素點(diǎn)的BZN接近中間值;

(3)內(nèi)部像素點(diǎn)的BZN最大,接近于N。

3.2 擴(kuò)散速度函數(shù)g(·)

g(BZNp(u,w,T))=

(11)

其中,BZNp(u,w,T)∈[0,N],g(BZNp(u,w,T))∈[0,1],如圖1所示。

圖1 擴(kuò)散速度隨BZN的變化情況

根據(jù)擴(kuò)散速度隨BZN的變化情況可知,擴(kuò)散速度函數(shù)在邊緣(BZN接近中間值)的擴(kuò)散速度最小,趨于零;在平滑區(qū)域(BZN最大)和噪聲點(diǎn)(BZN最小)的擴(kuò)散速度都比較大。這樣,在去噪的同時(shí)可以有效地保護(hù)邊緣。

3.3 控制保真項(xiàng)參數(shù)λ(·)的選取

結(jié)合BZN的性質(zhì),定義控制保真項(xiàng)參數(shù)λ(·)如下：

(12)

其中,λ(BZNp(u,w,T))∈[0,0.5]。

λ(BZNp(u,w,T)) 是單調(diào)遞增函數(shù),在噪聲點(diǎn)的值最小,在平滑區(qū)域的值最大,可以保證去噪后的圖像與原始圖像差異不大。

3.4 模型的提出

乘性噪聲的去除問(wèn)題就是從觀測(cè)圖像中恢復(fù)出原始圖像。本文考慮如下乘性噪聲模型：

u0(x)=u(x)v(x)

(13)

其中,x=(x1,x2)為圖像所在區(qū)域Ω內(nèi)的點(diǎn),u0、u、v分別表示觀測(cè)圖像、待恢復(fù)圖像以及噪聲。假定噪聲服從均值為1的Gamma分布，令f=logu0(x),z=logu(x),n=logv(x),式(13)可轉(zhuǎn)化為：

f=z+n

(14)

結(jié)合式(6)、式(7)和式(11),提出去噪模型如下：

zt=-(1-α)sign(Gσ*zη η)|z|+

αdiv(g(BZNp(u,w,T)z)+λ(logu0-z)

(15)

該模型的初始條件為：z0=logu0。

式(15)中,第一項(xiàng)是震動(dòng)濾波項(xiàng);第二項(xiàng)是擴(kuò)散項(xiàng);第三項(xiàng)是保真項(xiàng);α是尺度參數(shù);λ是控制保真項(xiàng)參數(shù),其取值按照式(12)計(jì)算。

4 模型求解

4.1 模型離散化

假定初始離散圖像u0的大小為m×m,令Δt為時(shí)間步長(zhǎng),空間采樣間隔Δx=Δy=h=1,w=2,T=30。在噪聲點(diǎn),α取值較大,接近于1;在邊緣處,α取值較小,接近于0;在平滑區(qū)域,α=0.5。對(duì)式(15)采用有限差分格式進(jìn)行離散化求解,其離散形式如下：

(16)

其中：

Dxx(zi,j)=(zi+1,j-2zi,j+zi-1,j)/h2；

Dyy(zi,j)=(zi,j+1-2zi,j+zi,j-1)/h2；

Dxy(zi,j)=(zi+1,j+1+zi-1,j-1-zi+1,j-1-zi-1,j+1)/4h2；

4.2 算法實(shí)現(xiàn)

經(jīng)過(guò)上述分析，給出震動(dòng)耦合擴(kuò)散的去噪算法如下：

Step1:設(shè)u0、u、v分別表示觀測(cè)圖像、待恢復(fù)圖像、噪聲。迭代次數(shù)為K,時(shí)間步長(zhǎng)Δt=0.2,置n=0,z0=log(u0),err=10-3;

Step2:用式(10)計(jì)算BZN;

Step3:用式(11)計(jì)算擴(kuò)散速度函數(shù)g(·);

Step4:用式(12)計(jì)算控制保真項(xiàng)參數(shù)λ(·);

Step5:用式(16)計(jì)算zn+1;

4.3 實(shí)驗(yàn)結(jié)果

利用上述算法對(duì)Lena圖像用Matlab編程進(jìn)行仿真實(shí)驗(yàn)。實(shí)驗(yàn)效果如圖2所示，圖中：(a)為原始圖像;(b)為加入服從均值為1的Gamma分布的含噪圖像;(c)為采用JY模型去噪后的圖像;(d)為采用HNW模型去噪后的圖像;(e)為本文去噪模型的圖像。通過(guò)比較這幾種算法的去噪效果圖可以看出:HNW模型的去噪效果要比JY模型的好些,但是會(huì)出現(xiàn)“階梯效應(yīng)”,而且會(huì)使邊緣模糊;本文去噪算法采用BZN來(lái)區(qū)分噪聲點(diǎn)、邊緣區(qū)域以及內(nèi)部平滑區(qū)域,可以在去噪的同時(shí)有效地保護(hù)邊緣。

從數(shù)據(jù)上比較算法的優(yōu)越性,采用峰值信噪比(PSNR)和均方誤差(MSE)來(lái)評(píng)價(jià)去噪后圖像效果。PSNR和MSE的計(jì)算公式如下,計(jì)算結(jié)果如表1所示。

其中,u0(m,n)表示原始圖像,大小為M×N;z(m,n)為去噪后的恢復(fù)圖像。

表1 不同模型去噪后的峰值信噪比(PSNR)和均方誤差(MSE)

5 結(jié) 語(yǔ)

本文分析了經(jīng)典的擴(kuò)散去噪模型的優(yōu)缺點(diǎn),將震動(dòng)濾波引入去除乘性噪聲模型。去噪模型由擴(kuò)散項(xiàng)、震動(dòng)濾波項(xiàng)和保真項(xiàng)構(gòu)成,文中把擴(kuò)散速度函數(shù)g(·)改進(jìn)為一個(gè)多項(xiàng)式函數(shù),也對(duì)控制保真項(xiàng)參數(shù)λ進(jìn)行了改進(jìn)。該去噪算法采用BZN來(lái)區(qū)分噪聲、邊緣以及平滑區(qū)域,在去噪的同時(shí)能有效地保護(hù)邊緣,同時(shí)也能抑制“階梯效應(yīng)”的產(chǎn)生,是一種有效的去噪算法。

參考文獻(xiàn)：

[1] 吳亮,胡云安.遙感圖像自動(dòng)道路提取方法綜述[J].自動(dòng)化學(xué)報(bào),2010,36(7):912-922.

[2] 霍春雷,程健,盧漢清,等.基于多尺度融合的對(duì)象級(jí)變化檢測(cè)新方法[J].自動(dòng)化學(xué)報(bào),2008,34(3):251-257.

[3] Rudin L I,Osher S,Fatemi E.Nonlinear Total Variation Based Noise Removal Algorithms[J].Physica D,1992,60:259-268.

[4] Rudin L,Lions P L,Osher S.Multiplicative Denoising and Debulurring:Theory and Algori-thms[C]//Geometric Level Set Methods in Imaging,Vision and Graphics.New York:Springer,2003：103-119.

[5] Aubert G,Aujol J.A Nonconvex Model to Remove Multiplicative Noise[C]//Proceedings of the 1st International Conference on Scale Space and Variational Methods in Computer Vision.Ischia,Italy:Springer,2007：68-79.

[6] Aubert G,Aujol J.A Variational Approach to Removing Multiplicative Noise[J].SIAM Journal on Applied Mathematics,2008,68(4):925-946.

[7] Jin Z M,Yang X P.Analysis of A New Variation Model for Multiplicative Noise Removal[J].Journal of Mathematical Analysis and Applications,2010,362(2):415-426.

[8] Huang Y M,Ng M K,Wen Y W.A New Total Variation Method for Multiplicative Noise Removal[J].SIAM Journal on Imaging Sciences,2009,2(1):20-40.

[9] Wu J,Chen T.PDE-Based Random-valued Impulse Noise Removal Based on New Class of Controlling Functions[J].IEEE Transactions on Image Processing,2011, 20(9):2428-2438.

[10] Osher S,Rudin L I.Feature-oriented Image Enhancement Using Shock Filters[J].SIAM Journal on Numerical Analysis,1990,27:919-40.

[11] Alvarez L,Mazorra L.Signal and Image Restoration Using Shock Filters and Anisotro-pic Diffusion[J].SIAM Journal on Numerical Analysis,1994,31:590-605.

[12] Perona P,Malik J.Scale-space and Edge Detection Using Anisotropic Diffusion[J].IEEE Transactions on Pattern Analysis and Machine Intelligence,1990,12(7):629-639.

[13] Catt F,Lions P-L,Morel J-M,et al.Image Selective Smoothing and Edge-detection by Nonlinear Diffusion[J].SIAM Journal on Numerical Analysis,1992,29(1):182-193.

[14] You Y L,Kaveh M.Fourth-order Partial Differential Equations for Noise Removal[J].IEEE Transactions on Image Processing,2000,9(10):1723-1730.

[15] Alvarez L,Lions P-L,Morel J-M.Image Selective Smoothing and Edge-detection by No-nlinear Diffusion.II[J].SIAM Journal on Numerical Analysis,1992,29(3):845-866.