反卷積復原算法

1、幾種圖像復原方法的對比一、Richardson-Lucy算法R-L算法是目前世界上應用最廣泛的函數(shù)恢復技術之一，它是一種迭代方法。MATLAB提供的deconvlucy（）函數(shù)還能夠用于實現(xiàn)復雜圖像重建的多種算法中，這些算法都基于Lucy-Richardson最大化可能性算法。R-L算法是一種迭代非線性復原算法，它是從最大似然公式推導出來的，圖像用泊松分布加以模型化的。當下面這個迭代收斂時模型的最大似然函數(shù)就可以得到一個令人滿意的方程：其中，*代表卷積，代表相關，代表未退化圖像的估計，g和h和以前定義一樣。在IPT中，L-R算法由名為deconvlucy的函數(shù)完成的。deconvlucy()函

2、數(shù)的調用格式：J=deconvlucy(I，PSF，NUMIT，DAMPAR，WEIGHT)。其中，I表示輸入圖像，PSF表示點擴散函數(shù)。其他參數(shù)都是可選參數(shù)：NUMIT表示算法的迭代次數(shù)，默認為10次；DAMPAR是一個標量，它指定了結果圖像與原圖像I之間的偏離閾值表，默認值為0(無衰減）；WEIGHT是一個與I同樣大小的數(shù)組，它為每一個像素分配一個權重來反映其重量，表示像素加權值，默認值為原始圖像的數(shù)值。圖像復原源代碼：% Deblurring Gray Images Using the Lucy-Richardson Algorithmclcclearclose allI=imread(

3、'E:lena512color.tif'); % 彩色圖像的像素為512*512I1=rgb2gray(I); % 灰度圖像的像素為512*512 % figure,imshow(I),title('Original color image');% figure,imshow(I1),title('Original gray image');I2=I1(1:2:end,1:2:end); % 圖像的像素設置為256*256figure,imshow(I2),title('Gray Image 256*256'); PSF = fs

4、pecial('gaussian',5,5); % 點擴散函數(shù)Blurred = imfilter(I2,PSF,'symmetric','conv'); figure;imshow(Blurred);title('Gaussian Blurred');V = 0.0001;BlurredNoisy = imnoise(Blurred,'gaussian',0,V);figure;imshow(BlurredNoisy);title('Blurred & Noisy');K=size(I2)

5、;WT=zeros(K);WT(5:end-4,5:end-4)=1;J1 = deconvlucy(BlurredNoisy,PSF); % H1 = deconvlucy(BlurredNoisy,PSF,5); % 迭代5次% H1_cell=deconvlucy(BlurredNoisy,PSF,5);% H2_cell=deconvlucy(H1_cell,PSF);% H2=im2uint8(H2_cell2);J2 = deconvlucy(BlurredNoisy,PSF,5,im2uint8(3*sqrt(V); % 迭代5次J3 =deconvlucy(BlurredNoi

6、sy,PSF,15,im2uint8(3*sqrt(V);% 迭代15次J4 =deconvlucy(BlurredNoisy,PSF,25,im2uint8(3*sqrt(V);% 迭代25次J5 =deconvlucy(BlurredNoisy,PSF,40,im2uint8(3*sqrt(V);% 迭代40次J6 =deconvlucy(BlurredNoisy,PSF,20,im2uint8(3*sqrt(V),WT);% 迭代20次,加WTJ7 = deconvlucy(BlurredNoisy,PSF,40,im2uint8(3*sqrt(V),WT); % 迭代40次,加WT%

7、figure, imshow(J1);title('J1:deconvlucy(A,PSF)');% figure, imshow(H1); title('H1:Restored Image NUMIT=5');% figure,imshow(H2),title('H2:Restored Image NUMIT=15');figure, imshow(J2);title('J2:deconvlucy(A,PSF,NUMIT=5,DAMPAR)');figure, imshow(J3);title('J3:deconvluc

8、y(A,PSF,NUMIT=15,DAMPAR)');figure, imshow(J4);title('J4:deconvlucy(A,PSF,NUMIT=25,DAMPAR)');figure, imshow(J5);title('J5:deconvlucy(A,PSF,NUMIT=40,DAMPAR)');figure, imshow(J6),title('J6:deconvlucy(A,PSF,NUMIT=20,DAMPAR,WEIGHT)');figure, imshow(J7),title('J7:deconvlucy

9、(A,PSF,NUMIT=40,DAMPAR,WEIGHT)');二、維納濾波維納濾波法是由Wiener首先提出的，在圖像復原領域，由于維納濾波計算量小，復原效果好，從而得到了廣泛的應用和發(fā)展。維納濾波最開始主要應用在一維信號處理里，取得了比較不錯的效果。之后，維納濾波法也用于二維信號處理中，也取得了比較好的效果。維納濾波器尋找一個統(tǒng)計誤差函數(shù)：最小的估計。E是期望值操作符，是未退化的圖像。該表達式在頻域可表示為其中，表示退化函數(shù)表示的復共軛表示噪聲的功率譜表示未退化圖像的功率譜比率稱為信噪功率比。在IPT中維納濾波使用函數(shù)deconvwnr來實現(xiàn)的。維納濾波能最佳復原的條件是要求已知

10、模糊的系統(tǒng)函數(shù)，噪聲功率譜密（或其自相關函數(shù)），原圖像功率譜密度（或其自相關函數(shù)）。但實際上，原圖像功率譜密度（或其自相關函數(shù)）一般難以獲知，再加上維納濾波是將圖像假設為平穩(wěn)隨機場的前提下的最佳濾波，而實際的圖像通常不能滿足此前提。因此維納濾波復原算法在實際中只能獲得次最佳實施，它更多的是具有理論價值，被用作度量其他算法性能優(yōu)劣的標桿。維納濾波復原函數(shù)deconvwnr()的調用格式：J=deconvwnr(I，PSF，NCORR，ICORR)其中，I表示輸入圖像，PSF表示點擴散函數(shù)，NSR(默認值為0)、NCORR和ICORR都是可選參數(shù)，分別表示信噪比、噪聲的自相關函數(shù)、原始圖像的自相關

11、函數(shù)。輸出參數(shù)J表示復原后的圖像。維納濾波復原源代碼：% 維納濾波在圖像復原中的應用clcclearclose allI=imread('pout.tif'); % 原始圖像noise=5*randn(size(I); % randn(1,lx)表示生成1*lx的矩陣，矩陣的每個元素都是隨機數(shù)noise=noise-min(min(noise); % randn(size(I)是返回一個和A有同樣維數(shù)大小的隨機數(shù)組J=double(I)+noise; R1=wiener2(J,10 10); % 未知噪聲R2=wiener2(J,10 10,noise); % 已知噪聲分布fi

12、guresubplot(2,2,1),imshow(uint8(I);title('原始圖像');subplot(2,2,2),imshow(uint8(J);title('退化圖像');subplot(2,2,3),imshow(uint8(R1);title('盲復原');subplot(2,2,4),imshow(uint8(R2);title('非盲復原');三、正則濾波另一個線性復原的方法稱為約束的最小二乘方濾波，在IPT中稱為正則濾波，并且通過函數(shù)deconvreg來實現(xiàn)。在最小二乘復原處理中，常常需要附加某種約束條件

13、。例如令Q為f的線性算子，那么最小二乘方復原的問題可以看成使形式為的函數(shù)，服從約束條件的最小化問題，這種有附加條件的極值問題可以用拉格朗日乘數(shù)法來處理。尋找一個，使下述準則函數(shù)為最?。菏街薪欣窭嗜障禂?shù)。通過指定不同的Q，可以得到不同的復原目標。實驗結果如下：正則濾波所用的源代碼：clcclearclose allI=imread('E:lena512color.tif');subplot(321), imshow(I),title('Original Image');I1=rgb2gray(I);subplot(322),imshow(I1),title(&#

14、39;Gray Image');PSF=fspecial('gaussian',7,10);V=.01;H=imfilter(I1,PSF);BlurredNoisy=imnoise(H,'gaussian',0,V);NOISEPOWER=V*prod(size(I1);J, LAGRA = deconvreg( BlurredNoisy, PSF,NOISEPOWER);K = deconvreg(BlurredNoisy, PSF,LAGRA/10);K0=deconvreg(BlurredNoisy, PSF,LAGRA*10);subplot(

15、323),imshow(BlurredNoisy);title('A=Blurred and Noisy');subplot(324),imshow(J);title('J LAGRA=deconvreg(A,PSF,NP)');subplot(325);imshow(K);title('deconvreg(A,PSF,0.1*LAGRA)');subplot(326);imshow(K0);title('deconvreg(A,PSF,10*LAGRA)');四、盲反卷積在圖像復原過程中，最困難的問題之一是，如何獲得PSF的恰當

16、估計。那些不以PSF為基礎的圖像復原方法統(tǒng)稱為盲去卷積。它以MLE為基礎的，即一種用被隨機噪聲所干擾的量進行估計的最優(yōu)化策略。工具箱通過函數(shù)deconvblind來執(zhí)行盲區(qū)卷積。實驗如下：% 盲反卷積圖像復原clcclearclose allI=imread('E:lena512color.tif');subplot(321),imshow(I),title('Original Image');I=rgb2gray(I);subplot(322),imshow(I),title('Gray Image');PSF=fspecial('ga

17、ussian',7,10); % PSF=7x7V=.0001;BlurredNoisy=imnoise(imfilter(I,PSF),'gaussian',0,V);BlurredNoisy=double(BlurredNoisy);WT=zeros(size(I);WT(5:end-4,5:end-4)=1; % 從第五行到倒數(shù)第五行，第五列到倒數(shù)第五列全部置為1INITPSF=ones(size(PSF); % INITPSF=ones(7x7)FUN=inline('PSF+P1','PSF','P1');J,P=deconvblind(BlurredNoisy,INITPSF,5,10*sqrt(V),WT,FUN,0); % 迭代5次K,P=deconvblind(BlurredNoisy,INITPSF,10,10*sqrt(V),WT,FUN,0); % 迭代10次L,P=deconvblind(BlurredNoisy,INITPSF,20,10*sqrt(V),WT,FU