嚴重噪聲污染下的模糊圖像復原方法

周 徽

(中國電子科技集團公司第二十八研究所，南京 210007)

數(shù)字圖像復原/恢復技術(shù)(以下簡稱復原技術(shù))是數(shù)字圖像處理的重要研究內(nèi)容，已在諸多領域取得了成功應用［1］。

由于現(xiàn)有的圖像系統(tǒng)都存在某種程度的缺陷，呈現(xiàn)出各種退化因素，因而進行圖像復原是必要的?？傮w上退化因素可以歸結(jié)為2類。①空間退化，又稱為模糊。成像系統(tǒng)中某些元件失常、圖像傳感器與成像目標之間的相對移動、攝像機失焦、不良天氣影響等因素都會造成圖像模糊。②點退化，通常指噪聲。最常見的點退化來自電子元件或數(shù)字化噪聲。

圖像退化過程可被模型化為一個退化函數(shù)和一個加性噪聲項共同作用于原始圖像f(x，y)，從而產(chǎn)生一幅退化的圖像g(x，y)，如圖1所示。根據(jù)該模型，退化圖像的數(shù)學描述見式(1)。

其中:g(x，y)代表退化圖像;f(x，y)代表未退化圖像;n(x，y)代表噪聲項;h(x，y)代表模糊算子，又稱點擴展函數(shù)(point spread function，PSF);符號“*”代表卷積。為了便于描述，式(1)改寫為矩陣形式

其中:f、g和n為NM×1維向量，由二維圖像按行或者列堆積而成;H是由PSF生成的MN×MN塊循環(huán)矩陣。

圖1 基本圖像退化/復原模型

理論上，在給定f(x，y)以及退化函數(shù)H和噪聲n(x，y)的一些先驗知識后，便可以獲得原始圖像的一個近似估計f^(x，y)。然而，式(2)具有病態(tài)特性，導致直接求解不容易實現(xiàn)，只能先設法對式(2)施加一定的約束，將其病態(tài)轉(zhuǎn)為良態(tài)后再設法求其最優(yōu)估計解。這個估計應能盡可能接近原始圖像。

針對該問題，人們提出過許多方法，如在已知g、H與n有關(guān)的先驗知識的前提下，去求原圖像的逼近解。在這些方法當中，Tikhonov正則化方法最具優(yōu)勢［2］。該方法是一個在L2空間的最小方差近似問題，能夠平衡精確性與平滑性［3］，因而得到廣泛研究。它將式(2)變?yōu)閹Ъs束條件的優(yōu)化解，即

其中:ψdata(f，g)表示數(shù)據(jù)逼近項;φreg(Qf)表示正則化項;Q為高通濾波算子;d為正則化參數(shù)，用以控制數(shù)據(jù)逼近項與正則化項的空間分布。理論上，在 L(α，f)最小化，即條件下的解就是原圖像f(x，y)的最優(yōu)逼近解。

在正則化方法應用中確定正則化參數(shù)是非常關(guān)鍵的一步。人們先后提出過許多關(guān)于如何選擇正則化參數(shù)的方法。Galatsanos等［4］提出了基于均方誤差(MSE)評判標準的正則化參數(shù)選擇和噪聲變量估計方法。Nakano等［5］提出用權(quán)矩陣的方式優(yōu)化正則化效果，而局部正則化優(yōu)化權(quán)矩陣從迭代過程獲得。Reeves［6］則給出了在限制性條件下，用交叉有效性準則判斷正則化最優(yōu)度的方法。Molina等［7］將等級貝葉斯方法應用到圖像復原中，提出了一個3步方法。該方法在每個步驟中完成1個正則參數(shù)的確定，降低了對噪聲先驗知識的依賴性。Leung［8］提出用L曲線作為選取正則化參數(shù)的依據(jù)，通過實驗驗證了L曲線的曲率與正則化參數(shù)優(yōu)異性之間的相關(guān)性，其結(jié)論是:當曲率最大化時，取得的參數(shù)是最優(yōu)的。Demoment［9］綜述了正則化方法的詳細內(nèi)容。

上面提到的正則化參數(shù)選取方法都具有全局性限制，即參數(shù)在算法開始之前就已經(jīng)確定，導致過分依賴于先驗信息，無法在算法迭代過程中自適應修正。與上面的方法不同，Katsaggelos等［10］提出用正則化函數(shù)取代全局正則化參數(shù)的方法，使算法能夠自適應地修改正則化參數(shù)，大幅提高了算法的效率。本文的工作就是充分利用該思想，提出一種改進的自適應加權(quán)正則化迭代圖像復原方法，增強算法對強噪聲存在時的適應性。

1 自適應正則化函數(shù)

用含有自變量f的正則化函數(shù)α(f)代替全局正則化參數(shù)α，建立目標平滑泛函式(4)。

目標泛函建立之后，正則化參數(shù)的選擇至關(guān)重要。自適應正則化函數(shù)α(f)的選取應滿足3條性質(zhì)［10］:① α(f)為平滑泛函的函數(shù);② α(f)與目標泛函之間需滿足一定的極值性質(zhì)，即當時，α(f)→0;‖g － Hf‖時，α(f)→∞;③ 在α(f)取值范圍內(nèi)不改變目標泛函的凸函數(shù)性質(zhì)。

根據(jù)上面性質(zhì)，假設存在正實數(shù)λ使得式(5)成立，則

正則化項用于懲罰解的粗糙性，然而過分的懲罰會造成解的過分平滑，抑制算法對高頻成分的恢復，從而導致復原結(jié)果在圖像強度變化劇烈處產(chǎn)生振鈴效應，影響圖像視覺效果。因此，在復原濾波時，必須將圖像局部細節(jié)特性考慮在內(nèi)，采取加權(quán)形式增強算法的局部適應能力，提高算法對振鈴的抑制能力，改善復原質(zhì)量［5］。

文獻［11］認為，圖像中小的梯度值對應噪聲，需對其進行平滑，而大的梯度值則對應圖像邊緣輪廓，應給予保護。根據(jù)該理論，定義加權(quán)正則化項為

其中W是由wij構(gòu)成的矩陣，wij定義為

加權(quán)后的自適應正則化函數(shù)如式(10)所示。

2 迭代復原算法及收斂性分析

圖像復原過程等價于式(4)定義的平滑函數(shù)逼近最優(yōu)解的過程。由于目標泛函的凸性，恰當?shù)剡x取正則化參數(shù)能夠使式(4)存在唯一的全局最優(yōu)解。該最優(yōu)解在平滑函數(shù)梯度為0時取得，即

將式(11)對f求導，得

因為α(fk)是第k次迭代結(jié)束后計算出的參數(shù)，用于第k+1次迭代，所以▽fα(f)=0，從而式(12)化簡為

建立逐步逼近的Van Cittert格式的迭代等式:

首先分析算法的收斂性，改寫式(14)為

對式(15)左右兩邊同時取范數(shù)，并根據(jù)矩陣范數(shù)與向量范數(shù)的協(xié)調(diào)性，有

反復使用這個不等式得

下面分析算法的收斂條件。因為‖A‖＜1等價于 ρ(A)＜ 1，有

其中ρ(·)定義為“·”的譜半徑。根據(jù)三角不等式，上式等價于

通過規(guī)范化H、C可使ρ(HTH)和ρ(CTC)均不大于1，從而算法收斂的充分條件演變?yōu)棣?f)＜1，即

算法流程描述:

步驟1 設置初始值f0=g，迭代計數(shù)值k=0，設定算法終止條件ε=10－7，開定時器。

步驟2 在第k步迭代中，自動修正正則化參數(shù)αk(fk)，并利用式(14)對圖像進行復原濾波。更新fk+1。

步驟4 關(guān)定時器，計算復原評價指標ISNR，顯示復原結(jié)果。

3 實驗與分析

選擇2組實驗驗證本文提出算法的有效性。第1組是樣條圖像，用來驗證算法對圖像邊緣恢復的效果;第2組是經(jīng)典的lena圖像，用來驗證算法對圖像細節(jié)恢復的效果。

實驗1 樣條圖像尺寸為158×300的灰度圖像，如圖2(a)所示。本實驗在PC機(CPU為P42.4GHz，內(nèi)存512 Mbytes)上進行，模糊算子H取為11×11的均勻模糊，隨機疊加信噪比為20 dB的高斯噪聲，C取2－D拉普拉斯算子。判據(jù)‖fk+1－fk‖/‖fk‖ ＜10－7作為算法終止迭代的條件，改進信噪比(improvement signal-to-noise ratio，ISNR)作為算法性能的評價指標，ISNR定義為

圖2 樣條圖像復原效果

圖2(a)、(c)、(e)分別為原始樣條圖像、退化圖像和復原圖像。圖2(b)、(d)、(f)分別為穿越圖2(a)、(c)、(e)中部的水平掃描線，其橫坐標為水平象素點，縱坐標為對應象素的亮度值。

由圖2(e)可以看出，相對于圖2(c)，算法恢復的視覺效果是比較理想的。圖2(f)給出的復原圖像中部的掃描線說明樣條交接的輪廓線基本得以復原，顯然較圖2(d)有明顯改善。實驗獲得的具體參數(shù)見表1和圖4。

實驗2 lena圖像尺寸為300×300的灰度圖像，如圖3所示。其他條件同實驗1。圖3(a)、(b)分別為退化的lena圖像和復原后的lena圖像。圖3(c)、(d)分別為圖3(a)、(b)的局部放大圖像。實驗獲得的具體參數(shù)見表1和圖4。

圖3 lena圖像恢復效果

表1 實驗結(jié)果參數(shù)

圖4 收斂誤差曲線

圖4給出2組實驗中每一迭代步驟對應的誤差判據(jù)變化曲線，其中圖4(a)是樣條圖像，圖4(b)是lena圖像，橫坐標對應迭代次數(shù)，縱坐標對應誤差判據(jù)值。由圖4可以看出，無論迭代初始條件如何，算法能夠以較少的迭代次數(shù)快速平穩(wěn)地趨向于最優(yōu)估計解。

4 結(jié)束語

本文提出一種改進的自適應加權(quán)正則化迭代圖像復原算法。該算法能夠自適應地選擇并自動修正正則化參數(shù)。算法對初始值不敏感，在每步迭代中，算法完成1次復原濾波，并更新正則化參數(shù)，確保結(jié)果能夠快速趨向于最優(yōu)。討論了算法的收斂性與收斂條件，通過建立顯式迭代復原算法尋找全局最優(yōu)解。該算法具有占用內(nèi)存小、計算速度快等優(yōu)點。無需任何噪聲信息作為先驗知識，算法通過觀測圖像本身自適應估計并調(diào)整正則化參數(shù)以適應強噪聲存在下的復原操作。用2組實驗驗證了算法的有效性。結(jié)果顯示:該算法能提高算法對強噪聲存在的適應性，抑制Gibbs振鈴波紋，具有較好的評價結(jié)果和較高的實時性。

［1］Banham M R，Katsaggelos A K.Digital image restoration［J］.Signal Processing，1997，14(2):24-41.

［2］Bouhamid A，Jbilou K.Sylvester Tikhonov-regularization methods in image restoration［J］.Journal of Computational and Applied Mathematics，2007，206(1):86-98.

［3］Leung C M，Lu W S.Optimal determination of regularization parameters and the stabilizing operator［J］.Signal Processing，1995，201:403-406.

［4］Galatsanos N P，Katsaggelos A K.Methods for choosing the regularization parameter and estimating the noise variance in image restoration and their relation［J］.IEEE Transactions on Image Processing，1992，1(3):322-336.

［5］Nakano K，Eguchi M，Toyota Y，et al.On regularization for image restoration problems from the viewpoint of a Bayesian information criterion ［C］//IECON-93.Hawaiian，USA:［s.n.］，1993:2257-2261.

［6］Reeves S J.Optimal regularized image restoration with constraints［C］//ICASSP-92.San Francisco，USA:［s.n.］，1992:301-304.

［7］Molina R，Katsaggelos A K，Mateos J.Bayesian and regularization methods for hyperparameter estimation in image restoration.Image Processing［J］.IEEE Transactions on Image Processing，1999，8(2):231-246.

［8］Leung C M，Lu W S.An L-curve approach to optimal determination of regularization parameter in image restoration.Electrical and Computer Engineering［J］.Electrical and Computer Engineering，1993，2:1021-1024.

［9］Demoment G.Image reconstruction and restoration:overview of common estimation structures and problems［J］.IEEE Transactions on Signal Processing，1989，37(12):2024-2036.

［10］Moon G K，Katsaggelos A K.General choice of the regularization functional in regularized image restoration［J］.IEEE Transactions on Image Processing，1995，4(5):594-602.

［11］Charbonnier N P，Blanc F L.Deterministic edge-preserving regularization in computed imaging［J］.IEEE Transactions on Image Processing，1997，6(2):298-311.