%Reads in a color image, motion blurs it with linear PSF of size (1XL) (length L),
%adds noise, then deblurs it with a Wiener filter with Tikhonov parameter lambda.
%Why do we need Tikhonov regularization? Try setting lambda=0.
clear;X=imread('toucan.jpg');X=double(X);[N,M,D]=size(X);%Know D=3.
lambda=0.1;sigma=10;L=48;M1=M+L-1;H=ones(1,L);%blurred image is NXM1.
Y(:,:,1)=conv2(X(:,:,1),H);Y(:,:,2)=conv2(X(:,:,2),H);Y(:,:,3)=conv2(X(:,:,3),H);%Blur image.
Y=Y+sigma*randn(N,M1,3);%Add 0-mean 2-D white Gaussian noise to blurred image.
FY=fft2(Y,N,M1);FH=fft2(H,N,M1);%zero-pad DFT of PSF to size of blurred image Y.
W=conj(FH)./(FH.*conj(FH)+lambda^2);%2-D Wiener filter wavenumber response.
FXHAT(:,:,1)=FY(:,:,1).*W;FXHAT(:,:,2)=FY(:,:,2).*W;FXHAT(:,:,3)=FY(:,:,3).*W;XHAT=real(ifft2(FXHAT));
XHAT=(XHAT-min(min(min(XHAT))))/(max(max(max(XHAT)))-min(min(min(XHAT))));%Scale for imagesc.
Y=(Y-min(min(min(Y))))/(max(max(max(Y)))-min(min(min(Y))));%Scale for imagesc.
figure,imagesc(Y),axis off,title('Blurred image')
figure,imagesc(XHAT),axis off,title('Reconstructed image')%This is zero-padded.