wrhwww

摘要：本文以VC++ 6.0為編程工具，講述了采取逆濾波和維納濾波兩種圖像恢復(fù)算法對退化圖像的恢復(fù)實現(xiàn)過程。

　　引言

　　圖像恢復(fù)技術(shù)是圖像處理領(lǐng)域一類重要的處理技術(shù)，與圖像增強等其他基本圖像處理技術(shù)類似，該技術(shù)也是以獲取視覺質(zhì)量得到某種程度改善為目的的，所不同的是圖像恢復(fù)過程需要根據(jù)指定的圖像退化模型來完成，根據(jù)這個退化模型對在某種情況下退化或惡化了的退化圖像進行恢復(fù)，以獲取到原始的、未經(jīng)過退化的原始圖像。換句話說，圖像恢復(fù)的處理過程實際是對退化圖像品質(zhì)的提升，并通過圖像品質(zhì)的提升來達到圖像在視覺上的改善。本文以VC++作為開發(fā)工具，講述了對退化圖像進行逆濾波和維納濾波處理算法。

　　逆濾波處理

　　對圖像進行恢復(fù)處理通常需要根據(jù)一定的圖像退化模型來進行，一個簡單的通用圖像退化模型可將圖像的退化過程模型化為一個作用在原始圖像f(x,y)上的退化系統(tǒng)H，作用結(jié)果與一個加性噪聲n(x,y)的聯(lián)合作用導(dǎo)致產(chǎn)生出了退化圖像g(x,y)，表現(xiàn)為數(shù)學(xué)形式為g(x,y)=H[f(x,y)]+n(x,y)。根據(jù)上述退化系統(tǒng)H可以從給定的退化圖像g(x,y)得到原始圖像f(x,y)的一個近似結(jié)果。逆濾波處理就是其中一種無約束恢復(fù)的圖像恢復(fù)技術(shù)，其恢復(fù)過程的數(shù)學(xué)形式可表示為F (u,v)=G(u,v)/H(u,v) (u,v=0,1,…,M-1),其中F(u,v)和G(u,v)分別為圖像f(x,y)和g(x,y)的頻域變換,H(u,v)可看作是一個濾波函數(shù)。由于圖像在退化過程中存在噪聲的干擾，因此通常情況下的濾波器往往不是正好的1/H(u,v),而是關(guān)于u和v的某個非線形的恢復(fù)轉(zhuǎn)移函數(shù)M(u,v)。經(jīng)過以上的分析，圖像的退化和恢復(fù)過程（模型）大致可用下圖來表示：


　　一種簡便的恢復(fù)方法是在選取恢復(fù)轉(zhuǎn)移函數(shù)M(u,v) 時，如果u2+v2≤w2,則取值1/H(u,v)，否則為1。這樣處理雖然簡單，但是恢復(fù)后的圖像往往存在較明顯的振鈴現(xiàn)象，通常為了消除振鈴現(xiàn)象，以H(u,v)的值作為判據(jù)，如不大于d(0
　　由于恢復(fù)過程需要在頻域進行，因此需要通過二維傅立葉變換將圖像由空域變換到頻域。二維的傅立葉變換較一維傅立葉變換要復(fù)雜的多，一般采取連續(xù)2次運用一維離散快速傅立葉變換的方法來實現(xiàn)，即先沿f(x,y)的每一個x對y求變換再乘以N得到F(x,v)，完成第一步變換。然后再將得到的F(x,v)沿f(x,v)的每一個v對x求變換即可得到f(x,y)的最終變換F(u,v)，這兩步的數(shù)學(xué)表達式如下：

F(x,v)=N*[(1/N)* f(x,y)exp[-j2πvy/N]] (v=0,1,……,N-1) F(u,v)=(1/N)* F(x,v)exp[-j2πux/N] (u,v=0,1,……,N-1)

　　類似也可以得出二維離散傅立葉變換逆變換用一維變換計算的表達式：

F(x,v)= F(u,v)exp[j2πux/N] (x,y=0,1,……,N-1) f(x,y)=(1/N)* F(x,v)exp[j2πvy/N]] (y=0,1,……,N-1)

　　在分布進行一維傅立葉變換時，多采用"蝴蝶圖"的快速算法（詳見信號處理方面資料），其核心算法如下：

int N=(int)pow(2,M); file://N：序列長度（2的整數(shù)次冪） ReverseOrder(A,N); file://對空間序列進行倒序 for(int i=1;i<=M;i++){ 　int b=(int)pow(2,(i-1)); 　for(int j=0;j<=(b-1);j++) { 　　float p=(float)(pow(2,(M-i))*j*2.0*PI/(float)N); 　　for(int k=j;k<=(N-1);){ 　　　float tr=(float)(A[k+b].Re*cos(p)+A[k+b].Im*sin(p)); file://計算復(fù)數(shù)運算A*U 　　　float ti=(float)(A[k+b].Im*cos(p)-A[k+b].Re*sin(p)); 　　　A[k+b].Re=A[k].Re-tr; file://復(fù)數(shù)運算A-tr 　　　A[k+b].Im=A[k].Im-ti; 　　　A[k].Re+=tr; file://復(fù)數(shù)運算A+tr 　　　A[k].Im+=ti; 　　　k+=b*2; 　　} 　} }

　　傅立葉逆變換的同傅立葉變換比較相似，只是在計算exp[j2πvy/N]時同正變換有符號的區(qū)別，以及在計算完成后逆變換需要將值除以N，因此不難寫出一維傅立葉逆變換的實現(xiàn)代碼。在進行二維傅立葉變換將圖像由空域變換到頻域之前，首先需要通過補0的手段將點數(shù)非2的整數(shù)次冪的非正方型網(wǎng)格采樣構(gòu)造為一個長寬均為2的整數(shù)次冪的正方型網(wǎng)格：

int WM=(int)(log(W)/log(2)+1.0f); file://計算圖像寬應(yīng)為2的多少次冪 int HM=(int)(log(H)/log(2)+1.0f); file://計算圖像高應(yīng)為2的多少次冪 WM=HM=max(WM,HM); file://取二者大值 int WN=(int)pow(2,WM); file://構(gòu)造網(wǎng)格寬度 int HN=(int)pow(2,HM); file://構(gòu)造網(wǎng)格高度 for{int i=0;i;for(int j=0;j 　　if(i 　　　U[i*WN*3+j].Re=D[i*W*3+j]; file://D為圖像序列 　　U[i*WN*3+j].Im=0.0f; 　　}else file://缺位補0 　　U[i*WN*3+j].Re=U[i*WN*3+j].Im=0.0f; 　} }

　　預(yù)處理完畢后，可對構(gòu)造網(wǎng)格的每一列分別進行一維快速傅立葉變換，并將結(jié)果存放在原位置，結(jié)果乘以N，完成第一步的轉(zhuǎn)換，求得F(x,v)：

for(i=0;i 　for(int j=0;j 　　UH[j].Re=U[j*WN*3+i].Re; 　　UH[j].Im=U[j*WN*3+i].Im; 　} 　DFT_FFT(UH,HM); file://對UH進行快速離散傅立葉變換 　for(j=0;j 　　U[j*WN*3+i].Re=HN*UH[j].Re; file://N=HN 　　U[j*WN*3+i].Im=HN*UH[j].Im; 　} }

　　隨即對構(gòu)造網(wǎng)格的每一行進行傅立葉變換，得到最終的變換結(jié)果F(u,v)：

for(i=0;i 　for(int k=0;k<3;k++){ file://對24位位圖的R、G、B三分量均各自進行變換 　　for(int j=0;j 　　　UW[j].Re=U[i*WN*3+j*3+k].Re; 　　　UW[j].Im=U[i*WN*3+j*3+k].Im; 　　} 　　DFT_FFT(UW,WM); file://對UW序列進行快速離散傅立葉變換 　　for(j=0;j 　　　U[i*WN*3+j*3+k].Re=UW[j].Re; 　　　U[i*WN*3+j*3+k].Im=UW[j].Im; 　　} 　} }

　　至于二維傅立葉逆變換則基本上是上述過程的逆過程，在此就不再贅述。根據(jù)逆濾波圖像恢復(fù)的設(shè)計方案，先通過前面的二維傅立葉變換將退化圖像g(x,y)從空域變換到頻域得到G(u,v)，然后在頻域經(jīng)過恢復(fù)轉(zhuǎn)移函數(shù)M(u,v)的恢復(fù)處理并經(jīng)過二維傅立葉逆變換將結(jié)果由頻域轉(zhuǎn)換回空域，就可得到經(jīng)過恢復(fù)處理的近似原始圖像：

…… dsp.DFT_2D_FFT(m_cpBuffer+54,m_nWidth,m_nHeight,U); file://進行二維傅立葉變換 for(int i=0;i 　for(int j=0;j 　　int k=(int)(j/3); 　　D1=(float)sqrt(i*i+k*k); 　　H=1.0f/(1+(D1/D0)*(D1/D0)); file://H(u,v)=1/(1+(u2+v2)/D02)) 　　if(H>0.45f){ file://閥值 d取0.45 　　　U[i*3*WN+j].Re/=H; file://在頻域與M(u,v)相乘 　　　U[i*3*WN+j].Im/=H; 　　}else{ 　　　U[i*3*WN+j].Re*=0.6f; file://如未超過閥值則M(u,v)取常數(shù)k=0.6 　　　U[i*3*WN+j].Im*=0.6f; 　　} 　} } dsp.DFT_2D_IFFT(m_cpBuffer+54,m_nWidth,m_nHeight,U); file://進行傅立葉逆變換

　　這里的逆濾波處理算法采用的是經(jīng)過改進的恢復(fù)轉(zhuǎn)移函數(shù)M(u,v)，因此恢復(fù)后的圖像不會出現(xiàn)振鈴現(xiàn)象。以標(biāo)準(zhǔn)檢測圖像Lina為處理對象應(yīng)用以上恢復(fù)處理算法，效果如下圖所示。其中間圖像為未經(jīng)過改進的簡單算法，在胳膊和臉部存在較明顯的振鈴現(xiàn)象，而采取了改進措施的圖像則沒有任何振鈴現(xiàn)象出現(xiàn)，圖像得到了較好的恢復(fù)。

維納濾波處理

　　維納（Wiener）濾波是對退化圖像進行恢復(fù)處理的另一種常用算法，是一種有約束的恢復(fù)處理方法，其采用的維納濾波器是一種最小均方誤差濾波器，其數(shù)學(xué)形式比較復(fù)雜：

F(u,v)=[(1/H(u,v))*(|H(u,v)|2)/(|H(u,v)|2+s*[Sn(u,v)/Sf(u,v)])]*G(u,v)

　　當(dāng)s為1時，上式就是普通的維納濾波；如果s為變量，則為參數(shù)維納濾波，如果沒有噪聲干擾，即Sn(u,v)=0時，上式實際就是前面的逆濾波。從其數(shù)學(xué)形式可以看出：維納濾波比逆濾波在對噪聲的處理方面要強一些。以上只是理論上的數(shù)學(xué)形式，在進行實際處理時，往往不知道噪聲函數(shù)Sn(u,v)和Sf(u,v)的分布情況，因此在實際應(yīng)用時多用下式進行近似處理：

F(u,v)=[(1/H(u,v))* (|H(u,v)|2)/(|H(u,v)|2+K)]*G(u,v)

　　其中K是一個預(yù)先設(shè)定的常數(shù)。由此可以寫出維納濾波的實現(xiàn)代碼：

…… float K=0.05f; file://預(yù)先設(shè)定常數(shù)K dsp.DFT_2D_FFT(m_cpBuffer+54,m_nWidth,m_nHeight,U); file://轉(zhuǎn)換到頻域 for(int i=0;i 　for(int j=0;j 　　int k=(int)(j/3); 　　D1=(float)sqrt(i*i+k*k); 　　float H=1.0f/(1+(D1/D0)*(D1/D0));//H(u,v)= 1/(1+(u2+v2)/D02)) 　　U[i*3*WN+j].Re=(U[i*3*WN+j].Re*H)/(H*H+K); file://維納濾波 　　U[i*3*WN+j].Im=(U[i*3*WN+j].Im*H)/(H*H+K); 　} } dsp.DFT_2D_IFFT(m_cpBuffer+54,m_nWidth,m_nHeight,U);//返回到空域

　　對經(jīng)過退化的Lina圖像應(yīng)用維納濾波處理，可得到如右圖所示的恢復(fù)效果圖。由于維納濾波在進行恢復(fù)時對噪聲進行了處理，因此其恢復(fù)效果要比逆濾波要好，尤其是退化圖像的噪聲干擾較強時效果更為明顯。


　　小結(jié)

　　本文對比較常用的兩種圖像恢復(fù)算法逆濾波和維納濾波的實現(xiàn)過程作了較為詳細(xì)的講述，通過對圖像質(zhì)量較低的退化圖像應(yīng)用上述算法可以使圖像質(zhì)量得到一定程度的改善，在視覺上可以得到較好的改觀。類似的圖像恢復(fù)算法還有有約束最小平方恢復(fù)算法等多種，應(yīng)視具體情況靈活選擇合適的算法以獲取最佳的恢復(fù)效果。本文所述程序在Windows 98下，由Microsoft Visual C++ 6.0編譯通過。