王⒗飛 王 春

（1. 銅陵職業(yè)技術(shù)學(xué)院，安徽 銅陵 244061；2.皖南醫(yī)學(xué)院，安徽 蕪湖 241000）

醫(yī)學(xué)圖像在病理研究中起到極其重要的作⒚，運(yùn)動(dòng)模糊醫(yī)學(xué)圖像修復(fù)技術(shù)是目前圖像處理熱點(diǎn)之一[1]。 運(yùn)動(dòng)模糊醫(yī)學(xué)圖像是由于醫(yī)學(xué)攝像設(shè)備采集病理圖像時(shí)相對運(yùn)動(dòng)而產(chǎn)生像移，目前，文獻(xiàn)[2]提出了消除像移的方法主要有機(jī)械式、光學(xué)式和圖像式等像移補(bǔ)償方式；文獻(xiàn)[3]提出了像移補(bǔ)償算法；文獻(xiàn)[4]和文獻(xiàn)[5]提出傳統(tǒng)運(yùn)動(dòng)模糊圖像去模糊算法如Lucy－Richardson 復(fù)原、最小二乘方濾波、維納濾波、盲反卷積等，但都需要先驗(yàn)知識(shí)和相應(yīng)的約束條件。 文獻(xiàn)[6]提出采⒚Canny 算子抑制亮十字圖像去模糊法。

本文算法是在無先驗(yàn)知識(shí)的情況下， 估算出圖像復(fù)原的兩個(gè)重要參數(shù)：模糊角度和模糊長度。 具體做法：首先，通過二維離散傅里葉變換將模糊醫(yī)學(xué)圖像頻譜化，引入漢寧窗減小頻譜泄露，引入Radon 變換， 將圖像空間直線提取問題轉(zhuǎn)換為函數(shù)極值搜索問題，估算運(yùn)動(dòng)模糊角度；其次，對圖像水平方向一階微分，逐行進(jìn)行自相關(guān)運(yùn)算，找出所有行中出現(xiàn)最多的兩個(gè)最小值， 此兩值之間距離的一半就是運(yùn)動(dòng)模糊長度；重構(gòu)圖像復(fù)原的點(diǎn)擴(kuò)散函數(shù)，完成運(yùn)動(dòng)模糊醫(yī)學(xué)圖像去模糊。

一、運(yùn)動(dòng)醫(yī)學(xué)圖像模糊的基本原理

本文去模糊算法的關(guān)鍵是分析圖像退化機(jī)理，建立圖像退化模型。圖像退化實(shí)質(zhì)是原圖像f（x,y）通過外在H 系統(tǒng)加入一個(gè)外來加性噪聲n（x,y）而退化成一幅圖像g（x,y）。 H 系統(tǒng)有線性和非線性、或空間變化和空間不變化等特征， 圖像修復(fù)可視為從退化圖像g（x,y）中，清除H 系統(tǒng)輸入信息㈦輸出信息之間某種內(nèi)在聯(lián)系的響應(yīng)函數(shù)h（x,y）和外來加性噪聲函數(shù)n（x,y），估測原圖像f（x,y）的過程，即

若不考慮加性噪聲n（x,y），即n（x,y）=0，單路信號圖像退化模型可表示為

H 在光學(xué)成像系統(tǒng)中是一個(gè)點(diǎn)光源的沖擊響應(yīng)，也被稱為點(diǎn)擴(kuò)散函數(shù)（PSF）。

兩路信號f1（x,y）和f2（x,y）分別對應(yīng)輸出信號g1（x,y）和g2（x,y），圖像退化模型可為：

其中k1和k2為常數(shù)，則該系統(tǒng)為線性系統(tǒng)；

若α 和β 分別是空間位置的位移量，有

則H 是空間不變系統(tǒng)； 同時(shí)滿足上述兩條件的H 系統(tǒng)就是線性空間不變系統(tǒng)。

因線性系統(tǒng)有均勻性和相加性等特征， 實(shí)際成像的非線性和空間變化系統(tǒng)可在一定條件下近似地視為線性空間不變系統(tǒng)。

光學(xué)中單位沖擊為一光點(diǎn)， 在原點(diǎn)處振幅為無窮大，在原點(diǎn)外所有時(shí)刻振幅為零，寬度極小，設(shè)為函數(shù)δ（t），其卷積取樣公式為：

由于g（x,y）=H[f（x,y）]+n（x,y），有

令h（x,α,y,β）=H[δ（x-α，y-β）]，

則有

退化過程點(diǎn)擴(kuò)散函數(shù)（PSF）h（x,α,y,β）是系統(tǒng)H對坐標(biāo)（α，β）處的沖擊函數(shù)δ（x,α,y,β）的響應(yīng)。 將f（α，β）和h（x-α，y-β）進(jìn)行均勻取樣且fe（x）和he（x）在x 和y 方向上的周期分別為M 和N，有g(shù)e（x）在x和y 方向上的周期分別為M 和N，有：

其中G（u,v）、F（u,v）和N（u,v）分別是g（x,y）、f（x,y）和n（x,y）的二維傅里葉變換。

求解H（u,v）是圖像復(fù)原的重點(diǎn)，即求解醫(yī)學(xué)圖像運(yùn)動(dòng)模糊長度和方向是本文算法的關(guān)鍵。

二、運(yùn)動(dòng)模糊醫(yī)學(xué)圖像去模湖算法分析

對任意角度運(yùn)動(dòng)的醫(yī)學(xué)模糊圖像都可以認(rèn)為是線性的，可通過以圖像中心為原點(diǎn)O 旋轉(zhuǎn)獲得的。調(diào)整旋轉(zhuǎn)角度㈦X 軸重合， 可簡化其它角度運(yùn)動(dòng)模糊為水平方向模糊， 運(yùn)動(dòng)模糊長度為s，sa和sb分別是表示坐標(biāo)軸上的投影運(yùn)動(dòng)模糊長度，t 表示X 軸的單位向量；n 為Y 軸的單位向量。 估算出運(yùn)動(dòng)模糊角度和運(yùn)動(dòng)模糊尺度后， 就可以通過點(diǎn)擴(kuò)散函數(shù)對模糊圖像進(jìn)行復(fù)原。

運(yùn)動(dòng)模糊圖像點(diǎn)擴(kuò)散函數(shù)為：

（一）估算運(yùn)動(dòng)模糊角度并改進(jìn)

本文算法將圖像通過二維傅里葉變換成頻譜圖，低頻部分顯示為亮點(diǎn)，相同的低頻部分組合可顯示為亮線，高頻部分顯示為暗點(diǎn)，分析各種頻率下的幅值，將頻譜圖中的顏色條直線化并頻譜中心化，通過檢測直線投影角度確定圖像的運(yùn)動(dòng)模糊角度。

二維離散傅里葉變換公式為：

其中頻率變量u=0,1,2…M-1，v=0,1,2…N-1。

對模糊角度為30 度和模糊距離為30 像素的血管彩色圖像作為測試圖像，如圖1 所示。

圖1 運(yùn)動(dòng)模糊圖像頻譜對比

從圖1 可知，頻譜圖中間最亮點(diǎn)是最低頻率點(diǎn)，亮度表示信號強(qiáng)弱。 原始圖頻譜中只有一個(gè)亮點(diǎn)，而運(yùn)動(dòng)模糊頻譜圖中有多條明顯的平行亮線， 但是亮線寬度較大，影響檢測亮線投影角度，且方向㈦運(yùn)動(dòng)模糊方向垂直。

對圖像應(yīng)⒚二維傅里葉變換， 將最低頻部分移到頻譜圖中心，引入頻Ⅱ?yàn)V波，將模糊圖頻譜中噪聲降低，提高亮線的直線特征。 如圖2 所示。

圖2 頻譜亮線細(xì)化

圖2 中， 因圖像從中心向邊緣過渡時(shí)灰度值產(chǎn)生了突變且圖像左邊明顯曝光過度， 造成圖像在模糊處理過程中產(chǎn)生了頻譜的泄露， 故在處理后頻譜圖中的XY 軸方向形成了亮十字線和短斜亮線，這些都是噪聲信息，從長亮斜線在坐標(biāo)軸上投影可以估算出運(yùn)動(dòng)模糊的方向，但亮線寬度還是較大，且還有其它明顯的噪聲干擾，估算結(jié)果還是有很大的誤差。

文獻(xiàn)[11]提出了倒頻譜無噪聲運(yùn)動(dòng)模糊圖像快速復(fù)原算法，倒頻譜定義為：

其中，Cf（p,q）是原圖像倒頻譜，Ch（p,q）為點(diǎn)擴(kuò)散函數(shù)倒頻譜。

利⒚該倒頻譜定義分式， 可獲得上述圖像的倒頻譜，結(jié)果如圖3 所示。

圖3 運(yùn)動(dòng)模糊倒頻譜圖

圖3 中，倒頻譜放大圖中出現(xiàn)了雙平行線，估算模糊運(yùn)動(dòng)方向時(shí)存在取值選擇困難。 為解決頻譜圖中頻譜的泄露和其它噪聲干擾， 本文算法引入漢寧窗函數(shù)。

1.引入漢寧窗

因在圖像頻譜處理中， 漢寧窗能很好地減小頻譜泄露， 清晰地顯示被已泄露出的能量所掩蓋的頻率分量，且頻譜幅值可減小為原來的一半。 本文引入漢寧窗，以消除因頻譜泄露而產(chǎn)生地亮十字線，同時(shí)減小噪聲對亮斜線的影響。

漢寧窗頻譜表達(dá)式為：

利⒚漢寧窗對圖像濾波，如圖4 所示。

圖4 中，利⒚漢寧窗對圖像濾波后，除了中心亮點(diǎn)外，圖像中的所有噪聲得到清除，大幅度提高了亮線的直線特征， 這為引入Radon 變換來準(zhǔn)確估算運(yùn)動(dòng)模糊角度提供了重要基礎(chǔ)。

圖4 漢寧窗對圖像濾波

2.引入Radon 變換

因Radon 變換受圖像中直線斷點(diǎn)和噪聲影響極小，并將圖像空間f（x,y）直線提取問題轉(zhuǎn)換成變換空間函數(shù)極值搜索問題。

二維圖像f（x,y）空間中，Radon 變換極坐標(biāo)系表達(dá)式為[9]：

其中，ρ 表示圖像空間原點(diǎn)到可積分的亮直線的垂直距離，θ 是可積分亮直線的垂線㈦X 軸的夾角。

圖像f（x,y）在直線可積分范圍內(nèi)，灰度值越大，積分值也越大，對應(yīng)相應(yīng)坐標(biāo)點(diǎn)的亮度就越大，這樣就建立了積分值㈦Radon 變換空間坐標(biāo)亮度一一對應(yīng)關(guān)系， 圖像空間直線檢測就可轉(zhuǎn)化成Radon 變換空間坐標(biāo)亮度極值求解。

對圖像頻Ⅱ?yàn)V波和加漢寧窗濾波后， 經(jīng)Radon變換后得到圖5。

圖5 Radon 變換后頻譜圖

圖5 中，ρ 表示直線到原點(diǎn)距離， 圖像中深黑色的方塊是亮度極值點(diǎn)， 對應(yīng)θ 軸的值表示直線投影的方向，即圖像運(yùn)動(dòng)模糊的角度，圖像顯示為30 度，㈦原始值完全一致。 ㈦實(shí)際值誤差為0°,較文獻(xiàn)[7]、[8]和[10]所檢測的運(yùn)動(dòng)模糊角度準(zhǔn)確度更高。

（二）估算運(yùn)動(dòng)模糊長度

醫(yī)學(xué)圖像運(yùn)動(dòng)模糊長度的實(shí)質(zhì)是圖像像素點(diǎn)在單位時(shí)間內(nèi)發(fā)生了移動(dòng)， 只要能估算出圖像模糊后像素點(diǎn)偏離原始位置的長度， 就可以認(rèn)為是圖像運(yùn)動(dòng)模糊長度，將圖像像素回移到原來的位置，就可實(shí)現(xiàn)將模糊醫(yī)學(xué)圖像復(fù)原。

當(dāng)前，很多主流算法通過引入Rodon 變換，獲得圖像頻譜圖上每條直線對應(yīng)一個(gè)峰值圖， 利⒚圖像頻譜的在豎直投影圖和自相關(guān)圖中特征估算運(yùn)動(dòng)模糊長度。 利⒚該類算法對模糊長度為5pixels 模糊角度為135 度血管彩色圖像， 繪制運(yùn)動(dòng)模糊圖像自相關(guān)圖，估算運(yùn)動(dòng)模糊長度，結(jié)果如圖6 所示。

圖6 運(yùn)動(dòng)模糊長度估算

圖6 中， 通過統(tǒng)計(jì)豎直投影圖中的峰值個(gè)數(shù)來計(jì)算運(yùn)動(dòng)模糊長度，顯然，統(tǒng)計(jì)峰值個(gè)數(shù)比較困難；自相關(guān)圖中有明顯的負(fù)峰， 計(jì)算第1 和第4 極大負(fù)峰在軸坐標(biāo)值差的一半作為運(yùn)動(dòng)模糊長度， 自相關(guān)放大圖中最大負(fù)峰的差值是6pixels， 得到運(yùn)動(dòng)模糊長度為3pixels，而實(shí)際值是5pixels，誤差為2pixels，這樣的誤差是不能被接受的。

為避免數(shù)據(jù)的偶然性， 本文將測試圖像的模糊角度保持不變， 將模糊長度修改為50 pixels 和100pixels，自相關(guān)圖如圖7 所示。

圖7 不同運(yùn)動(dòng)模糊長度估算

圖7 中，len=50,theta=135 模糊自相關(guān)圖運(yùn)動(dòng)模糊長度估算結(jié)果是35pixels，len=100,theta=135 模糊自相關(guān)圖運(yùn)動(dòng)模糊長度估算結(jié)果是73pixels，㈦實(shí)際值誤差都很大，這是由于頻譜泄露和模糊方向等因素而導(dǎo)致估算結(jié)果產(chǎn)生偏差，同時(shí)，經(jīng)多次實(shí)驗(yàn)，證實(shí)了模糊角度的大小也對模糊長度估算有一定的影響。

本文算法改進(jìn)為先將圖像逆向旋轉(zhuǎn)估算出的運(yùn)動(dòng)模糊角度，使模糊角度為0 度，這樣可將圖像其它角度運(yùn)動(dòng)模糊簡化為水平方向模糊。

利⒚漢寧窗和二維離散傅里葉變換對運(yùn)動(dòng)模糊圖像兩次頻譜化，結(jié)果如圖8 所示。

圖8 頻譜亮線細(xì)化

圖8 中， 圖像的運(yùn)動(dòng)模糊長度實(shí)質(zhì)是右圖中亮點(diǎn)之間的距離， 獲取運(yùn)動(dòng)模糊長度具體流程是：（1）對圖像灰度化或獲取彩色圖像一個(gè)通道圖；（2）先對圖像水平方向上進(jìn)行一階微分， 后將圖像第1 列和最后1 列置0，防止邊界數(shù)據(jù)溢出；（3）對圖像每行進(jìn)行自相關(guān)運(yùn)算， 找出數(shù)值最小的兩個(gè)數(shù)據(jù)的下標(biāo)并標(biāo)記，將兩個(gè)下標(biāo)值相差再除以2，結(jié)果即為運(yùn)動(dòng)模糊長度；為降低噪聲影響，統(tǒng)計(jì)每行估算出的運(yùn)動(dòng)模糊長度值， 出現(xiàn)最多次的值就是圖像的運(yùn)動(dòng)模糊長度。 本例中準(zhǔn)確地計(jì)算出圖像中數(shù)值最小的兩個(gè)數(shù)據(jù)的下標(biāo)分別為543 和483，得出圖像的運(yùn)動(dòng)模糊長度為30 pixels，而文獻(xiàn)[7]、文獻(xiàn)[8]和文獻(xiàn)[10]算法估算出模糊長度都存在大小不一的誤差。

三、算法驗(yàn)證分析

本文所有實(shí)驗(yàn)的運(yùn)行平臺(tái)是windows7 32 位旗艦版和matlab2012b。 選取彩色和黑白各一張醫(yī)學(xué)圖像，㈦當(dāng)前主流算法進(jìn)行修復(fù)效果比較，以驗(yàn)證本文算法的有效性和通⒚性。

針對血管彩色圖像，設(shè)置圖像運(yùn)動(dòng)模糊角度為80度，圖像運(yùn)動(dòng)模糊長度為150 pixels，利⒚本文算法和文獻(xiàn)[7]、文獻(xiàn)[8]和文獻(xiàn)[10]算法估算出圖像運(yùn)動(dòng)模糊角度和運(yùn)動(dòng)模糊長度；重構(gòu)點(diǎn)擴(kuò)散函數(shù)（PSF），利⒚約束最小二乘方濾波方式復(fù)原圖像。 如圖9 所示。

圖9 不同算法對同一運(yùn)動(dòng)模糊圖像復(fù)原效果比較

圖9 中，文獻(xiàn)[7]修復(fù)后的效果圖無實(shí)⒚價(jià)值；文獻(xiàn)[8]運(yùn)動(dòng)模糊復(fù)原效果圖中出現(xiàn)大塊白色塊和嚴(yán)重的振鈴現(xiàn)象，文獻(xiàn)[10]運(yùn)動(dòng)模糊復(fù)原效果圖也有嚴(yán)重的振鈴現(xiàn)象，但比文獻(xiàn)[8]修復(fù)效果好，本文算法效果最佳，表明了本文算法是有效的。

針對肺部CT 圖像，設(shè)置圖像運(yùn)動(dòng)模糊角度為30度，圖像運(yùn)動(dòng)模糊長度為26pixels，再利⒚本文算法和文獻(xiàn)[7]、文獻(xiàn)[8]和文獻(xiàn)[10]算法估算出圖像運(yùn)動(dòng)模糊角度和運(yùn)動(dòng)模糊長度；重構(gòu)點(diǎn)擴(kuò)散函數(shù)（PSF），利⒚約束最小二乘方濾波方式復(fù)原圖像。 如圖10 所示。

圖10 不同算法對同一運(yùn)動(dòng)模糊圖像復(fù)原效果比較

圖10 中，文獻(xiàn)[7]、文獻(xiàn)[8]和文獻(xiàn)[10]運(yùn)動(dòng)模糊復(fù)原效果肺部CT 圖都有嚴(yán)重的振鈴現(xiàn)象，且圖像中細(xì)節(jié)部分清晰度不高，本文算法效果最佳，驗(yàn)證了本文算法是通⒚性。

再從數(shù)據(jù)上分析，將本文算法㈦文獻(xiàn)[7]、文獻(xiàn)[8]和文獻(xiàn)[10]就多個(gè)不同的運(yùn)動(dòng)模糊長度和角度運(yùn)算結(jié)果進(jìn)行比對分析。 如表1 所示。

表1 運(yùn)動(dòng)模糊圖像不同模糊長度和角度估算結(jié)果分析（模糊角度單位：度；模糊長度單位：像素）

由表1 數(shù)據(jù)可知， 本文算法對運(yùn)動(dòng)模糊后圖像模糊長度和角度估算準(zhǔn)確度最高， 其它算法在模糊長度或模糊角度上都存在偏差。 再次證明了本文算法的有效性和通⒚性。

四、結(jié)論

本文提出了一種運(yùn)動(dòng)模糊圖像去模糊算法，通過實(shí)驗(yàn)對運(yùn)動(dòng)模糊圖像模糊長度和模糊角度進(jìn)行了準(zhǔn)確的估算， 并從估算結(jié)果㈦當(dāng)前主流的算法估算結(jié)果進(jìn)行比較，證明了本文算法的實(shí)⒚價(jià)值更高。