王 巧，喬 雙

（東北師范大學(xué)物理學(xué)院，吉林 長(zhǎng)春 130024）

Richardson－Lucy（RL）算法是目前最為廣泛使用的圖像復(fù)原方法之一，是一種基于貝葉斯分析的迭代算法.RL算法假設(shè)圖像服從Poission分布，采用最大似然法進(jìn)行估計(jì)［1］.該算法在噪聲影響可忽略或較小的情況下具有唯一解.對(duì)于實(shí)際應(yīng)用中常見的低信噪比圖像，RL算法存在放大噪聲的缺陷，難以獲得較好的復(fù)原效果［2］.調(diào)制核算法是通過分析，估計(jì)梯度穩(wěn)定獲取圖像的局部結(jié)構(gòu)，并根據(jù)這個(gè)結(jié)構(gòu)信息去決定一個(gè)準(zhǔn)確的調(diào)制核的形狀和尺寸［3］.在平坦的區(qū)域，調(diào)制核的尺度是寬的，而且基本上是各向同性的，表明在這個(gè)區(qū)域內(nèi)部沒有強(qiáng)的方向性結(jié)構(gòu).在邊緣區(qū)域，調(diào)制核的形狀描述邊緣輪廓，并且核的值基本上表明像素的強(qiáng)度與感興趣的像素的相似性.在包含小規(guī)模圖像細(xì)節(jié)的區(qū)域，相應(yīng)的調(diào)制核也縮小到一個(gè)小區(qū)域［4］.得到了確定的調(diào)制核后，對(duì)退化圖像進(jìn)行處理，可以得到去噪圖像.由于Richardson－Lucy存在放大噪音的缺點(diǎn)，而調(diào)制核算法能夠?qū)υ胍暨M(jìn)行有效控制，因而本文提出將2種算法相結(jié)合，圖像復(fù)原的效果得到了明顯改善.

1 RL算法

RL算法是一種基于貝葉斯迭代運(yùn)算的圖像復(fù)原技術(shù)，迭代中加入了一些圖像的先驗(yàn)知識(shí)作為約束條件，這種算法和以往的基于傅里葉變換的圖像復(fù)原算法不同，它假定觀察圖像，點(diǎn)擴(kuò)散函數(shù)和原始圖像服從泊松概率統(tǒng)計(jì)模型，并最終收斂于泊松統(tǒng)計(jì)的最大似然解［5］.

1.1 RL算法的基本原理

RL算法是從最大似然中引出來的，它不是簡(jiǎn)單地對(duì)退化圖像作求逆運(yùn)算，而是通過對(duì)復(fù)原圖像的逐步迭代來逼近真實(shí)圖像［6］.

設(shè)真實(shí)圖像信號(hào)f及其觀測(cè)成像結(jié)果g為N×M 的實(shí)值矩陣，令H∶＝［hi，j］是r×s的模糊核，其中

通常人們認(rèn)定模糊算子是線性的，同時(shí)滿足平移不變性，即在觀測(cè)得到的降質(zhì)圖像中，其像素值都是由同一個(gè)模糊核產(chǎn)生的.于是圖像降質(zhì)過程可以解釋為

這里P表示針對(duì)H＊f的Poisson噪音污染， 表示卷積運(yùn)算.

在實(shí)際成像過程中，我們只能獲得降質(zhì)圖像g，于是我們將g看成是一個(gè)服從Poisson分布的隨機(jī)變量，因此有

將所有的像素視為彼此獨(dú)立分布的，并代入（3）式，從而得到似然函數(shù)

對(duì)（4）式兩端取對(duì)數(shù)得

于是我們可以定義E（f）＝－logP（g｜f），f的最大似然估計(jì)fML則成為下列極值問題的求解

這里

其中l(wèi)表示為N×M的矩陣，〈〉表示RN×M的內(nèi)積計(jì)算.

對(duì)E（f）中的f求導(dǎo)，并令其導(dǎo)數(shù)等于零，于是可推得

這里HT表示為矩陣H的共軛轉(zhuǎn)置，而上式中的表示A和B矩陣中對(duì)應(yīng)分量的逐點(diǎn)除法.由（1）式中

所以顯然可有HT＊l＝1.于是將f視為不動(dòng)點(diǎn)，可以給出如下RL算法

這里A·b表示A和B矩陣中對(duì)應(yīng)分量的逐點(diǎn)乘法.當(dāng)k→∞時(shí)，fk→f.

1.2 RL算法迭代次數(shù)的影響

用Matlab進(jìn)行仿真實(shí)驗(yàn)，所用圖像采用長(zhǎng)度為6，方向?yàn)?°的運(yùn)動(dòng)模糊，噪聲均值為0，方差為0.0001的高斯白噪聲的lena圖像，如圖1所示.

圖1 高斯白噪聲的lena圖像

使用RL方法在不同迭代次數(shù)下復(fù)原圖1（b）的結(jié)果見圖2.其中k為迭代次數(shù)，SNR為信噪比.

1.3 RL算法對(duì)不同模糊程度圖像的復(fù)原效果

圖3所示是不同模糊程度的圖像.使用RL算法對(duì)于不同模糊程度圖像迭代20次的結(jié)果見圖4.

從實(shí)驗(yàn)結(jié)果來看，RL算法能夠很好地為圖像去模糊、去噪，理論上迭代算法的迭代次數(shù)越多，圖像就會(huì)越精確，但由于實(shí)際的退化圖像中存在噪聲，迭代的次數(shù)過多會(huì)放大噪聲，使復(fù)原圖像的信噪比降低［7］.這個(gè)問題若得不到解決，對(duì)于實(shí)際應(yīng)用中常見的低信噪圖像，RL算法難以獲得較好的復(fù)原效果.

圖2 不同迭代次數(shù)復(fù)原圖1（b）的圖像

圖3 不同模糊程度的含噪圖像

圖4 對(duì)不同模糊程度的含噪圖像的復(fù)原

2 調(diào)制核算法

調(diào)制核算法能夠從模糊含噪圖像中獲得圖像的結(jié)構(gòu)，該算法的核心是根據(jù)圖像的某個(gè)像素的局部結(jié)構(gòu)確定該像素處的調(diào)制核的尺寸和大小，當(dāng)獲得確定的調(diào)制核后與含噪圖像卷積處理，從而能夠有效地抑制噪聲.

2.1 調(diào)制核算法的基本原理

假設(shè)一個(gè)像素位于位置xi＝［xi，yi］，它的Steering－Kernel（簡(jiǎn)記為SK）的數(shù)學(xué)表示為

其中：xi代表一個(gè)已經(jīng)給定的位置在以xi為中心的SK的窗口內(nèi)；Cl是從在圍繞xl的分析窗口內(nèi)的灰度的集合估計(jì)出的協(xié)方差矩陣.從定義上看，SK就像一個(gè)簡(jiǎn)單的高斯函數(shù)，但是最終的SK結(jié)果不只是簡(jiǎn)單的橢圓外形，這主要是因?yàn)閱为?dú)的協(xié)方差矩陣Cl是在每一個(gè)像素位置處被估計(jì)后使用的.

對(duì)于以xi為中心的窗口wl的局部灰度矩陣被定義為

Gx（xm）和Gy（xm）表示在［xm，ym］的灰度.在wl區(qū)域內(nèi)其主導(dǎo)方向V1和其垂直方向V2能夠通過計(jì)算D的奇異值分解被估計(jì)為

這里的奇異值V1≥s2≥0代表V1和V2各自方向上的能量.矩陣Cl可以通過（11）式估計(jì)出來.

當(dāng)調(diào)制核獲得后，該算法可表示為

2.2 用Matlab編程檢驗(yàn)調(diào)制核去噪效果

圖5（a）為原始的lena圖像，圖5（b）為模糊加噪后的圖像，圖5（c）為使用SK算法去噪后的圖像.從實(shí)驗(yàn)結(jié)果可以看出，調(diào)制核算法能夠有效地抑制圖像中的強(qiáng)噪音.

圖5 原始lena、模糊加噪和SK去噪效果對(duì)比

3 調(diào)制核算法與RL算法結(jié)合

RL算法在迭代次數(shù)比較低的情況下，隨著迭代次數(shù)的增加，提高圖像的信噪比逐漸提高，但當(dāng)?shù)螖?shù)增加到一定程度時(shí)，迭代次數(shù)的增加會(huì)放大噪音，信噪比反而下降.調(diào)制核算法能夠根據(jù)圖像的局部結(jié)構(gòu)很好地抑制了噪聲，能夠提高圖像的信噪比.因此，先對(duì)模糊含噪圖像運(yùn)用調(diào)制核算法去噪，提高圖像的信噪比，然后運(yùn)用Richardson－Lucy迭代算法復(fù)原圖像.圖6是幾種方法用于中子圖像去模糊、去噪的結(jié)果.

從實(shí)驗(yàn)結(jié)果可以看出，RL算法由于存在噪聲放大的缺陷，視覺效果不理想.將調(diào)制核與Richardson－Lucy算法相結(jié)合后，圖像復(fù)原的效果得到顯著改善.

圖6 幾種方法應(yīng)用于中子圖像的結(jié)果

4 結(jié)論

Richardson－Lucy算法通過反復(fù)迭代能夠很好地復(fù)原信噪比高的模糊圖像.調(diào)制核算法能根據(jù)圖像的局部結(jié)構(gòu)抑制圖像中存在的強(qiáng)噪聲，將調(diào)制核算法與圖像復(fù)原算法相結(jié)合，即先對(duì)模糊含噪圖像運(yùn)用調(diào)制核算法去噪，提高圖像的信噪比，再運(yùn)用Richardson－Lucy算法復(fù)原圖像能夠獲得很好的復(fù)原視覺效果.與單純使用某種復(fù)原算法相比較，圖像的信噪比明顯增強(qiáng)，視覺質(zhì)量明顯提高.

［1］陳云龍，王平，王鵬.基于L－R非線性迭代的降質(zhì)圖像復(fù)原算法［J］.計(jì)算機(jī)工程，2010，36（4）：202－205.

［2］TAKEDA H，F(xiàn)ARSIU S，MILANFAR P.Kernel regression forimage processing and reconstruction［J］.IEEE Transactions on Image Processing，2007，16（2）：349－366.

［3］ACHIM A，TSAKALIDES P，BEZERIANOS A.SAR Image denoising via bayesian wavelet shrinkage based on heavy－tailed modeling［J］.IEEE Trans on Geoscience and Remote Sensing，2003，1（8）：1773－1784.

［4］LUCY L J.An iterative technique for the rectification of observed distributions astron［J］.Astronomical Journal，1974，79（6）：745.

［5］閆河，閆衛(wèi)軍，李唯唯.基于Lucy－Richardson算法的圖像復(fù)原［J］.計(jì)算機(jī)工程，2010，36（15）：204－206.

［6］MINH N D，VETTERLI M.The contourlet transform：an efficient directional multisolution image representation［J］.IEEE Trans On Image Processing，2005，14（12）：2091－2106.

［7］MASSCHAELE B，DIERICK M，HOOREBEKE L V，et al.Neutron CT enhancement by iterative de－blurring of neutron transmission images［J］.Nucl Instrand Meth A，2005，542（1）：361－366.

［8］RICHARDSON W H.Bayesian－based iterative method of image restoration［J］.J Opt Soc Am，1972，62（1）：55－59.