洪景新，陳 栩

（廈門大學(xué)信息科學(xué)與技術(shù)學(xué)院，福建 廈門361005）

由于受到成像系統(tǒng)本身的缺陷、成像過程中傳播介質(zhì)的雜質(zhì)，以及相對運(yùn)動(dòng)的變化等因素的影響，圖像信息的采集過程難免會對圖像信息造成某些失真和不同程度的降質(zhì)，這一現(xiàn)象稱為圖像退化.利用數(shù)字圖像處理技術(shù)，對退化圖像加以改善并復(fù)原原圖像信息，這個(gè)過程就稱為圖像復(fù)原.如今眾多的應(yīng)用領(lǐng)域?qū)D像質(zhì)量的要求都相當(dāng)高，尤其對精度和清晰度的要求，因此圖像的復(fù)原問題愈加具有意義.

復(fù)原技術(shù)是基于將退化過程模型化并采用逆過程處理模糊圖像，實(shí)現(xiàn)圖像的復(fù)原.逆濾波算法［1］是最早用于圖象復(fù)原的一種標(biāo)準(zhǔn)技術(shù)，該算法容易受到噪聲的干擾，噪聲較大時(shí)，圖像恢復(fù)的效果很差.隨后，Helstrom［2］采用最小均方誤差估計(jì)方法，提出了基于維納濾波器的濾波算法，該方法只能針對某個(gè)具體圖像設(shè)計(jì)最優(yōu)濾波器，對于其他圖像則不一定能達(dá)到較好效果.Hunt和Andrews［3］對逆濾波、維納濾波進(jìn)行了對比研究，并提出了約束最小二乘濾波方法，這是一種基于線性表達(dá)的恢復(fù)方法，對各種退化圖像的恢復(fù)具有較好的適應(yīng)性.近幾年，也出現(xiàn)了將稀疏表示［4］運(yùn)用到圖像復(fù)原中的研究，使模糊圖像的恢復(fù)效果不斷改善.

本文首先闡述圖像復(fù)原的一般退化模型與勻速運(yùn)動(dòng)的特性，從中推導(dǎo)出勻速運(yùn)動(dòng)模糊的退化矩陣.然后，從求解大型線性方程組的角度，分析常見的約束最小二乘濾波算法存在的問題，說明l1正則化在圖像復(fù)原上的優(yōu)勢，提出一種基于l1正則化的復(fù)原算法，并利用基追蹤算法對正則化問題進(jìn)行求解.最后對不同模糊尺度的圖像進(jìn)行仿真實(shí)驗(yàn)，大量的實(shí)測數(shù)據(jù)和圖像表明，本文算法對較寬范圍的模糊尺度都有著很好的復(fù)原效果，同時(shí)能有效抑制振鈴效應(yīng).

1 模糊圖像退化模型

1.1 一般退化模型

一般退化過程可以建模為一個(gè)退化函數(shù)和一個(gè)加性噪聲項(xiàng)的級聯(lián)［1］，其退化模型可以表示為：

其中，y，x∈Rm分別為模糊圖像和原圖像，H∈Rm×m是退化過程的表示矩陣，n∈Rm為退化過程中的加性噪聲.不考慮加性噪聲時(shí)，退化模型可簡化為：

此時(shí)，模糊圖像y∈Rm可看作是原圖像x∈Rm在H下的線性投影，而圖像去模糊的過程就是在退化函數(shù)H∈Rm×m已知或可測量的前提下，從模糊圖像y∈Rm中求解未知的原圖像x∈Rm.因此，去模糊的過程等價(jià)于求解線性逆問題，也可理解為從模糊圖像重構(gòu)出原圖像的過程.

1.2 勻速運(yùn)動(dòng)模糊圖像退化模型

運(yùn)動(dòng)模糊，是圖像退化的一個(gè)主要形式，是由于在圖像信息采集時(shí)，目標(biāo)物體與成像裝置的相對位移而產(chǎn)生的一種退化問題，例如在運(yùn)動(dòng)物體上（飛機(jī)、宇宙飛行器等）拍攝外界的照片或拍攝高速運(yùn)動(dòng)物體的照片均存在運(yùn)動(dòng)模糊問題.

勻速運(yùn)動(dòng)模糊是最基本的運(yùn)動(dòng)退化，各類復(fù)雜的運(yùn)動(dòng)退化模型可分解為多個(gè)勻速運(yùn)動(dòng)退化模型.假設(shè)成像時(shí)間T內(nèi)，目標(biāo)物體與成像裝置以勻速v做相對運(yùn)動(dòng)，選取圖像中一行像素點(diǎn)g（i），i＝0，1，…，n，將成像時(shí)間T劃為L個(gè)等長的單位時(shí)間t，且在0時(shí)刻，令物體上點(diǎn)f（0）投影在像素點(diǎn)g（0）上，若目標(biāo)物體每單位時(shí)間移動(dòng)一個(gè)像素，則kt時(shí)，物體上點(diǎn)f（0）的投影位置為像素點(diǎn)g（kt）.由此不難看出，每個(gè)像素點(diǎn)接受的信息為成像時(shí)間內(nèi)物體移動(dòng)范圍的各點(diǎn)的加權(quán)平均，即：

其中，L即模糊尺度.令模糊路徑上的點(diǎn)擴(kuò)散函數(shù)為：

則式（3）可寫成離散卷積的形式：

若用矩陣形式表達(dá)，即可得到退化矩陣：

2 l1－min圖像復(fù)原算法

2.1 約束最小二乘算法

求解線性方程的過程就是求解最小二乘問題.常見的最小二乘算法［5］是求解線性逆問題的經(jīng)典算法，其中心思想是通過最小化估算數(shù)據(jù)與實(shí)際數(shù)據(jù)之間的誤差的平方來尋找問題的最優(yōu)解，其函數(shù)表達(dá)式為：

在圖像復(fù)原中，矩陣H∈Rm×m一般是病態(tài)的，最小二乘的解對誤差（舍入誤差或其他誤差）很敏感，導(dǎo)致求得的解往往不是理想的，進(jìn)而使得復(fù)原效果很差，最小二乘也因此變得沒有意義.為了解決這個(gè)問題，可通過正則化來約束最終的解.

正則化的基本思想是用一族與不適定問題有著近似解的適定問題來逼近原問題.選擇合適的正則化方法，可確保我們得到理想的解.約束最小二乘算法中通常采用Tikhonov正則化，對式（7）進(jìn)行Tikhonov正則化可得：

式（8）即為常見的約束最小二乘算法.式中第一項(xiàng)為殘差約束項(xiàng)，用以控制解的逼近程度，第二項(xiàng)為正則約束項(xiàng)，用以約束信號的能量.正則參數(shù)λ＞0用于平衡精確度與噪聲敏感度.L是微分算子，一般選用拉普拉斯算子.約束最小二乘算法屬于基于l1范數(shù)的正則化，實(shí)踐表明，當(dāng)方程y＝Hx中的矩陣H與向量y有擾動(dòng)時(shí)，約束最小二乘法的解往往不穩(wěn)健.

2.2 l1 正則化

l1正則化是采用l1范數(shù)作為了約束條件，來求解最小二乘問題的正則化方法.因此用l1范數(shù)作為能量約束取代式（8）中的正則約束項(xiàng)，則可以得到：

其中，l1正則約束項(xiàng)是稀疏約束條件.l1范數(shù)指是向量中各元素的絕對值之和，若個(gè)別元素出現(xiàn)擾動(dòng)，對整體影響相比l2范數(shù)較小.因此，不同于基于l2的正則化形式，基于l1的正則化對異常值較不敏感［6］.在圖像處理中，色塊之間往往有尖銳的邊緣，會出現(xiàn)突變的異常值，因此基于l1的正則化更適用于圖像處理.但同時(shí)，由于l1范數(shù)的稀疏約束條件，使得求解出來的x∈Rm變稀疏.而對于圖像這類自然信號來說，幾乎所有的像素灰度值都是非零的，所以要對圖像進(jìn)行稀疏變換.

圖像的稀疏表示是基于 Mallat等［7］提出的信號在過完備庫上的分解.根據(jù)稀疏表示理論，若存在一個(gè)過完備字典W，使圖像或信號可以表示為：

若向量α∈Rm中只有很少的大系數(shù)，則稱圖像或信號x∈Rm是可稀疏的.若α∈Rm只有K個(gè)非零元素，則稱α∈Rm為圖像或信號的K稀疏表示.因此，式（9）可轉(zhuǎn)化為：

常見的稀疏變換有傅里葉變換、離散余弦變換（DCT）、小波變換等.本文采用DCT變換對圖像進(jìn)行稀疏表示.

3 基追蹤算法

基于l1正則化的線性規(guī)劃問題的常見求解方法有基追蹤（BP）算法和迭代閾值算法（IST）等，本文選用BP算法進(jìn)行求解.

BP方法［8］是Chen等提出的一種信號稀疏表示法.它尋求從過完備的函數(shù)（基）集合中得到信號的最稀疏的表示，即用盡可能少的基盡可能精確地表示原信號，從而獲得信號的內(nèi)在本質(zhì)特性.

BP算法的思想是通過極小化l1范數(shù)來獲得最稀疏的解.即求解下式：

由文獻(xiàn)［9］可知，極小l1范數(shù)最優(yōu)問題與線性規(guī)劃問題等價(jià).可通過等價(jià)變換m?2n；x?（u，v）；c?（1，1）；A?（Φ，－Φ）；b?s，使其重構(gòu)成如式（13）所示的標(biāo)準(zhǔn)線性規(guī)劃問題：

其中，x∈Rm是一個(gè)重新定義的向量.因此，BP問題的解可以通過求解線性規(guī)劃問題來獲得.線性規(guī)劃問題的通常解法有單純形法、內(nèi)點(diǎn)法.本文采用基于內(nèi)點(diǎn)法的原－對偶障礙算法（primal－dual log－barrier algorithm）求解線性規(guī)劃問題：

其中，γ和δ為正則化系數(shù)，在不同的實(shí)際問題，對其進(jìn)行相應(yīng)的初始化.對于正則化形式如式（14）的線性規(guī)劃問題，原－對偶障礙算法的求解流程圖如下：

圖1 原－對偶障礙算法流程Fig.1 The flow chart of primal－dual log－barrier algorithm

原－對偶障礙算法的具體算法步驟如下：

Step 1：設(shè)置參數(shù)：可行性容差FeaTol、對偶間隙容差PDGapTol和正則化系數(shù)γ，δ.

Step 2：設(shè)置邊界：x＞0，y，z＞0，μ＞0，初始化k＝0.

Step 3：循環(huán)：

（i）t＝c＋γ2x－z－ATy，r＝b－Ax－δ2y，v＝μe－Zx，D＝（X－1Z＋γ2I）－1，其中，X，Z是分別由x，z構(gòu)成的對角矩陣，e是單位向量.

（ii）從（ADAT＋δ2I）Δy＝r－AD（X－1v－t）中求解Δy，并計(jì)算

Δx＝DATΔy＋D（X－1v－t），

Δz＝X－1v－X－1ZΔx.

（iii）計(jì)算原步長和對偶步長ρp，ρd并更新變量

ρp＝0.99×max｛ρ：x＋ρΔx≥0｝，

ρd＝0.99×max｛ρ：z＋ρΔz≥0｝，

x＝x＋ρpΔx，y＝y(tǒng)＋ρdΔy，z＝z＋ρdΔz，

μ＝（1－min（ρp，ρd，0.99））μ.

（iv）k＋＋.

Step 4：終止條件：

4 實(shí)驗(yàn)設(shè)計(jì)與結(jié)果分析

為了驗(yàn)證l1－min圖像復(fù)原算法對勻速運(yùn)動(dòng)模糊圖像的去模糊效果，本文從視覺效果和客觀指標(biāo)評價(jià)兩方面設(shè)計(jì)實(shí)驗(yàn)，并選用經(jīng)典維納濾波算法和約束最小二乘算法進(jìn)行比較.實(shí)驗(yàn)中，圖像選擇分辨率為256×256的標(biāo)準(zhǔn)灰度圖片，根據(jù)勻速運(yùn)動(dòng)圖像退化模型，生成退化矩陣分別對圖像進(jìn)行退化處理.模糊尺度L＝50的退化效果如圖2所示：

圖2 實(shí)驗(yàn)圖像和水平運(yùn)動(dòng)模糊退化圖像Fig.2 Original image and horizontal motion－blurred image

4.1 視覺效果

圖3為標(biāo)準(zhǔn)Lena和Peppers圖像在模糊尺度L＝50下3種算法的復(fù)原圖像.圖4為標(biāo)準(zhǔn)Cameraman圖像在不同模糊尺度下本文算法的復(fù)原圖像.從圖3實(shí)驗(yàn)結(jié)果可以清楚看到，3種復(fù)原算法都使圖片質(zhì)量相比退化圖像得到一定的改善.顯然，維納濾波圖像效果最差，圖像模糊嚴(yán)重，可以看到很多重影.最小二乘算法與l1－min算法的恢復(fù)圖像較為清晰.不過，從放大圖像可看出，最小二乘算法的恢復(fù)結(jié)果帶有振鈴現(xiàn)象，而本文算法恢復(fù)的圖像在完整地保存了邊界信息的同時(shí)幾乎看不出振鈴現(xiàn)象.從圖4實(shí)驗(yàn)結(jié)果可以看出，本文算法在不同模糊尺度下，都有良好的復(fù)原效果且隨模糊尺度的增長，復(fù)原效果只有少量下降.

圖3 不同算法的復(fù)原效果對比，模糊尺度L＝50Fig.3 Restoration comparison of different algorithms with blur length of 50

圖4 不同模糊尺度下本文算法復(fù)原效果對比Fig.4 Restoration comparison of the proposed algorithm with different blur length

4.2 客觀指標(biāo)評價(jià)

客觀指標(biāo)評價(jià)選取灰度平均梯度值（gray mean grads，GMG）［10］和結(jié)構(gòu)相似指數(shù)（structural similarity index measurement，SSIM）［11］作為衡量復(fù)原圖像的標(biāo)準(zhǔn).GMG屬于無參評價(jià)，它能較好反映復(fù)原圖像的對比度和紋理變化特征，其值越大表示圖像越清晰，圖像質(zhì)量越好.SSIM是有參評價(jià)，它可以衡量復(fù)原圖像與原圖像的相似度，其值越高表示2幅圖像越相似.

圖5為標(biāo)準(zhǔn)的Peppers圖片，在不同退化程度下，3種算法的復(fù)原圖像的評價(jià)指標(biāo)對比.其中，本文算法與最小二乘算法默認(rèn)迭代20次.圖中的實(shí)測數(shù)據(jù)可以看出，3種算法下復(fù)原圖像相比退化圖像的客觀指標(biāo)都有所提高，l1－min圖像復(fù)原算法在各退化程度下的兩種指標(biāo)都優(yōu)于其他算法.在不同模糊尺度下結(jié)構(gòu)相似指數(shù)SSIM表明l1－min圖像復(fù)原算法的復(fù)原圖像在結(jié)構(gòu)上比其他算法更相似性于原圖.而在GMG上，本文算法復(fù)原圖像的清晰度的優(yōu)勢更大.其中維納濾波算法隨模糊尺度的增長GMG值逐漸降低，而本文算法與最小二乘算法隨模糊尺度的增長相對比較穩(wěn)健，但l1－min算法明顯維持在更好的水平上.

圖5 實(shí)驗(yàn)中不同算法的評價(jià)指標(biāo)Fig.5 Evaluation criteria of different algorithms in experiment

5 結(jié) 論

本文在最小二乘算法的基礎(chǔ)上，用l1范數(shù)代替l2范數(shù)，再引入BP算法精確求解l1正則化的線性規(guī)劃問題，提出一種基于l1－min的勻速運(yùn)動(dòng)模糊圖像復(fù)原算法.算法對較寬的模糊范圍下的退化圖像都可以獲得一個(gè)較清晰的圖像復(fù)原效果，且能有效抑制約束最小二乘算法中的振鈴現(xiàn)象.但在個(gè)別特殊模糊尺度下，算法恢復(fù)出來的圖像的質(zhì)量突降，如圖4中模糊尺度80下的復(fù)原圖像中出現(xiàn)的波紋現(xiàn)象.因此認(rèn)為下一步工作應(yīng)探討和研究這些特殊模糊尺度導(dǎo)致圖像質(zhì)量下降的原因，爭取改善圖像質(zhì)量.

［1］Gonzalez R C，Woods R E.岡薩雷斯數(shù)字圖像處理［M］.北京：電子工業(yè)出版社，2003.

［2］Helstrom C W.Image restoration by the method of least squares［J］.JOSA，1967，57（3）：297－303.

［3］Andrews H C，Hunt B R.Digital image restoration［M］.New York：Englewood Cliffs，1977.

［4］Zhang H C，Yang J C，Zhang Y N，et al.Sparse representation based blind image deblurring［C］∥IEEE International Conference on Multimedia and Expo（ICME）.Barcelona，Spain：IEEE，2011：1－6.

［5］Bjorck A.Numerical methods for least squares problems［M］.Philadelphia：Society for Industrial and Applied Mathematics（SIAM），1996.

［6］Beck A，Teboulle M.A fast iterative shrinkage－thresholding algorithm for linear inverse problems［J］.SIAM Journal on Imaging Sciences，2009，2（1）：183－202.

［7］Mallat S G，Zhang Z.Matching pursuits with time－frequency dictionaries［J］.IEEE Transactions on Signal Processing，1993，41（12）：3397－3415.

［8］Chen S S，Donoho D L，Saunders M A.Atomic decomposition by basis pursuit［J］.SIAM Journal on Scientific Computing，1998，20（1）：33－61.

［9］Bloomfield P，Steiger W L.Least absolute deviations：theory，applications，and algorithms［M］.Boston：Birkh?user，1983.

［10］袁萬立.模糊圖像復(fù)原及評價(jià)方法的研究［D］.無錫：江南大學(xué)，2012.

［11］Wang Z，Bovic A C，Sheikh H R，et al.Image quality assessment：from error visibility to structural similarity［J］.IEEE Transactions on Image Processing，2004，13（4）：600－612.