劉洪琛,劉朝霞,張 龍
中央民族大學(xué) 理學(xué)院,北京100081
圖像去噪是圖像處理研究領(lǐng)域中經(jīng)常研究的問題[1]。為了在去除噪聲的同時更好地保持圖像的清晰度、紋理特征等有用的信息,大量的去噪方法相繼出現(xiàn)[2]。其中,1992 年Rudin 等提出的全變差正則化[3](Rudin Osher Fatemi,ROF)模型是迄今最為成功的圖像復(fù)原模型之一。全變差正則化模型在圖像去噪領(lǐng)域得到了廣泛的研究與應(yīng)用,但該模型有一個很大的缺點(diǎn),易在相對平滑區(qū)域產(chǎn)生“階梯效應(yīng)”[4];基于L2范數(shù)的調(diào)和模型能有效地在平滑區(qū)域去除噪聲,但處理圖像的邊緣易模糊[5]。2007 年,艾必剛等[6]利用分?jǐn)?shù)階微積分的定義,推導(dǎo)出分?jǐn)?shù)階微分梯度算子的頻率濾波器。該濾波器對圖像增強(qiáng)的效果較好,在增強(qiáng)細(xì)節(jié)的同時,能夠有效去除噪聲,但同時也增加了計算復(fù)雜程度;2016 年,石玉英等[7]提出了基于L1和L2混合保真項的圖像恢復(fù)模型,充分發(fā)揮了L1和L2保真項的優(yōu)點(diǎn),有效地去除了高斯-椒鹽混合噪聲并保留圖像的有效信息,但還是存在較大程度的階梯效應(yīng);2018 年,馬紅強(qiáng)等[8]提出了基于改進(jìn)棧式稀疏去噪自編碼器的圖像去噪模型,降噪效果良好。但是需要數(shù)據(jù)集訓(xùn)練模型,增加了訓(xùn)練成本;同年,尹芳等[9]提出了新的小波域Wener濾波器和Perona-Malik融合去噪算法,具有良好的去噪和抑制偽吉布斯的能力,但是對于計算時間的控制有待提升。
常見噪聲有高斯噪聲、脈沖噪聲、泊松噪聲等。高斯噪聲是指它的概率密度函數(shù)服從高斯分布(即正態(tài)分布)的一類噪聲,包括起伏噪聲、熱噪聲等;泊松噪聲是一種噪聲分布符合泊松分布模型的與信號有關(guān)的噪聲,包含在各種應(yīng)用中,如射線照相、熒光顯微鏡、正電子發(fā)射斷層掃描(Positron Emission Tomography,PET)、光學(xué)納米檢測和天文成像應(yīng)用[5]?;謴?fù)被泊松噪聲破壞的模糊圖像是困難的。工業(yè)領(lǐng)域中,X射線圖像檢測作為一種無損、無接觸和高分辨率的缺陷檢測方法,X 射線在成像過程中,會受到不同類型噪聲污染。X射線圖像的噪聲不能簡單用單一的噪聲模型來表示,它表現(xiàn)為高斯-泊松混合噪聲,因此,去除高斯-泊松混合噪聲具有很廣泛的應(yīng)用場景[10]。
為了有效抑制高斯-泊松混合噪聲,獲得更好的視覺效果,本文提出基于L2和Kullback-Leibler(KL)散度作為保真項(KL保真項)與兩個正則項相結(jié)合的圖像恢復(fù)的變分模型,并用增廣拉格朗日算法進(jìn)行數(shù)值實(shí)現(xiàn)。實(shí)驗(yàn)結(jié)果表明,針對目前的幾種去除高斯-泊松混合噪聲方法的局限性,所提出的新模型不僅保留了圖像的邊緣信息,減少了圖像的“階梯效應(yīng)”,增強(qiáng)了圖像的視覺效果,同時有效降低了均方根誤差,得到了更高的峰值信噪比,CPU 的計算時間更短,并且不需要針對特定的訓(xùn)練集進(jìn)行訓(xùn)練,提高了算法的效率,說明了算法的必要性和優(yōu)越性。
不失一般性,用M×N 矩陣表示一幅灰度圖像,V=RM×N,離散的梯度算子?:V →Q,其中Q=V×V ,f=Ku+n,f ∈V 是觀察到的圖像,u ∈V 是真實(shí)圖像,n ∈V 為噪聲算子,K:V →V 為模糊算子,M ,N 是圖像大小。
1977年,Tikhonov等[11]提出調(diào)和模型:
此模型能在在平滑區(qū)有效地去除噪聲,但是該模型朝各個方向的擴(kuò)散能力相同,容易模糊圖像的邊緣輪廓等細(xì)節(jié)信息。
基于L2保真項的全變差圖像恢復(fù)模型(ROF 模型),是圖像去噪中很成功的模型之一[3]。ROF模型為:
在模型中,保真項體現(xiàn)出所恢復(fù)的圖像與原圖之間的估計殘差,估計殘差越大,保真項的值越大,反之亦然。為保證函數(shù)為凸函數(shù),保真項的指數(shù)取值范圍是[1,2],取值2 時即為L2保真項。在各向同性擴(kuò)散的過程中L2保真項針對不同的估計殘差賦予其相應(yīng)的權(quán)重[12-13],保證了在去除高斯噪聲時的魯棒性[14]。
此模型能很好保存圖像的邊緣特征,但是在處理圖像相對光滑的區(qū)域部分會產(chǎn)生“階梯效應(yīng)”。
根據(jù)泊松噪聲分布的特點(diǎn),在文獻(xiàn)[5,15]中使用Kullback-Leibler散度作為保真項(KL保真項)的全變差圖像恢復(fù)模型(Total Variation Kullback-Leibler,TV-KL模型)用于解決這個問題。TV-KL模型為:
TV-KL 模型比標(biāo)準(zhǔn)的期望最大化(EM)重建(其中 不應(yīng)用TV 正則項)[16]表現(xiàn)得更穩(wěn)定和穩(wěn)健,并且比ROF 模型去除Poisson噪聲的效果更好,但也會產(chǎn)生“階梯效應(yīng)”。
由于‖ ‖?uL1會產(chǎn)生階梯效應(yīng),為了克服這個問題,文獻(xiàn)[17]中提出了一個改進(jìn)的模型——MS(Mumford-Shah)模型。文獻(xiàn)[7]中提出了一個求解MS模型的近似模型:
使用KL保真項的全變差圖像恢復(fù)模型能有效地去除Poisson 噪聲。基于模型式(4),相應(yīng)去除Poisson 噪聲的模型為:
文獻(xiàn)[10]中為了去除高斯-泊松混合噪聲提出了基于混合Kullback-Leibler散度作為保真項和L2保真項的全變差模型,該模型為:
此模型中L2保真項能有效去除高斯噪聲,基于混合KL保真項能有效去除泊松噪聲,增強(qiáng)了去噪性能。
綜合考慮調(diào)和模型、TV-KL模型、改進(jìn)MS模型、保真項混合模型的優(yōu)點(diǎn),并且受文獻(xiàn)[6,18]的啟發(fā),提出了對高斯-泊松混合噪聲去噪的圖像恢復(fù)模型:
其中,α ≥0 是平衡正則項和保真項的參數(shù),r1≥0,r2≥0用來調(diào)節(jié)高斯噪聲和泊松噪聲在總數(shù)據(jù)保真項中的比例,依據(jù)文獻(xiàn)[10],比例系數(shù)可根據(jù)對噪聲的先驗(yàn)知識或采用統(tǒng)計方法來確定。在此模型中,保真項。
從式(7)可知,當(dāng)r1>0,r2>0,α=0,調(diào)整為模型(6),其可以較好地去除高斯-泊松混合噪聲,并且保持邊緣等細(xì)節(jié)信息[19];當(dāng)r1=0,r2>0,α >0 時,調(diào)整為模型(5),可較好地去除高斯噪聲,并且在一定程度上克服階梯效應(yīng)。
模型(7)和模型(5)不同的是,保真項采用L2和KL保真項的混合項,提高了去除高斯噪聲的能力,一定程度上有效避免了階梯效應(yīng);模型(7)和模型(6)不同的是,正則項采用的∫Ω||?u2和∫Ω||?u 的混合正則項,提高了在平滑區(qū)去除噪聲的能力。因此,通過融合思想而提出的新模型充分體現(xiàn)了KL保真項和L2保真項的特點(diǎn),集合了模型(5)(6)中的優(yōu)點(diǎn),可以更為有效地去除高斯-泊松混合噪聲,并且優(yōu)化了模型(5)(6)的不足,在一定程度上避免了階梯效應(yīng),保持了邊緣細(xì)節(jié)信息。在本文中,將提出的模型(7)記為新模型。
增廣拉格朗日方法[10]的基本思想是引入拉格朗日乘子和懲罰項將帶約束的問題轉(zhuǎn)換成無約束問題進(jìn)行求解。本章使用增廣拉格朗日方法求解新模型。由于?u 的不可微性和非線性,因此需要引入新的變量[20]。
下面引入兩個新的變量p ∈Q 和z ∈V,上述式(7)等價于解下面這個帶約束泛函的極值問題:
使用增廣拉格朗日方法將上述帶約束問題轉(zhuǎn)換成下列不帶約束的問題,定義增廣拉格朗日泛函如下:
其中拉格朗日乘子λp∈Q,λz∈V ,正則化參數(shù)rp>0,rz>0。
求解問題式(3)轉(zhuǎn)化為求解增廣拉格朗日泛函的鞍點(diǎn)問題:
仿照文獻(xiàn)[20],可以類似地證明鞍點(diǎn)問題(10)至少有一個解。對于上述鞍點(diǎn)問題,使用迭代算法解決,算法如下:
(1)初始化u0=0,p0=0,z0=0 和
(3)更新拉格朗日乘子:
u-子問題。給定p,z,有:
p-子問題。給定u,z,有:
z-子問題。給定u,p,有:
下面使用交替求極小的算法求解上述三個子問題:
(1)解決u-子問題(12)
由邊界條件和極小化條件[19-20],式(12)對應(yīng)的Euler-Lagrange方程為:
使用快速傅里葉變換(Fast Fourier Transformation,F(xiàn)FT)求解上述線性方程,定義F(u) 為u 的傅里葉變換,通過傅里葉變換可計算u:
(2)解決p-子問題(13)式(13)簡化為:
由邊界條件和極小化條件[19-20],式(14)對應(yīng)的Euler-Lagrange方程為:
參考文獻(xiàn)[18]中對模型(15)算法收斂性的證明,以及結(jié)合文獻(xiàn)[22]中對增廣拉格朗日算法應(yīng)用于TV模型的收斂準(zhǔn)確性證明,可以證明提出的模型(7)在增廣拉格朗日算法中是嚴(yán)格收斂的,準(zhǔn)確性較高。
圖1 測試圖像
為了驗(yàn)證新模型的有效性,仿真實(shí)驗(yàn)均在AMD A8-3500M APU1.5 GHz CPU,24 GB RAM,Windows 7操作系統(tǒng)且運(yùn)行環(huán)境為MATLAB 8.4.0.150421(R2014b)實(shí)驗(yàn)環(huán)境下實(shí)現(xiàn)。首先,給出模型參數(shù)的設(shè)置,參數(shù)拉格朗日乘子λp,λz是在算法迭代過程中引入的,它的選取對去噪結(jié)果圖像的影響非常小,均取值為0,參考文獻(xiàn)[10]、[23]與[24]中選取參數(shù)的方法,根據(jù)對噪聲的先驗(yàn)知識,通過重復(fù)實(shí)驗(yàn),選取r1=20 ,r2=10 ,α=2,rp=20,rz=4。在此實(shí)驗(yàn)中,對測試圖像使用Matlab代碼“imnoise(u,‘gaussian’,0,r)”添加方差為r 的高斯噪聲(Gaussian Noise,GN),使用Matlab代碼“imnoise(u,‘Poisson’)(根據(jù)圖像數(shù)據(jù)數(shù)據(jù)生成泊松噪聲,而不是向數(shù)據(jù)添加人為噪聲)添加泊松噪聲(Poisson Noise,PN),其中,u 是原始圖像,為了簡單起見,G(0.01)表示對測試圖像添加方差為0.01的高斯噪聲;P表示對測試圖像添加泊松噪聲:G(0.02)+P表示對測試圖像添加方差為0.02的高斯噪聲和泊松噪聲。圖1中的圖像是從測試圖像中選取的6 幅圖像,Pepper 圖像、Man 圖像、Cameraman 圖像、Tablet 圖像、Lena 圖像大小均為256×256,Brain 圖像大小為709×591。為了客觀評價圖像去噪模型的性能,引入了峰值信噪比(Peak Signal to Noise Ratio,PNSR)和結(jié)構(gòu)相似度(Structural Similarity Index,SSIM)兩個客觀評價指標(biāo)。PSNR越大且SSIM越大,算法的去噪效果越好。其定義分別如下:
其中,f 是去噪后圖像,u 是原始圖像,M ,N 是圖像大小,μf和μu分別是f 和u 的均值,σf和σu分別是f和u 的方差。實(shí)驗(yàn)中為了簡介,所選取的對比模型有使用Kullback-Leibler 散度作為保真項的全變差圖像恢復(fù)模型(TV-KL)、改進(jìn)MS模型(Mix Regularization Term,MRT),以及保真項混合模型(Mix Fidelity Term,MFT)。
表1和表2給出了5幅測試圖像在添加不同比例噪聲時的峰值信噪比(PSNR)和結(jié)構(gòu)相似度(SSIM)的比較,從表1 和表2 中可得,在添加方差較小的高斯-泊松混合噪聲,新模型和MFT 模型的PSNR 與SSIM 比較接近,并且高于TV-KL模型和MRT模型。但是,隨著添加方差增大的高斯-泊松混合噪聲時,四種模型對應(yīng)的PSNR都比較低,分析其原因,當(dāng)混合噪聲中高斯噪聲方差較大時,圖像中的邊緣、紋理等幾何結(jié)構(gòu)信息丟失很多。在添加高斯噪聲方差增大過程中,總體上新模型的PSNR 的值比其他模型大,SSIM 的值也大于其他模型,由于新模型結(jié)合了調(diào)和模型有效處理圖像平滑區(qū)域與TV-KL 模型有效保存圖像邊緣信息的優(yōu)勢,并且L2保真項能有效去除高斯噪聲,基于混合Kullback-Leibler散度作為保真項能有效去除泊松噪聲。
表1 去噪圖像PSNR對比 dB
表2 去噪圖像SSIM對比
圖2 比較了在添加不同高斯-泊松混合噪聲下Pepper 圖像的恢復(fù)結(jié)果。第一行到第四行,對測試圖像Pepper 分別添加G(0.01)+P、G(0.02)+P、G(0.03)+P、G(0.04)+P的混合噪聲。Pepper圖像的特點(diǎn)主要包含平滑和邊緣部分,紋理較少。根據(jù)圖像的視覺效果,針對添加不同方差的混合噪聲,MFT 模型的圖像恢復(fù)能力相對較弱,在添加混合噪聲中高斯噪聲方差較大時,去噪后得到的圖像比較模糊,同時圖像的結(jié)構(gòu)和邊緣信息也未能恢復(fù)。MRT模型對圖像平滑部分的去噪效果較好,不足的是當(dāng)圖像添加高斯噪聲方差增大時邊緣細(xì)節(jié)信息丟失較多。TV-KL模型去噪能力很強(qiáng),但是易產(chǎn)生“階梯效應(yīng)”。而提出的新模型不僅去噪能力很強(qiáng),圖像的結(jié)構(gòu)和邊緣信息得到有效的恢復(fù),而且在一定程度上程度上抑制了“階梯效應(yīng)”。
對于如Man 圖像、Cameraman 圖像這類圖像紋理、邊緣等幾何結(jié)構(gòu)較多的圖像,分別添加G(0.01)+P、G(0.02)+P、G(0.03)+P、G(0.04)+P 的混合噪聲。選擇其中一部分去噪后的圖像,如圖3和圖4所示,其中圖3的第一行添加G(0.01)+P 混合噪聲,第二行添加了G(0.03)+P 混合噪聲,整體而言,當(dāng)添加混合噪聲中高斯噪聲方差較小時,新模型不僅去除了大部分噪聲,輪廓更加簡潔,紋理信息得到更好的保持,有很好的視覺效果,并且有較高的評價指標(biāo)。雖然MRT 模型也有較高的評價指標(biāo),但是得到結(jié)果圖像過于光滑,產(chǎn)生了“階梯效應(yīng)”。
對于如Lena圖像、Tablet這類含有豐富紋理和邊緣信息的圖像,圖5~7 中,四種方法得到的圖像恢復(fù)結(jié)果和圖2基本一致。從局部放大圖像可看到,新模型去除高斯-泊松混合噪聲的效果更好,噪聲斑點(diǎn)較少(圖5(i)),MFT 模型雖然能在一定程度上去除高斯-泊松混合噪聲(圖5(h)),但是當(dāng)混合噪聲中高斯噪聲方差較大,MFT模型處理的結(jié)果圖像變得模糊,圖像分辨率降低,并且邊緣特征信息也不能被保持(圖6(h))。新模型可以有效地去除混合噪聲并且較好地恢復(fù)了圖像的紋理結(jié)構(gòu)信息(圖5(i)),較好地保留了圖像的邊緣信息,當(dāng)混合噪聲中高斯噪聲方差較大時,新模型處理的效果也比較好(圖6(i)和圖7(i))。MRT 模型去除混合噪聲效果也比較好(圖5(g)),但是當(dāng)混合噪聲中高斯噪聲方差較大時,MRT模型處理的結(jié)果圖過于光滑(圖6(g)和圖7(g))。
圖2 添加不同高斯-泊松混合噪聲的Pepper圖像恢復(fù)效果比較
圖3 添加G(0.01)+P,G(0.03)+P噪聲的Man圖像恢復(fù)效果比較
圖4 添加G(0.04)+P混合噪聲的Cameraman圖像恢復(fù)效果比較
圖5 添加G(0.01)+P混合噪聲的Lena圖像恢復(fù)圖像及局部放大圖像比較
圖6 添加G(0.03)+P混合噪聲的Lena圖像恢復(fù)圖像及局部放大圖像比較
圖7 添加G(0.02)+P混合噪聲的Tablet圖像恢復(fù)圖像及局部放大圖像比較
最后,通過CPU 時間來比較幾種模型的時間復(fù)雜度,表3統(tǒng)計了5幅圖像在添加G(0.01)+P、G(0.02)+P、G(0.03)+P、G(0.04)+P 的混合噪聲下的平均運(yùn)行時間。由表3可知,新模型運(yùn)行時間最短,其次是TV-KL,然后是MRT模型,MFT模型的運(yùn)行時間最長。因此,從去噪的運(yùn)行效率上來看,新模型相對于TV-KL 模型、MRT模型、MFT模型有比較強(qiáng)的優(yōu)勢。
表3 不同混合噪聲下圖像恢復(fù)運(yùn)行時間對比s
為了有效抑制高斯-泊松混合噪聲,獲得更好的視覺效果,本文在綜合考慮TV-KL模型、MRT模型和MFT的優(yōu)點(diǎn),提出基于L2和KL保真項與兩個正則項相結(jié)合的圖像恢復(fù)的變分模型,并用增廣拉格朗日方法對模型進(jìn)行化簡,化為三個子問題進(jìn)行快速求解。為了加速計算,使用FFT 來解決u 子問題。對添加不同高斯-泊松混合噪聲的測試圖像進(jìn)行數(shù)值模擬實(shí)驗(yàn),實(shí)驗(yàn)表明提出的新模型有效去除了混合噪聲,恢復(fù)了圖像的紋理結(jié)構(gòu)信息,圖像的邊緣特征被更好地保持,同時模型的算法也具有較高的執(zhí)行效率,并且一定程度上抑制了“階梯效應(yīng)”。另一方面,本文提供的測試圖片和實(shí)際生活需要處理的圖片還存在一定的差距,因此,在后續(xù)的研究工作中,會尋找更加高效的、實(shí)現(xiàn)效果更好的去除高斯-泊松混合噪聲的圖像恢復(fù)模型。