clear
close all

%% read in the image
B=imread('dalai_lama.jpg');
A=rgb2gray(B);

figure
imshow(mat2gray(A))
title(['Original b/w'])

%% choose how many terms we want to keep of this SVD
[U,S,V]=svd(double(A));
s=diag(S);

n = 20;
Ur = U(:, 1:n);
Vr = V(:, 1:n);
sr = s(1:n);
Sr = diag(sr);
Ar = Ur*Sr*Vr';

%Plot the singular values
figure
subplot(2,1,1)
plot(s,'linewidth',3)
axis tight
title('Singular values')
grid
subplot(2,1,2)
plot(log(s),'linewidth',3)
title('Singular values (log)')
grid

%% Compare reduced and original images
figure
imshow(mat2gray(Ar))
title(['n=' num2str(n)])

figure
subplot(2,1,1)
imshow(mat2gray(Ar))
title(['n=' num2str(n)])
subplot(2,1,2)
imshow(mat2gray(A))
title(['n=' num2str(length(s)) ' (full)'])
