基于路徑提取的旋轉模糊圖像復原算法

羅院紅，許廷發，朱振福

(1.北京理工大學 光電學院 光電成像技術與系統教育部重點實驗室，北京100081;2.北京環境特性研究所 目標與環境光學特征國防科技重點實驗室，北京100854)

0 引言

隨著傳感器成像技術的發展，越來越多的光學成像傳感器安裝在運動平臺上。當光學傳感器隨著旋轉平臺作高速旋轉運動時，在傳感器成像的曝光時間內，由于光學傳感器與場景之間存在較大的相對運動，所獲取的圖像將會產生嚴重的旋轉模糊。模糊圖像給后續的視頻圖像跟蹤處理，包括目標檢測、識別和跟蹤帶來很大的困難，因此，有效地對旋轉模糊圖像進行復原是視頻圖像跟蹤的一項重要工作。

國內外針對旋轉模糊圖像的復原問題已開展了一些相關的研究。Sawchuk［1-2］使用極坐標變換的方法，將二維直角坐標系下的圖像變換到極坐標下，這樣直角坐標系下的圖像旋轉模糊在極坐標下就變成旋轉角θ 方向的直線運動模糊，缺點是要進行兩次坐標變換和兩次灰度插值，運算量很大而且容易丟失原始圖像信息。洪漢玉等［3-4］提取一種基于Bresenham 路徑提取的最小二乘復原算法，該算法避免了坐標變換的繁瑣運算，但是最小二乘算法運算量仍然較大，且容易產生高頻振蕩，不利于實時實現。

為了避免大量的坐標變換和灰度插值運算，本文提出一種基于模糊路徑的圖像復原方法。這種方法通過數據分離，將二維圖像旋轉模糊復原轉換為一維的反卷積復原。本方法能夠快速地復原提取后的一維數據;能夠解決旋轉角度較大時的圖像復原問題，且速度較快，能達到實時性的要求。

1 旋轉模糊退化模型

旋轉圖像模糊不同于直線運動模糊，旋轉模糊是沿路徑的徑向運動模糊。其模糊長度隨半徑的變換而變化，即離旋轉軸心越遠，模糊長度越長。原始圖像的象素灰度值沿以旋轉軸心為圓心的一系列同心圓路徑在像平面內進行模糊，稱這一系列同心圓路徑為模糊路徑。

電荷耦合元件(CCD)成像的過程其實就是像素灰度的積分過程，而圖像模糊就可看成是場景中運動軌跡上若干個點的能量在CCD 像平面上的積累，設清晰圖像為f(x，y)，旋轉模糊圖像為g(x，y)，則在曝光時間T 內，實際采集的模糊圖像g(x，y)與清晰圖像f(x，y)的關系可表述為

采用極坐標表示為

為了便于離散化處理，將弧長l =rθ 來代替角度θ 坐標，并令s=rωt，ar=rωT，則有

將(3)式沿模糊路徑分離并進行空間離散化后可得

式中:gr(i)為模糊圖像的像素灰度值;fr(i)為清晰圖像像素灰度值;Nr為圓周上像素的總數。這樣就將二維圖像的旋轉模糊問題轉化為一維的直線運動模糊問題。gr(i)即是對應半徑為r 的圖像數據，其中ar為沿模糊路徑的模糊長度。

觀察(4)式，將fr(i)周期延拓，即有fr(i)=fr(i +Nr);并定義對應于半徑為r 時的一維點擴散函數可將其寫成離散卷積的形式:

寫成矩陣形式，有g=Cf，其中C 是點擴散函數hr(i)的循環矩陣。引入噪聲項，則基于模糊圖像的退化模型為

這樣就將二維的旋轉模糊復原問題轉化為模糊路徑上一維數據的復原，而同一模糊路徑上的模糊是空間不變的，即同半徑上的圖像像素模糊長度相同，具有相同的點擴散函數。因此，為了復原圖像，可以逐一地沿圓周模糊路徑反卷積去模糊。

2 模糊路徑提取

由第1 節的理論推導可知，旋轉模糊復原首先要提取模糊路徑上的象元灰度值。傳感器成像是將場景中的能量信息記錄在二維的離散像平面上，離散平面上的像素點并不在標準的圓周上，如圖1(a)所示。一種可行的方法是用雙線性差值將離散直角坐標上的像素值轉化到圓周上，但這種方法難免導入復雜的運算，計算量大。因此本文使用最小偏差圓弧插補法來提取模糊路徑，提取的模糊路徑為近似圓周，見圖1(b).

本文結合最小偏差圓弧插補算法和圖像旋轉模糊原理［5-6］，用插補法提取1/4 圓弧的模糊路徑，根據對稱性求得其他3 個象限模糊路徑。提取的路徑如圖1(b)所示，并非理想的標準圓，但插補法提取的路徑上的點與理想圓弧的誤差小于0.5 個像素點，而圖像灰度值是高度相關法的，所以模糊路徑上的各像素點與理想圓周上點的像素值非常接近。最小偏差插補法只需要在步進時進行簡單的誤差判斷，避免了三角函數、開方等復雜的運算，大大減少了計算量。

最小偏差插補法在給定運動模式下能夠提取任意軌跡的路徑，包括直線、圓弧和非圓弧曲線，所以在圖像旋轉模糊路徑為非正規圓時，也可以用最小偏差插補法提取模糊路徑。也就是在旋轉軸心不穩定的情況下，最小偏差插補法提取路徑仍然有效，仍可以在模糊路徑上復原圖像，這種能力是極坐標變化法所不具備的。

3 旋轉模糊圖像復原模型

3.1 維納濾波復原模型

對于上述的一維圖像模糊數據可用一維的維納濾波復原?；诰矸e定理，對(5)式作傅里葉變換，轉換為頻域，有

圖像噪聲可以看作是隨機過程，圖像復原其實就是一個基于統計估計的數據復原過程。維納濾波就是找一個原圖像f 的估計值，使得它們之間的均方誤差最小，即使e2=E{(f-)2}最小。由此可推導，維納濾波得到一維復原數據的頻譜［7-8］:

圖1 基于最小偏差圓弧插補的模糊路徑提取方法Fig.1 Trail extraction based on circular minimum deviation interpolator

式中:k(u)=Sn(u)/Sf(u)，Sn和Sf分別表示噪聲和未退化圖像的功率譜，由于未退化圖像的功率譜無法精確獲取，所以常采用信噪比作為參數，k(u)=1/SNR，SNR =(Sg-Sn)/Sn，Sg為模糊圖像的功率譜，因此通過計算模糊圖像和噪聲的功率譜，即可獲維納濾波參數。實際模糊圖像復原過程中，綜合考慮算法的實時性和有效性，對于同一傳感器的圖像，往往在初始運算過程中通過計算某一幅圖像的信噪比，進而獲取維納濾波參數k，后續的圖像復原都使用該濾波參數進行復原，即采用常量維納濾波復原圖像。

3.2 對角加載算法復原

3.2.1 對角加載算法原理

對角加載算法主要是源于對逆濾波的改進。退化后的一維圖像數據可表示為g =Cf，C 就是退化矩陣，為一循環矩陣。無約束的代數法復原可得，復原數據C-1g，對循環矩陣C 對角化后可推導得，，這就是逆濾波。

由于矩陣C 具有零特征值或小特征值，C 具有零特征值時C 不可逆，具有小特征值時，求逆會帶來嚴重的震蕩。所以逆濾波會產生嚴重的病態問題。一種解決矩陣病態性的方法是對角加載技術，修正循環矩陣C［9］。修正后的矩陣為C' = C +ΔλI，可證明C'的特征值變為ξ(k)= λ(k)+ Δλ，λ(k)是原循環矩陣C 的特征值。

證明1 C'的特征值E(C')=E(C +ΔλI)=E(C)+E(ΔλI)=E(C)+ΔλE(I)，由于E(C)=λ(k)，所以E(C')=λ(k)+Δλ. 即修正后循環矩陣的特征值為λ(k)+Δλ，且C'只是改變原矩陣C對角元素的值，所以修正后的矩陣仍然為循環矩陣。

修正前的循環矩陣為

修正后的矩陣為

顯然，C'仍然是循環矩陣，所以仍然可以對角化。應用修正后的估計退化模型就可以得到f 的估值=(C+ΔλI)-1g. 其中h =［h(0)，h(1)，…，h(M -1)］是點擴散函數，比較修正前后的循環矩陣可知，經過對角加載后的點擴散函數變為h' =［h(0)+Δλ，h(1)，…，h(M-1)］，設其頻譜為H'(u). 加載后的濾波矩陣變為

微量加載量Δλ 的選取是自適應的過程，經多次仿真實驗，取Δλ =1/L 可得到較好的復原效果，L為圓弧長度，所以越往圖像邊緣，加載量越小。

3.2.2 對角加載算法的改進——梯度加載算法

但是另一方面，對角加載會影響矩陣的小特征量，使得小特征量趨向等于加載量λ，導致高頻分量加強，反映在圖像上就是復原圖像引起一定的震蕩，導致抗噪能力降低。如果在保持低復雜度的情況下，引入新的加載，在修復病態性的同時，提高算法的抗噪能力和魯棒性，將是一舉兩得的解決方案。在此實驗過程中可以發現，由于退化函數模型的特殊性，在引入對角加載的過程中，將產生針對像素點的單向平滑效應。針對此特點，只需要再引入相應有方向性的簡單約束，即可實現對圖像的平滑，而梯度算子正是這樣的一種方案。梯度加載C″ =C +Δλ·Δ，即由梯度矩陣取代原來的單位矩陣。前面分析，逆濾波是一種無約束的復原，對角加載只是改善了矩陣的病態性，并沒有給任何約束，仍然是一種無約束的復原。在f(x)的一階導數即梯度計梯度算子為 ，引入梯Δ度算子就相當于給簡單的方向性約束，即可實現圖像的平滑。

梯度矩陣為

即修正后的矩陣為

C″仍然是循環矩陣，比較可得，梯度算子加載后的點擴散函數變為

其頻譜為H″(u)，剃度加載后f 的估計值

加載后的濾波矩陣變為

通過梯度加載量的引入，既可以避免圖像復原的退化矩陣的病態性，又可以實現圖像的平滑。相對極坐標變換法和約束最小二乘法，梯度加載算法的計算量大大減少，梯度加載復原只需要進行一次退化函數的傅里葉變換，可以用快速傅里葉變換算法計算，避免了前面算法的大量乘法運算和多次傅里葉變換計算。計算量的降低對于彈上系統等實時性要求很高的系統具有重大的意義。

不管是維納濾波還是對角加載濾波復原，都是首先沿模糊路徑快速提取模糊像素的灰度值g(i)，然后在頻域內進行濾波得到估計值的頻譜(u)，再快速傅里葉逆變換得到(i)，這樣從小到大改變半徑，逐步恢復各模糊路徑上的圖像數據。由于模糊路徑上的點不能全部映射直角坐標系下的像素點，所以恢復后的圖像會出現空穴點(即未恢復的像素點)，可以使用中值濾波或鄰域求均值的方法消除空穴點。

4 實驗結果及分析

4.1 復原效果評價指標

除了直觀的視覺圖像復原效果，本文采用峰值信噪比作為客觀評價指標來評價圖像復原前后的效果。峰值信噪比PSNR 定義為

式中:Ii分別表示清晰圖像和復原圖像的灰度值，M、N 為圖像的長和寬。

4.2 極坐標下的圖像復原過程

圖2為基于極坐標變換(CRT)的圖像復原過程，圖2(a)為旋轉角度10°的模糊圖像，圖2(b)為模糊圖像在極坐標下的變換，圖2(c)為復原后的清晰極坐標圖像，圖2(d)為清晰極坐標圖像變換到直角坐標系下的圖像，是最終的復原圖像。由于作了兩次的坐標變換和兩次灰度插值，在離散的極坐標變換過程中，圖像中心和圖像邊緣采用相同的量化步長，導致圖像中心處信息冗余而圖像邊緣處信息丟失，越靠近邊緣信息丟失越多，復原效果也就越差;在極坐標逆變換過程再次插值運算，導致信息進一步丟失，所以復原效果不佳，且有一定的振蕩效應。復原圖像的算法效果(峰值信噪比PSNR)見表1.

圖2 基于CRT 的圖像復原過程Fig.2 The restoration process of rotary motion-blurred image based on CRT

表1 算法效果(PSNR)比較Tab.1 The effects of different algorithms dB

4.3 基于模糊路徑提取的圖像復原效果

與CRT 算法不同，基于模糊路徑的復原方法，是沿模糊路徑由小到大逐個去除其空間可變模糊。為了證明算法在大角度旋轉模糊條件下的適應性，在不同的旋轉角度下進行基于模糊路徑提取的維納濾波仿真實驗，表2為仿真圖像復原效果，可見即使在旋轉模糊角度較大的情況下，本文算法仍然得到較好的復原效果。這主要是因為采用了沿模糊路徑的空間自適應算法，不同于CRT 算法使用相同的點擴散函數，算法根據不同的模糊路徑長度采用不同的點擴散函數，因而能夠較好地復原模糊路徑較大時的旋轉模糊圖像。

表2 在不同旋轉角度下基于模糊路徑提取的復原效果圖Tab.2 Comparison of blurred images and restored images at different rotary angle

噪聲是影響圖像復原效果的一個重要因素，為了證明本文算法的魯棒性，模糊退化圖像疊加不同信噪比高斯噪聲進行仿真實驗，模糊角度為10°，分別疊加40 dB、30 dB 和20 dB 的加性高斯噪聲，復原效果如表3所示。在信噪比下降的情況下，復原圖像的信噪比也逐漸下降，且圖像邊緣產生了一定的振鈴效應，但仍然具有較好的復原效果。

表3 不同信噪比下的復原效果Tab.3 Restored results under different SNRs

按照(17)式計算各種算法復原圖像的峰值信噪比PNSR，見表1.

極坐標變換的方法對圖像的改善效果較差，復原圖像的信噪比明顯低于基于模糊路徑提取的方法。維納濾波由于引入平滑效應，復原圖像的視覺效果好于對角加載，但其銳度較低，峰值信噪比不如對角加載和梯度加載;對角加載銳度較高，但是其高頻加強的特性也使得噪聲得到放大，且邊緣處有震蕩條紋，主要是直接引入對角加載無約束造成的;而引入約束的梯度加載復原圖則明顯地降低了震蕩效應，且有較好的銳度和視覺效果，集合了維納濾波和對角加載的優點。

4.4 算法處理時間比較

所有算法都在Intel 酷睿i3 3.1 GHZ，VC2005的環境下仿真，不同尺寸的圖像使用不同算法的時間比較，如表4所示，算法處理時間隨圖像尺寸的增大呈指數級增長，基于模糊路徑的方法(維納濾波、對角加載和梯度加載方法)處理時間相差不大，但遠小于傳統的極坐標變換方法，且圖像越大差距越明顯。

表4 算法處理時間比較Tab.4 Processing time of different algorithms s

4.5 真實圖像復原效果

為了驗證算法的實際可行性，利用實驗室架構的旋轉模糊圖像平臺采集模糊圖像，旋轉平臺的角速度ω 可調，且CCD 的曝光時間t 已知，相乘可得到曝光時間內圖像的旋轉角度，即模糊角度θ =ωt.如圖3所示，圖3(a)為實際CCD 采集的旋轉模糊圖像，圖3(b)為其復原效果圖。但實際圖像的復原效果不如仿真效果，分析其原因，主要有兩點:1)轉臺旋轉時其旋轉軸心并不完全穩定，因而模糊路徑并不是模型中假設的標準圓;2)實際采集的圖像有噪聲的影響，且圖像信噪比無法精確得知，影響了圖像反卷積復原的效果。

5 結論

圖3 實際模糊圖像的復原效果圖Fig.3 The restored results of real blurred images

本文對徑向可變模糊圖像復原問題進行了研究。不同于以往的CRT 方法，本文不需要進行坐標變換的繁瑣運算，而是用最小偏差插補算法提取模糊路徑上的像素灰度值。針對提取后的一維模糊像素值，提出了維納濾波法、對角加載和梯度加載的方法。進行了仿真實驗，驗證了算法復原效果，有良好的抗噪能力;對實際采集的旋轉模糊圖像復原，也得到了較好的復原效果。

References)

［1］ Sawchuk A A，Space variant image restoration by coordinate transformations［J］.JOSA，1974，64(2):138 -144.

［2］ Sawchuk A A.Space variant image motion degradation and restoration［J］.Proceedings of IEEE，1972，60(7):854 -861.

［3］ 洪漢玉.成像探測系統圖像復原算法研究［D］. 武漢:華中科技大學，2004:122 -125.HONG Han-yu. Research on image restoration algorithms in imaging detection system［D］. Wuhan:Huazhong University of Science and Techonology，2004:122 -125. (in Chinese)

［4］ 洪漢玉，張天序.非零邊界旋轉運動模糊圖像的恢復算法［J］.中國圖像圖形學報，2004，9(3):265 -274.HONG Han-yu，ZHANG Tian-xu. Restoration algorithm for rotational motion blurred images with non-zeros boundary［J］. Image and Graphics，2004，9(3):265 -274. (in Chinese)

［5］ 趙巍.數控系統的插補算法及減速控制方法研究［D］. 天津:天津大學，2004:23 -31.ZHAO Wei.Research on interpolation algorithm and accelerationdeceleration controller in CNC system［D］. Tianjin:Tianjin University，2004:23 -31. (in Chinese)

［6］ Wright W E. Parallelization of Bresenham’s line and circle algorithms［J］. Proceedings of IEEE Computer graphics and Applications，1990，9(6):60 -67.

［7］ Heltrom C W. Image restoration by the method of least squares［J］.JOSA，1967，57(3):297 -303.

［8］ 鄒謀炎. 反卷積和信號復原［M］. 北京:國防工業出版社，2001:91 -111.ZOU Mou-yan. Deconvolution and signal recovery［M］. Beijing:National Defense Industry Press，2001:91 -111. (in Chinese)

［9］ 張杰，廖桂生，王玨. 對角加載對信號源數檢測性能的改善［J］.電子學報，2004，32(12):2094 -2097.ZHANG Jie，LIAO Gui-sheng，WANG Jue. Performance improvement of source number detection using diagonal loading［J］. Acta Electronica Sinica，2004，32(12):2094 -2097. (in Chinese)