圖像復原中的模糊參數(shù)估計

王偉，鄭津津，劉星，周洪軍，沈連婠

(1.中國科學技術大學精密機械與精密儀器系，安徽合肥 230027;2.中國科學技術大學國家同步輻射實驗室，安徽合肥 230029)

運動模糊大多數(shù)情況下是由相機抖動造成的，很多情況下沒有足夠的光照來避免使用長曝光時間，不可避免地得到拍攝結(jié)果不滿意的模糊圖像.從單獨一張運動模糊圖像恢復出一張非模糊圖像在數(shù)字圖像處理領域一直是個基礎性研究問題.如果模糊核(或者叫點擴散函數(shù))是線性移不變的話，問題就簡化成了圖像去卷積.圖像去卷積又可以分成盲恢復和非盲恢復.在非盲恢復的情況下，運動模糊核認為是已知或者已經(jīng)計算得到的，剩下的任務就是估計出原始圖像.傳統(tǒng)方法有維納提出的維納濾波［1］、Richardson 和 Lucy共同提出的 Lucy-Richardson濾波［2］，它們雖然年代久遠，但是由于簡單實用，仍然在圖像恢復中廣泛運用.盲恢復情況下，問題變得更加病態(tài)化，因為模糊核和原始圖像都是未知的.自然圖像的結(jié)構(gòu)和模糊核形狀的多樣性使得先驗分布很容易就過擬合或者欠擬合.盲去卷積是個更具挑戰(zhàn)性和病態(tài)的難題，因為模糊核未知.一些技術利用額外的輸入能夠使問題更易處理，比如多幅圖像.Rav-Acha利用2幅運動模糊圖像的信息［3］，Yuan用一副含噪和一副模糊圖像來優(yōu)化低光照條件下的拍攝效果;其他去模糊系統(tǒng)充分利用了額外特定的硬件設備［4］，Ben-Ezra和Nyar通過一個低分辨率的相機和一個高分辨率的相機來估計模糊核［5］;Rasker拍動相機的快門，使其在曝光時一會兒開一會兒關，來使高空間和高頻率的信息損失減少［6］.這些方法中，目標的運動路徑必須由用戶指定.

在直線運動模糊圖像參數(shù)確定上，目前主要分為頻域法和空域法兩大類.頻域法中主要方法是由Canon提出的檢測模糊圖像頻譜圖中黑色條紋個數(shù)的方法［7］;空域法中較有影響力的是Yitzhaky提出的計算模糊圖像的自相關函數(shù)法［8］.本文所做工作就是估計出運動模糊參數(shù)，使盲恢復問題轉(zhuǎn)化為非盲恢復來解決.

在旋轉(zhuǎn)運動模糊參數(shù)的確定上，國內(nèi)外的研究并不多，一般是事先默認已獲得了退化參數(shù)，然后利用bresenham畫圓算法提取像素點并進行重排列［9］，再運用直線運動模糊的參數(shù)估計方式進行處理.

本文在分析運動模糊特性的基礎上，針對不同運動模糊方式，提出了由粗及精的自適應參數(shù)估計方法，為后續(xù)圖像復原工作提供模糊參數(shù).

1 退化圖像的數(shù)學原理

1.1 直線運動模糊退化原理

對于運動模糊圖像來說，它每一點的灰度值是在運動參數(shù)所確定的路徑上各景物點灰度值的線性迭加.根據(jù)運動模糊降值的特點，通過一定算法構(gòu)造出與退化原理相似的恢復模型，使得模糊圖像盡量清晰，還原出那些被模糊掉的有用信息，無論在航空航天、道路交通、醫(yī)學圖像、工業(yè)控制等領域都是有重要意義的.首先要分析的是圖像降值退化的過程，圖1為圖像退化過程的示意圖.

圖1 圖像退化過程模型Fig.1 Model of image degrading

假定成像系統(tǒng)是線性位移不變系統(tǒng)，則獲取的圖像g(x，y)表示為

式中:f(x，y)表示理想沒有退化的圖像，g(x，y)表示退化后的圖像，h(x，y)表示系統(tǒng)退化核函數(shù)(點擴散函數(shù))，*表示卷積，n(x，y)表示加性噪聲.

對式(1)兩邊進行傅里葉變換，得

式中:H(u，v)稱為系統(tǒng)的傳遞函數(shù)，從頻率域角度看，它使圖像退化，因而反映了成像系統(tǒng)的性能.

先討論直線運動模糊狀況:設圖像f(x，y)在平面上相對運動，x0(t)、y0(t)分別是物體在x、y方向上的運動分量，t是運動時間，T是曝光時間，則退化圖像g(x，y)是T時間內(nèi)的積分:

為方便討論，假設僅是x方向的運動，且令曝光時間T內(nèi)的總移動量為a，則

所以有

相應的系統(tǒng)退化函數(shù)為

這表明，系統(tǒng)退化函數(shù)在運動方向上為sinc函數(shù)，具有帶狀調(diào)制的外觀，h(x)和H(u)的函數(shù)圖像如圖2所示.

圖2 勻速直線運動模糊點擴散函數(shù)和傳遞函數(shù)Fig.2 PSF and transfer function of line motion blur

由圖2可以看出，運動模糊的傳遞函數(shù)H(u)在u=n/a處為零，其中n為整數(shù).

1.2 旋轉(zhuǎn)運動模糊退化原理

再來分析一下旋轉(zhuǎn)運動模糊的退化機理.若只考慮某一條圓弧上的圖像像素模糊退化模型，該圓弧即為一模糊路徑，若該圓弧半徑為r，曝光時間為T，旋轉(zhuǎn)角速度為ω，則模糊角度θ=ωT，則退化圖像g(x，y)是T時間內(nèi)的積分:

表示為極坐標形式為

式中:x0(t)=rcos(ωt)，y0(t)=rsin(ωt).

令 l=rθ，s=rωt，并將 r表示為下標形式，則有

整理得

令 ar=2πT，則

這就是旋轉(zhuǎn)模糊運動的數(shù)學模型.

2 算法描述

運動模糊參數(shù)估計的關鍵是建立正確的模糊系統(tǒng)模型.

2.1 直線運動模糊參數(shù)估計算法

模糊方向和模糊尺度是直線運動模糊的2個重要參數(shù).

由于傳輸函數(shù)零點的存在，可推知模糊圖像的頻譜上有一系列的平行暗條紋，并且條紋的位置一一對應于各零點.由圖像的運動分析可知運動方向垂直于這些暗條紋.Paul Hough于1962年提出的Hough變換屬于特征提取技術，最初只是用于二值圖像直線檢測，后來擴展到圓等任意形狀的檢測.Hough變換的基本原理在于利用點與線的對偶性，將原始圖像空間中的給定曲線變?yōu)閰?shù)空間的一個點.這樣就將原始圖像中給定曲線的檢測問題轉(zhuǎn)化為尋找參數(shù)空間中的峰值問題［10］.

運用邊緣檢測算子可以檢測出灰度或者二值圖像的邊界，在經(jīng)過實驗比較以后，發(fā)現(xiàn)canny算子能夠很好地檢測出運動模糊頻譜的邊緣，它檢測的方法是使用2個閾值來尋找圖像梯度的局部最大值［11］.

基于上述分析，對原圖像進行傅里葉變換所得的圖像，取其模值的對數(shù).由于傅里葉變換的疊加性，同樣可見在垂直于θ的方向有2條亮線，對圖像進行canny算子邊緣提取，得到二值化的2條直線如圖3所示.再利用Hough變換檢測直線的方向即可得運動模糊圖像的相對運動方向θ，兩直線間的距離即為2倍的模糊尺度 s.求得的點擴散函數(shù)h(x，y)如圖4 所示.

圖3 Canny算子邊緣提取結(jié)果Fig.3 Edge extracted by canny operator

圖4 直線運動模糊的點擴散函數(shù)Fig.4 PSF path of linear motion blur parameter

鄒謀炎在文獻［12］中提出誤差－參數(shù)分析法，當參數(shù)從大到小變化時，誤差參數(shù)曲線在真正參數(shù)附近的變化率顯著變小，據(jù)此可得參數(shù)的估計值.

根據(jù)上述2種思想，可以根據(jù)Hough變換的方法來確定退化參數(shù)的大概范圍，然后再利用誤差－參數(shù)法來進一步確定，根據(jù)圖像恢復質(zhì)量的優(yōu)劣來確定最優(yōu)解，以下是程序算法步驟:

1)首先對直線運動模糊圖像進行對數(shù)傅里葉變換得到G1;

2)用canny算子提取G1邊緣得到G2，通過Hough變換檢測G2圖中直線的方向，即為模糊方向θ;

3)對G2進行θ角度旋轉(zhuǎn)，然后將所有像素向水平軸上投影，得到曲線相鄰極小值之間的距離即為條紋寬度s;

4)將s±5作為搜索范圍，選定迭代步長.for i=1:10

由參數(shù)α產(chǎn)生點擴散函數(shù)h，用維納濾波復原算法，由h和觀測圖像y確定圖像x，

5)作出E-α曲線，由此判定實際點擴散函數(shù)的參數(shù)值α，進而獲得實際的點擴散函數(shù).

2.2 旋轉(zhuǎn)運動模糊參數(shù)估計算法

旋轉(zhuǎn)運動模糊的2個參數(shù)為旋轉(zhuǎn)運動中心和旋轉(zhuǎn)運動角度.同樣地，根據(jù)邊緣提取的方法，提取旋轉(zhuǎn)運動模糊的模糊圓弧，然后根據(jù)3種不同的方法［13］來找圓弧，估計出旋轉(zhuǎn)運動的圓心，并用聚類方法確定最后的估計結(jié)果，最后根據(jù)bresenham畫圓結(jié)合極坐標轉(zhuǎn)換方法算出旋轉(zhuǎn)運動角度，從而確定退化模糊參數(shù).

具體操作如下:

1)首先對旋轉(zhuǎn)運動模糊圖像進行邊緣提取.

2)挑出規(guī)則的候選圓弧段，用三點找圓心法、Hough變換法和最小二乘擬合3種方法.分別對每一段候選圓弧進行擬合并求出圓心參數(shù)，最后對所有得到的參數(shù)進行聚類并用類心作為旋轉(zhuǎn)中心的估計結(jié)果.

3)將圖像由極坐標系轉(zhuǎn)換為直角坐標系，利用bresenham畫圓算法提取像素點，對每行相對空間不變模糊的情況下進行直線運動模糊估計確定s，s/r作為θ的初始值，根據(jù)θ來確定退化的點擴散函數(shù).

4)將θ±5作為搜索范圍，選定迭代步長0.5°.for i=1:10

由參數(shù)θ產(chǎn)生點擴散函數(shù)h，用維納濾波復原算法，由h和觀測圖像y確定圖像x，

5)作出E-α曲線，由此判定實際點擴散函數(shù)的參數(shù)值α，進而獲得實際的點擴散函數(shù).

3 實驗結(jié)果與分析

3.1 直線運動模糊參數(shù)估計仿真實驗

為了驗證算法的效果，選取一張256×256的8位灰度級的lena圖像作為測試圖像.直線運動模糊選取的參數(shù)為模糊角度為70°，模糊尺度為30個像素.圖5顯示的是實際模糊參數(shù)為(30，70°)的直線運動模糊Lena圖像及其頻譜.圖6是不同直線運動模糊參數(shù)恢復出的圖像的PSNR值曲面圖.

圖5 直線運動模糊圖像及其對數(shù)頻譜Fig.5 Linear motion blur image with frequency spectrum

圖6 不同直線運動模糊參數(shù)恢復出的PSNR曲面Fig.6 PSNR curved surface restored by different linear motion blur parameters

表1顯示的是針對不同直線運動模糊參數(shù)，本文算法和Yizhaky提出的傳統(tǒng)頻域方法識別得到的估計結(jié)果對比.

表1 直線運動模糊參數(shù)隨機抽測樣本Table 1 Results of parameter detection

圖7為利用本算法檢測出模糊參數(shù)后進行維納濾波恢復后的圖像.

圖7 s=30，q=70運動模糊及其恢復圖像Fig.7 Linear motion blur image with s=30，q=70

3.2 旋轉(zhuǎn)運動模糊參數(shù)估計仿真實驗

旋轉(zhuǎn)運動模糊選取的仿真實驗參數(shù)為圓心在圖片中心即(128，128)，旋轉(zhuǎn)角度為 20°.

圖8是實際模糊參數(shù)為((128，128)，20°)的旋轉(zhuǎn)運動模糊圖像及其頻譜.圖9是canny算子進行邊緣提取后的結(jié)果.圖10是不同旋轉(zhuǎn)角度恢復出的圖像的PSNR值曲線圖(20°對應最高PSNR，即為估計出的旋轉(zhuǎn)運動角度).表2給出的是不同擬合方法得到的旋轉(zhuǎn)圓心坐標.圖11是利用表2估計出的旋轉(zhuǎn)運動圓心坐標和圖10估計出的旋轉(zhuǎn)角度，最后按極坐標重新排列得到的恢復結(jié)果圖.

圖8 旋轉(zhuǎn)運動模糊圖像及其頻譜Fig.8 Rotation blur image and its frequency spectrum

圖9 Canny算子提取結(jié)果Fig.9 Edge extracted by canny operateor

圖10 不同旋轉(zhuǎn)運動模糊參數(shù)恢復出的PSNR曲線Fig.10 PSNR curve corresponding to different rotate motion blur parameters

表2 旋轉(zhuǎn)運動模糊圓心隨機抽測樣本Table 2 Results of parameter detection

圖 11 ((128，128)，20°)運動模糊及其恢復圖像Fig.11 Rotation blur image with((128，128)，20°)

以上2組測試數(shù)據(jù)表明，本文算法在大多數(shù)情形下都能準確估計出運動模糊參數(shù)，相比傳統(tǒng)意義上的頻域方法效果更好，而且步驟也并不復雜多少.對于大模糊尺度下如何進行恢復是個病態(tài)問題，原因是大模糊尺度(如45個像素以上模糊)會使圖像丟失很多信息，諸如邊界像素等.

4 結(jié)束語

通過對不同運動形式的運動模糊圖像的退化模型進行數(shù)學分析，根據(jù)不同運動模糊圖像的特點，提出了自適應的模糊參數(shù)檢測方法.通過對多幅圖像不同運動模糊參數(shù)測試的實驗證明，該方法可以便捷有效地估計出中小尺度運動模糊圖像的模糊參數(shù)，能真實反映實際運動模糊參數(shù)，并且恢復效果良好.要在大模糊尺度下進行參數(shù)估計，只需改變參數(shù)變化范圍.但不可避免地，由于大模糊尺度下的運動模糊，導致信息的過多損失，不僅計算量會增加，而且估計出的模糊參數(shù)結(jié)果遠遠不如小模糊尺度準確，有興趣的讀者可以利用本文算法一試，這是以后作者要改進和解決的一個方向.

［1］WIENER N.Extrapolation，interpolation，and smoothing of stationary time series［M］.New York:Wiley，1949:2-15.

［2］WILLIAM H R.Bayesian-based iterative method of image restoration［J］.Journal of the Optical Society of America，1972，62(1):55-59.

［3］RAVACHA A，ANDPELEG S.Two motion blurred images are better than one［J］.Pattern Recognition Letters，2005，26:311-317.

［4］YUAN L，SUN J，QUAN L，et al.Image deblurring with blurred noisy image pairs［J］.ACM Trans Graphics，2007，26(3):1-10.

［5］BENEZRA M，ANDNAYAR，S.K.Motion-based motion deblurring［J］.IEEE Transaction on Pattern Analysis and Machine Intelligence，2004，26(6):689-698.

［6］RASKAR R，AGRAWAL A，ANDTUMBLIN J.Coded exposure photography:motion deblurring using fluttered shutter［J］.ACM Transactions on Graphics，2006，25(3):95-804.

［7］HUNT B R，CANNON T M.Nonstationary assumptions of Gaussian models of images［J］.IEEE Transactions on Systems，Man，and Cybern，1976，2:876-882.

［8］YITZHAKY Y，KOPEIKA N S.Identification of blur parameters from motion blurred images［J］.Graphical Models and Image Processing，1997，59(5):321-332.

［9］JACK E.BRESENHAM.Algorithm for computer control of a digital plotter［J］.IBM Systems Journal，1965，4(1):25-30.

［10］DUDA R O，HART P E.Use of the hough transformation to detect lines and curves in pictures［J］.Communications of the ACM，1972，15:11-15.

［11］CANNY J.A computational approach to edge detection［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，1986，8(6):679-698.

［12］鄒謀炎.反卷積和信號復原［M］.北京:國防工業(yè)出版社，2001:200-202.

［13］ALDRICH J.Doing least squares:perspectives from Gauss and Yule［J］.International Statistical Review，1998，66(1):61-81.