基于路徑提取的旋轉(zhuǎn)模糊圖像復(fù)原算法

羅院紅，許廷發(fā)，朱振福

(1.北京理工大學(xué) 光電學(xué)院 光電成像技術(shù)與系統(tǒng)教育部重點(diǎn)實(shí)驗(yàn)室，北京100081;2.北京環(huán)境特性研究所 目標(biāo)與環(huán)境光學(xué)特征國防科技重點(diǎn)實(shí)驗(yàn)室，北京100854)

0 引言

隨著傳感器成像技術(shù)的發(fā)展，越來越多的光學(xué)成像傳感器安裝在運(yùn)動(dòng)平臺(tái)上。當(dāng)光學(xué)傳感器隨著旋轉(zhuǎn)平臺(tái)作高速旋轉(zhuǎn)運(yùn)動(dòng)時(shí)，在傳感器成像的曝光時(shí)間內(nèi)，由于光學(xué)傳感器與場景之間存在較大的相對運(yùn)動(dòng)，所獲取的圖像將會(huì)產(chǎn)生嚴(yán)重的旋轉(zhuǎn)模糊。模糊圖像給后續(xù)的視頻圖像跟蹤處理，包括目標(biāo)檢測、識(shí)別和跟蹤帶來很大的困難，因此，有效地對旋轉(zhuǎn)模糊圖像進(jìn)行復(fù)原是視頻圖像跟蹤的一項(xiàng)重要工作。

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

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

1 旋轉(zhuǎn)模糊退化模型

旋轉(zhuǎn)圖像模糊不同于直線運(yùn)動(dòng)模糊，旋轉(zhuǎn)模糊是沿路徑的徑向運(yùn)動(dòng)模糊。其模糊長度隨半徑的變換而變化，即離旋轉(zhuǎn)軸心越遠(yuǎn)，模糊長度越長。原始圖像的象素灰度值沿以旋轉(zhuǎn)軸心為圓心的一系列同心圓路徑在像平面內(nèi)進(jìn)行模糊，稱這一系列同心圓路徑為模糊路徑。

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

采用極坐標(biāo)表示為

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

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

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

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

寫成矩陣形式，有g(shù)=Cf，其中C 是點(diǎn)擴(kuò)散函數(shù)hr(i)的循環(huán)矩陣。引入噪聲項(xiàng)，則基于模糊圖像的退化模型為

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

2 模糊路徑提取

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

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

最小偏差插補(bǔ)法在給定運(yùn)動(dòng)模式下能夠提取任意軌跡的路徑，包括直線、圓弧和非圓弧曲線，所以在圖像旋轉(zhuǎn)模糊路徑為非正規(guī)圓時(shí)，也可以用最小偏差插補(bǔ)法提取模糊路徑。也就是在旋轉(zhuǎn)軸心不穩(wěn)定的情況下，最小偏差插補(bǔ)法提取路徑仍然有效，仍可以在模糊路徑上復(fù)原圖像，這種能力是極坐標(biāo)變化法所不具備的。

3 旋轉(zhuǎn)模糊圖像復(fù)原模型

3.1 維納濾波復(fù)原模型

對于上述的一維圖像模糊數(shù)據(jù)可用一維的維納濾波復(fù)原。基于卷積定理，對(5)式作傅里葉變換，轉(zhuǎn)換為頻域，有

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

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

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

3.2 對角加載算法復(fù)原

3.2.1 對角加載算法原理

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

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

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

修正前的循環(huán)矩陣為

修正后的矩陣為

顯然，C'仍然是循環(huán)矩陣，所以仍然可以對角化。應(yīng)用修正后的估計(jì)退化模型就可以得到f 的估值=(C+ΔλI)-1g. 其中h =［h(0)，h(1)，…，h(M -1)］是點(diǎn)擴(kuò)散函數(shù)，比較修正前后的循環(huán)矩陣可知，經(jīng)過對角加載后的點(diǎn)擴(kuò)散函數(shù)變?yōu)閔' =［h(0)+Δλ，h(1)，…，h(M-1)］，設(shè)其頻譜為H'(u). 加載后的濾波矩陣變?yōu)?/p>

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

3.2.2 對角加載算法的改進(jìn)——梯度加載算法

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

梯度矩陣為

即修正后的矩陣為

C″仍然是循環(huán)矩陣，比較可得，梯度算子加載后的點(diǎn)擴(kuò)散函數(shù)變?yōu)?/p>

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

加載后的濾波矩陣變?yōu)?/p>

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

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

4 實(shí)驗(yàn)結(jié)果及分析

4.1 復(fù)原效果評價(jià)指標(biāo)

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

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

4.2 極坐標(biāo)下的圖像復(fù)原過程

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

圖2 基于CRT 的圖像復(fù)原過程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 基于模糊路徑提取的圖像復(fù)原效果

與CRT 算法不同，基于模糊路徑的復(fù)原方法，是沿模糊路徑由小到大逐個(gè)去除其空間可變模糊。為了證明算法在大角度旋轉(zhuǎn)模糊條件下的適應(yīng)性，在不同的旋轉(zhuǎn)角度下進(jìn)行基于模糊路徑提取的維納濾波仿真實(shí)驗(yàn)，表2為仿真圖像復(fù)原效果，可見即使在旋轉(zhuǎn)模糊角度較大的情況下，本文算法仍然得到較好的復(fù)原效果。這主要是因?yàn)椴捎昧搜啬：窂降目臻g自適應(yīng)算法，不同于CRT 算法使用相同的點(diǎn)擴(kuò)散函數(shù)，算法根據(jù)不同的模糊路徑長度采用不同的點(diǎn)擴(kuò)散函數(shù)，因而能夠較好地復(fù)原模糊路徑較大時(shí)的旋轉(zhuǎn)模糊圖像。

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

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

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

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

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

4.4 算法處理時(shí)間比較

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

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

4.5 真實(shí)圖像復(fù)原效果

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

5 結(jié)論

圖3 實(shí)際模糊圖像的復(fù)原效果圖Fig.3 The restored results of real blurred images

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

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］ 洪漢玉.成像探測系統(tǒng)圖像復(fù)原算法研究［D］. 武漢:華中科技大學(xué)，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］ 洪漢玉，張?zhí)煨?非零邊界旋轉(zhuǎn)運(yùn)動(dòng)模糊圖像的恢復(fù)算法［J］.中國圖像圖形學(xué)報(bào)，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］ 趙巍.數(shù)控系統(tǒng)的插補(bǔ)算法及減速控制方法研究［D］. 天津:天津大學(xué)，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］ 鄒謀炎. 反卷積和信號(hào)復(fù)原［M］. 北京:國防工業(yè)出版社，2001:91 -111.ZOU Mou-yan. Deconvolution and signal recovery［M］. Beijing:National Defense Industry Press，2001:91 -111. (in Chinese)

［9］ 張杰，廖桂生，王玨. 對角加載對信號(hào)源數(shù)檢測性能的改善［J］.電子學(xué)報(bào)，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)