張宗良，顏永豐

(西北農(nóng)林科技大學(xué)信息工程學(xué)院，陜西楊凌712100)

0 引言

圖像去噪具有極高的理論和實(shí)用價(jià)值，近年來數(shù)字圖像的急劇增多更加突顯了圖像去噪的重要意義。傳統(tǒng)的去噪方法通常先假設(shè)圖像具有某種良好的性質(zhì)，如:具有較小的總變分［1-2］，或者在小波基中稀疏，然后通過最小化與該性質(zhì)對(duì)應(yīng)的目標(biāo)函數(shù)來去噪。

最近幾年利用圖像塊之間的相似性來去噪取得了極大的成功。文獻(xiàn)［3］最早研究了基于塊處理的圖像去噪，并提出了非局部均值 (NLM)算法。該算法通過加權(quán)累加非局部塊來估計(jì)當(dāng)前所求塊。權(quán)值由非局部塊與當(dāng)前所求塊的相似度決定。在非局部塊方法的基礎(chǔ)上，文獻(xiàn)［4］提出了一種全新的算法BM3D。該算法先用塊匹配算法找到那些與當(dāng)前塊相似的塊，再將這些相似塊組成一個(gè)柱體。由于該柱體具有極大的冗余性，對(duì)該柱體進(jìn)行三維離散余弦變換、閾值截?cái)?、余弦逆變換后便可去噪。文獻(xiàn)［5］的算法與BM3D算法的不同之處是它對(duì)柱體進(jìn)行的是PCA變換。

上述基于塊相似的去噪簡單有效，但是只有少量的文章［6-9］給出其背后的數(shù)學(xué)解釋。總結(jié)這些文章我們知道無噪圖像中的相似塊會(huì)在圖像塊空間中形成低維流形，而在帶噪圖像中低維流形受到噪聲污染而變形。復(fù)原這些低維流形便可去噪。本文據(jù)此提出了一個(gè)多步迭代的去噪算法。與大多數(shù)類BM3D算法的兩步迭代相比，本文算法的多步迭代更具有一般性。

1 圖像塊空間中的低維流形

本文只考慮灰度數(shù)字圖像。一幅圖像是從平面區(qū)域到像素域的映射u:D→［0，255］ R，其中D Z2是圖像所在的平面區(qū)域。ux，y∈［0，255］是位置 (x，y)∈D上的像素值，［0，255］是像素的取值范圍即“像素域”。記位置(x，y)的鄰域?yàn)镹x，y。本文涉及到的都是長方形或方形的鄰域。有兩種常見的對(duì)(x，y)取鄰域的方法，一種是(x，y)在鄰域的左上角，另一種是(x，y)在鄰域的中心。例如若采用第一種取法，則(x，y)的大小為I×J的鄰域Nx，y={(i，j)}，i∈［x，x+I－1］∩Z ，j∈［y，y+J－1］∩Z 。圖像塊是尺寸遠(yuǎn)小于原圖像的子圖像，即原圖像上的一個(gè)小塊。位置(x，y)上的圖像塊是指原圖像在鄰域Nx，y上的一幅子圖像px，y=uNx，y。這里有個(gè)細(xì)節(jié)問題，當(dāng) (x，y)靠近圖像的邊界時(shí)，鄰域Nx，y有可能超出圖像所在的平面區(qū)域即Nx，yD。此時(shí)可對(duì)圖像進(jìn)行擴(kuò)展，本文采用的擴(kuò)展方法是使擴(kuò)展部分與原圖像沿邊界對(duì)稱。接下來我們進(jìn)一步的記由圖像塊px，y的像素按行列順序排列而成的向量為塊向量px，y，稱塊向量所在的空間為塊空間。

下面我們考察圖像塊向量在塊空間中的分布是什么樣的。真實(shí)自然的圖像具有很多良好的性質(zhì)。例如圖像中含有大量相同或相似的塊，這些相同或相似的塊向量將分布在塊空間中一個(gè)半徑很小的球中。如圖1中的標(biāo)為1、2的塊分別與標(biāo)為R1、R2的塊相似。還有些圖像塊由于光照強(qiáng)度的漸變等原因而有規(guī)律的變化，這些圖像塊向量也將按一定的規(guī)律分布在塊空間中。由此可見，特定區(qū)域內(nèi)的所有圖像塊向量將在塊空間中形成若干個(gè)流形［6-8］。這里的“流形”可以理解為是由相似塊或有規(guī)律的塊組成的一個(gè)有規(guī)律可循的組織或結(jié)構(gòu)。圖2是一個(gè)二維塊分布示意圖，在該圖中，塊空間是二維的即每個(gè)塊只含二個(gè)像素。從左圖可看到若干條直線 (近似)，每條直線可看成一個(gè)低維的流形。右圖是對(duì)應(yīng)的受到噪聲污染后的流形，可以看到噪聲的影響直觀上是使直線變粗?？梢?，只要我們能從受噪聲污染的流形 (圖1b和圖2b)中恢復(fù)出原本低維的流形(圖1a和圖2a)，再將恢復(fù)后流形中的點(diǎn)即圖像塊放回到該圖像塊原來所在的位置便能達(dá)到去噪效果。

接下來看看真實(shí)圖像中的情況。同樣在塊大小等于2(即一個(gè)塊只含左右兩個(gè)像素)的情況下畫出塊分布圖。圖4(a)、(b)、(c)分別是圖3中實(shí)線條框出的甲、乙、丙三個(gè)區(qū)域及其中的圖像塊在塊空間中的分布圖。左邊的是無噪的情況，右邊是加入方差為5的高斯白噪聲的情況。從圖4可以看出塊向量在塊空間中的分布確實(shí)有一定的規(guī)律，或聚成小圓，或近似排成直線。從圖4還能明顯看出噪聲對(duì)分布圖的影響，或使小圓變大，或使直線變粗。由此驗(yàn)證了低維流形假設(shè)的正確性。以上是圖像塊大小等于2時(shí)的情況，我們相信圖像塊更大時(shí)在塊空間也會(huì)呈現(xiàn)出一定的規(guī)律。

2 算法描述

輸入:含噪圖像u，塊尺寸參數(shù)n，搜索范圍尺寸參數(shù)M，距離閾值T，高斯函數(shù)平滑度參數(shù)h1、h2，間隔參數(shù)step，初始化偏移圖像ur=u。

第一步:塊匹配。對(duì)偏移圖像中的每個(gè)圖像塊 (記為當(dāng)前塊)，在一定的搜索范圍內(nèi)找到與之相似的塊。即視那些與當(dāng)前塊距離小于一定閾值的塊為相似的塊。形式化的描述塊匹配如下。簡潔起見，以下用k表示位置(x，y)。令分別表示帶噪圖像u和偏移圖像ur中所有塊的集合，其中分別是將帶噪圖像和偏移圖像中第k塊 (圖像塊以k為中心，大小為(2n+1)×(2n+1))按列堆疊起來的向量，K是塊的總個(gè)數(shù)。特別的我們稱為位置k上的參考?jí)K，簡稱參考?jí)K。接下來對(duì)每個(gè)位置k，搜索所有與相似的塊q。理想情況下，搜索范圍應(yīng)該是整幅帶噪圖像即Ω。但是這樣的搜索代價(jià)太大，本文將搜索范圍縮小到以位置k為中心的(2M+1)×(2M+1)大小的區(qū)域，并記這個(gè)區(qū)域內(nèi)塊的集合為Ωk。易知Ωk為Ω的子集。從而所有與相似的塊的集合，即所在的流形可定義為:

圖4 圖3中甲乙丙區(qū)域及其中的塊分布

其中T是給定的距離閾值，距離小于該閾值的兩塊被認(rèn)為是相似的。是塊p與塊q的距離。還有一點(diǎn)需要說明的是，在實(shí)際中不必對(duì)每個(gè)塊都進(jìn)行塊匹配，也就是可以隔幾行隔幾列取塊進(jìn)行塊匹配。本文中行間隔和列間隔取相等的值，記為step。

第三步:加權(quán)累加?，F(xiàn)在需要將復(fù)原后低維流形中的圖像塊放回到它原來的位置。由于搜索范圍的重疊和圖像塊的重疊，這一步將需要兩次的加權(quán)累加。在前面一步中我們得到了去噪后的低維流形對(duì)應(yīng)的低秩矩陣Ak，現(xiàn)在將Ak中的塊加上偏移量放回到原來的位置。放回的時(shí)候還有講究。由于在第一步中搜索范圍Ωk有重疊，因此將會(huì)出現(xiàn)同一位置上的塊被多次處理的情況，即該位置上會(huì)有多個(gè)處理結(jié)果。我們需要對(duì)這些結(jié)果加權(quán)累加以得到該位置上圖像塊的去噪結(jié)果。令給定位置k上的塊被處理過J次，第j次的結(jié)果為則第j次對(duì)應(yīng)的權(quán)為

其中參數(shù)h1用來控制高斯函數(shù)的平滑度。得到權(quán)值后，進(jìn)行累加，則位置k上的塊的去噪結(jié)果k為

類似的，參數(shù)h2用來控制高斯函數(shù)的平滑度表示取絕對(duì)值。最終去噪圖像在位置k處的像素值為

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

本文采用兩種客觀標(biāo)準(zhǔn)來評(píng)價(jià)去噪效果。第一種是峰值信噪比

式中:u0——原始無噪圖像——無噪圖像的像素總數(shù)?！ピ牒蟮膱D像，·F——Frobenius范數(shù)。PSNR越大去噪效果越好。另一種是SSIM指數(shù)［13］，同樣越大越好。為簡單起見，在實(shí)驗(yàn)中固定參數(shù)n=1，M=10，h1=1，h2=25，step=3，T=2σ，σ為所加噪聲的方差。本文的噪聲都是與原始無噪圖像不相關(guān)的加性高斯白噪聲。

圖5展示了對(duì)噪聲方差為20的House圖像的去噪過程，圖像長寬都為256，進(jìn)行了4次迭代。從該圖可以看出經(jīng)過迭代后去噪效果逐漸提高。

圖6是對(duì)噪聲方差為20的House和Lena含噪圖像的去噪對(duì)比圖，圖像長寬都為256。其中圖6(c)、(f)是本文算法迭代5次時(shí)的結(jié)果。由圖可以看出本文所提方法確實(shí)能達(dá)到很好的去噪效果，視覺效果已與當(dāng)前去噪效果最好的算法 BM3D［14］差別不大。

表1和表2分別是對(duì)加入不同方差噪聲后的House和Lena圖像的去噪結(jié)果，圖像長寬都為128，迭代次數(shù)不一。表中數(shù)據(jù)上一行是PSNR值，下一行是SSIM值。初步結(jié)果表明本文所提算法去噪效果良好，與BM3D算法差距較小。值得一提的是前面所用的實(shí)驗(yàn)參數(shù)還不是最優(yōu)的，如何選擇好的參數(shù)是我們進(jìn)一步研究的內(nèi)容。

表1 圖像“House”去噪

表2 圖像"Lena"去噪

4 討論

前面提到過，傳統(tǒng)方法通過最優(yōu)化含約束項(xiàng)的目標(biāo)函數(shù)來實(shí)現(xiàn)去噪。本文的方法也有類似的形式，約束是無噪的圖像塊空間中具有低維流形。令分別是圖像塊空間中無噪流形和被噪聲污染流形的集合，其中M是流形的個(gè)數(shù)。則要得到降噪的圖像塊可通過解下面的最優(yōu)化問題

其中Dm和Am分別是重新確定坐標(biāo)原點(diǎn)后對(duì)應(yīng)的矩陣，R(Am)是Am的秩，p指定范數(shù)的類型，λm是調(diào)節(jié)參數(shù)。進(jìn)一步的，我認(rèn)為NLM和BM3D算法的成功可以看成是低維流形假設(shè)的特殊情況:復(fù)原塊空間中的0維流形，即他們?cè)噲D將受噪聲污染變成的小球復(fù)原成一個(gè)點(diǎn)。

5 結(jié)束語

本文進(jìn)一步分析了基于塊處理的圖像去噪。前人的低維流形假設(shè)指出:真實(shí)無噪圖像在塊空間中普遍存在低維流形，而噪聲的污染使這些低維流形變形。本文用實(shí)際例子給出了低維流形的幾何解釋。并根據(jù)該假設(shè)給出了通過復(fù)原塊空間中的低維流形來去噪的一個(gè)迭代算法。在算法上給出了偏移圖像的概念。初步實(shí)驗(yàn)表明本算法能達(dá)到良好的去噪效果，與當(dāng)今最好的去噪算法有可比性。本算法還有很大改進(jìn)的空間，未來的工作主要是尋找更好的塊匹配和低維流形恢復(fù)方法，考慮在其他領(lǐng)域的應(yīng)用如圖像的超分辨率重建等。

［1］Beck A，Teboulle M.Fast gradient-based algorithms for constrained total variation image denoising and deblurring problems［J］.IEEE Transactions on Image Processing，2009，18(11):2419-2434.

［2］Louchet C，Moisan L.Total variation as a local filter［J］.SIAM Journal on Imaging Sciences，2011，4(2):651-694.

［3］Buades A，Coll B，Morel J M.A review of image denoising algorithms，with a new one ［J］.SIAM Journal on Multiscale Modeling and Simulation，2005，4(2):490-530.

［4］Dabov K，F(xiàn)oi A，Katkovnik V，et al.Image denoising by sparse3-D transform-domain collaborative filtering［J］.IEEE Transactions on Image Processing，2007，16(8):2080-2095.

［5］ZHANG L，DONG W，ZHANG D，et al.Two-stage image denoising by principal component analysis with local pixel grouping ［J］.Pattern Recognition，2010，43(4):1531-1549.

［6］Peyre G.Manifold models for signals and images［J］.Computer Vision and Image Understandin，2009，113(2):1077-3142.

［7］Niyogi P，Smale S，Weinberger S.A topological view of unsupervised learning from noisy data ［J］.SIAM Journal on Computing，2011，40(3):646-663.

［8］NI J，Turaga P，Patel V，et al.Example-driven manifold priors for image deconvolution ［J］.IEEE Transactions on Image Processing，2011，20(11):3086-3096.

［9］Wetzler A，Kimmel R.Efficient beltrami flow in patch-space［J］.Scale Space and Variational Methods in Computer Vision，2012(6667):134-143.

［10］Wright J，Ganesh A，Rao S，et al.Robust principal component analysis:Exact recovery of corrupted low-rank matrices via convex optimization［C］//Vancouver:Twenty-Fourth Annual Conference on Neural Information Processing Systems，2009:1-44.

［11］LIN Z，CHEN M，MA Y.The augmented lagrange multiplier method for exact recovery of corrupted low-rank matrices［R］.Urbana-Champaign:University of Illinois at Urbana-Champaign Technical Report，2010:1-23.

［12］LIN Z，CHEN M，MA Y.Inexact augmented lagrange multiplier(ALM)method［CP/OL］.［2012-05-04］.http://perception.csl.uiuc.edu/matrix-rank/sample_code.html.

［13］WANG Z，Bovik A C，Sheikh H R，et al.Image quality assessment:From error visibility to structural similarity［J］.IEEE Transcactions on Image Processing，2004，13(4):600-612.

［14］Dabov K.Block-matching and 3D filtering(BM3D)algorithm and its extensions［CP/OL］.［2012-05-04］.http://www.cs.tut.fi/～ foi/GCF-BM3D/.