ch8- 圖象復原2.ppt

1、逆濾波比較簡單，但沒有清楚地說明如何處理噪聲，而維納濾波綜合了退化函數(shù)和噪聲統(tǒng)計特性兩個方面進行復原處理.,5 維納濾波器,逆濾波方法不能完全恢復原始信號f(x,y)，而只能求出f(x,y)的一個估計值,希望找到一種方法，在有噪聲條件下，從退化圖像 g(x,y)復原出f(x,y)的估計值，該估計值符合一定的準則。,維納濾波 (Wiener filtering)=最小均方差濾波,維納濾波是最常用的圖像恢復方法,基于維納濾波的圖像恢復方法是1967年提出的,C.W. Helstrom, “Image restoration by the method of lest sqaures,” Journ

2、al of the Optical Scoiety of America, vol.57, no.3, pp.297-303, 1967.,C.W.Helstrom, This weeks citation classic, 1982 1967-1982年SCI引用超過125次.,N.Wiener, “The extrapolation, interpolation and smoothing of stationary time series”, New York: Wiely, 1949.,在均方誤差值最小的準則下得到的 稱為對f(x,y)的最小二乘方估計。,按照該準則得到的濾波器叫維納濾

3、波器。,目標：使得復原后圖像 與原始圖像 的均方 誤差最?。?因此維納濾波器又稱為最小均方差濾波器,線性濾波：尋找點擴散函數(shù)hw(x,y)，使得,則有,由Andrews和Hunt推導滿足這一要求的傳遞函數(shù)為：,這里， 是成像系統(tǒng)傳遞函數(shù)H(u,v)的復共軛； Sn(u,v) 是噪聲功率譜； Sf (u,v)是輸入圖像的功率譜。,或者：,計算退化圖像g(x,y)的二維Fourier變換G(u,v),Wiener濾波的過程:,計算理想圖像的頻譜估計,計算點擴展函數(shù)h(x,y)的二維Fourier變換H(u,v),計算退化 圖像和噪聲的功率譜Sf(u,v),Sn(u,v),計算濾波器HW(u,v),

4、求反Fourier變換,這一方法有如下特點： （1）當H（u,v）0或幅值很小時，分母不為零，不會造成嚴重的運算誤差。 （2）在信噪比高的圖像中，即Sn(u,v)Sf(u,v),維納濾波復原法特點：,如果沒有噪聲，就成為逆濾波,(3) 當理想圖像功率譜Sf (u,v)0)時 ,表明我們不可能從全是噪聲的圖像中恢復出任何有意義的信號,（4）往往未退化圖像的功率譜Sf (u,v)難以知道，用下式近似表示：,維納濾波函數(shù)：,測試維納濾波效果：,逆濾波和維納濾波的比較,維納濾波的結果非常接近原始圖像,比逆濾波要好,逆濾波和維納濾波的比較,(a)運動模糊及均值為0方差為650的加性高斯噪聲污染的圖像 (

5、b) 逆濾波的結果 (c) 維納濾波的結果 (d)-(f) 噪聲幅度的方差比(a)小1個數(shù)量級 (g)-(i) 噪聲幅度的方差比(a)小5個數(shù)量級,逆濾波(inverse filtering),線性代數(shù)恢復,圖像恢復方法,頻域圖像恢復,維納濾波(Wiener filtering),代數(shù)逆濾波,奇異值偽逆濾波。,線性代數(shù)恢復在1970s提出,B.R.Hunt, “A matrix theory proof the discrete convolution theorem,” IEEE Trans. Audio Electroacoust , vol.19, no.4, pp.285-288, 1

6、971.,B.R.Hunt, “The application of constrained least squares estimation to image restoration by digital computer,” IEEE Trans. Computers, vol.22, no.9, pp.805-812, 1973.,約束最小二乘濾波,Matrix-vector形式,投影法。,6 代數(shù)逆濾波,無約束代數(shù)逆濾波 (unconstrained restoration),已知退化圖像的向量形式g和退化矩陣H,則無約束逆濾波恢復的的圖像為:,結果,設恢復的圖像為 ，如果不考慮噪聲，

7、用 表示恢復誤差。則代數(shù)逆濾波的目的是最小化目標函數(shù),即,推導:,使導數(shù)為零,無約束代數(shù)逆濾波 (unconstrained restoration),約束最小二乘濾波 (constrained least squares restoration),令Q為f的線性矩陣算子，約束最小二乘恢復問題就是在滿足約束條件 下，使 最小化的問題。,例如選擇Q為有限差分矩陣，使得二階差分的能量最小,約束最小二乘濾波 (constrained least squares restoration),f (x, y)在(x, y)處的二階微分,用拉格朗日法求,微分，,可以用來調節(jié)以滿足約束條件。,約束最小二乘濾波

8、 (constrained least squares restoration),推導,clear F=checkerboard(8); subplot(2,2,1),imshow(F,); title(a) 原始圖像,FontSize,12) PSF=fspecial(motion,7,45); MF=imfilter(F,PSF,circular); noise=imnoise(zeros(size(F),gaussian,0,0.001); MEN=double(MF)+noise; subplot(2,2,2),imshow(MEN,); title(b) 運動模糊+高斯噪聲后的圖像,

9、FontSize,10) NSR=sum(noise(:).2)/sum(MEN(:).2); subplot(2,2,3),imshow(deconvwnr(MEN,PSF),); %matlab自帶維納濾波函數(shù)，若噪聲功率譜設為0，則為逆濾波 title(c) 逆濾波后的圖像,FontSize,10) subplot(2,2,4),imshow(deconvwnr(MEN,PSF,NSR),);%維納濾波 title(d) 維納濾波后的圖像,FontSize,10),Matlab約束濾波函數(shù)：deconvreg,奇異值（SVD）偽逆濾波 ( SVD=Singular Value Decom

10、position),把退化矩陣H進行SVD分解，并得到H的廣義逆矩陣。從而可以把逆濾波表示成迭代形式，通過控制迭代次數(shù)，可以在恢復效果較好時停止迭代，提供了一種人機交互機會。,基本思想,其中U的列向量是AAT的特征值， V的列向量是ATA的特征值； 是對角矩陣，對角線上的值i稱為H的奇異值。,退化矩陣H的SVD分解為,于是恢復圖像為,則退化矩陣H的廣義逆矩陣為,寫成迭代形式,7 圖象灰度校正和幾何畸變校正,一、灰度校正 假設理想圖象為f（x,y），由于灰度失真因 子D（x,y）的影響，實際得到的圖象為g（x,y）， g（x,y）=D（x,y）f（x,y） 灰度校正是要從畸變的圖象g（x,y）中

11、復原原 始圖象f（x,y）。一種最直接的方法是用光密度計 測量出被拍攝景物中的某一部分區(qū)域S內真實的灰 （亮）度數(shù)值C，而對應的圖象灰度為gc（x,y）， 則 代入圖象復原方程則有,三、幾何畸變校正 在圖象獲取或顯示過程中可能產生圖象的幾何失真， 如下圖所示,從上圖可以看出，幾何畸變是將無失真坐標系中函 數(shù)f（x,y）變換到另外一個坐標上，例如，原先在（x,y）點上的象素（灰度）變化到（u,v），在圖象上 反映有些位置被擠壓，而另一些位置被擴張。我們希望 找到這兩個坐標系之間的關系。,幾何基準圖像的坐標系統(tǒng)用(x, y)來表示 需要校正的圖像的坐標系統(tǒng)用(x, y)表示,設兩個圖像坐標系統(tǒng)之間

12、的關系用解析式表示,通常h1(x,y)和h2(x,y)用多項式來表示：,通常用線性畸變來近似較小的幾何畸變,更精確一些可以用二次型來近似,若基準圖像為f(x,y)，畸變圖像為g(x,y)，對于景物上的同一個點，假定其灰度不變，則,若已知h1(x,y)和h2(x,y),通常用線性畸變來近似較小的幾何畸變,零級內插,雙線性插值,則：,h1(x,y)和h2(x,y)未知,通常用已知的多對對應點來確定系數(shù)a, b,線性畸變,可由基準圖找出三個點(r1, s1), (r2, s2), (r3, s3)與畸變圖像上三個點(u1, v1), (u2, v2), (u3, v3)一一對應。,將對應點代入，有：

13、,解聯(lián)立方程組，得出6個系數(shù)。,此時,二次型畸變,有12個未知量，需要6對已知對應點,代入上式,同樣有,解方程組，得到ai, bi. 12個系數(shù)。,記做向量矩陣形式:,通常實際應用中，會取多余的對應點對，這時A不是方陣而是高矩陣，這時矩陣的逆用廣義逆矩陣來求解。高矩陣的廣義逆矩陣為,在廣告制作和計算機動畫中常常要使物體變形。,內插方法除了零級內插和雙線性內插外，還有B樣條插值和sinc插值函數(shù)內插,幾何畸變復原的一套方法也可以用于使圖像畸變的工作中。,clear f=checkerboard(24); subplot(2,2,1),imshow(f) title(a) 原始圖像,FontSiz

14、e,10) s=0.7; theta=pi/6; T=s*cos(theta) s*sin(theta) 0; -s*sin(theta) s*cos(theta) 0; 0 0 1; tform=maketform(affine,T); g1=imtransform(f,tform,nearest); subplot(2,2,2),imshow(g1) title(b) 最鄰近插值變換,FontSize,10) g2=imtransform(f,tform); subplot(2,2,3),imshow(g2) title(c) 雙線性插值變換,FontSize,10) g3=imtrans

15、form(f,tform,Fillvalue,0.5); subplot(2,2,4),imshow(g3) title(e) 修改c的背景為灰色,FontSize,10),clear f=imread(Fig0515(a)(base-with-control-points).tif); subplot(2,2,1),imshow(f) title(a) 原始圖像,FontSize,10) inputpoints=83 81;450 56;43 293;249 392;436 442; outputpoints=68 66;375 47;42 286;275 434;523 532; tform=cp2tform(inputpoints,outputpoints,projective); g1=imtransform(f,tform); subplot(2,2,2),imshow(g1) title(b) 幾何失真圖像,FontSize,10) tform=cp2tform(outputpoints,inputpoints,projective); g2=imtransform(g1,tf