%
% Optical Imaging and Spectroscopy
%
% David J. Brady
% Duke University
% www.opticalimaging.org
%
% figures 8.5, 8.6 and 8.5
%
figure(1);
set(gcf,'color','white');
colormap gray;
%dog=double(imread('cameraman.tif'));
load flowers;
dog=dog/max(max(dog));
subplot(3,1,1);imagesc(dog);axis image; axis off;
title('(a) Original');
sigma=0.01;
N=4;
Hl=toeplitz([1/N zeros(1,255)],[ones(1, N)/N,zeros(1,256-N)]);
Hr=Hl;
dData=Hl*dog*Hr'+sigma*randn(size(dog));
dNoise=pinv(Hl)*dData*pinv(Hr)';
subplot(3,1,2);imagesc(dNoise);axis image; axis off;
title(['(b) Least squares reconstruction mse=' num2str(mean(mean((dNoise-dog).^2)),'%3.1f')]);
[al bl cl]=svd(Hl);[ar br cr]=svd(Hr);
rl=125;
rr=125;
hlpi=cl(:,1:rl)*pinv(bl(1:rl,1:rl))*al(:,1:rl)';
dNoise=hlpi*dData*hlpi';
subplot(3,1,3);imagesc(dNoise,[0 1]);axis image; axis off;
title(['(c) Truncated SVD mse=' num2str(mean(mean((dNoise-dog).^2)),'%5.2e')]);
figure(3);
set(gcf,'color','white');
colormap gray;
subplot(2,1,1);plot(diag(br),'-k');axis([0 256 0 1]);title('Singular value spectrum');
subplot(2,1,2);imagesc(ar);axis off;title('Singular vectors');
figure(5);
set(gcf,'color','white');
colormap gray;
subplot(1,2,1);imagesc(dog(100:150,100:150));axis off; axis image;
subplot(1,2,2);imagesc(dNoise(100:150,100:150));axis off; axis image;
%