clear all;
%number of points
N=100;


%noise amplitude
e=.5;

%domain
x=linspace(0,2*pi,N);
y=linspace(0,pi,N);
[X,Y]=meshgrid(x,y);

%pause

%****
%test field F(x,y)
a1=.5;
F=sin(X).*sin(2*Y)+a1*sin(5*X).*cos(7*Y);

cmax=.95*max(abs(F(:)));


figure(1),clf;set(gcf,'Renderer','Painters');
subplot(311)
contourf(X,Y,F);caxis([-cmax cmax]);colorbar;
xlabel('$x$','Interpreter','latex','FontSize',18)
ylabel('$y$','Interpreter','latex','FontSize',18)
title('initial field $F(x,y)$','Interpreter','latex','FontSize',18)


% noise
xi=(2*rand(N)-ones(N));
xi=randn(N);

% Fn = F + noise
Fn=F+e*xi;

figure(1);
subplot(312)
contourf(X,Y,Fn);colorbar;
xlabel('$x$','Interpreter','latex','FontSize',18)
ylabel('$y$','Interpreter','latex','FontSize',18)
title('$F(x,y)$ + noise','Interpreter','latex','FontSize',18)

%%

% SVD analysis of F+noise
[Un,S,Vn]=svd(Fn);
Sdn=diag(S);

h=figure(2)
%subplot(212)
Nind=[1:10];
plot(Nind,log10(Sdn(Nind)),'or')
xlabel('mode mumber'),ylabel('log_{10}( \Sigma )')
hold on


% SVD analysis of F
[U,S,V]=svd(F);
Sd=diag(S);

figure(2)
%subplot(211)
Nind=[1:10];
plot(Nind,log10(Sd(Nind)),'*b')
xlabel('mode mumber','Interpreter','latex','FontSize',18);
ylabel('$\log_{10}( \Sigma )$','Interpreter','latex','FontSize',18)
grid on;
hold off
h=gca;
set(h,'XTick',[0:10])
h1=legend(h,'singular values of $F$','singular values of $F$+noise$~~~~~~~~~$');
set(h1,'Location','NorthOutside','Interpreter','latex','FontSize',14);

%%
%reconstruction

Fr=Sdn(1)*Un(:,1)*Vn(:,1)'+1*Sdn(2)*Un(:,2)*Vn(:,2)'+0*Sdn(3)*Un(:,3)*Vn(:,3)';

figure(1)
subplot(313)
contourf(X,Y,Fr);caxis([-cmax cmax]);colorbar;
xlabel('$x$','Interpreter','latex','FontSize',18)
ylabel('$y$','Interpreter','latex','FontSize',18)
title('reconstruction of $F(x,y)$ using first 2 svd modes','Interpreter','latex','FontSize',18)
hold off


%%

Fd1=Sdn(1)*Un(:,1)*Vn(:,1)';
Fd2=Sdn(2)*Un(:,2)*Vn(:,2)';
Fd3=Sdn(3)*Un(:,3)*Vn(:,3)';
%[norm(Fd1-sin(X).*sin(2*Y))/norm(Fd1),norm(Fd2-a1*sin(5*X)*cos(7*Y))/norm(Fd2)]


cmax=.95*max(max(abs(Fd1(:))),max(abs(Fd2(:))));

figure(4),clf;set(gcf,'Renderer','Painters');
subplot(311)
contourf(X,Y,Fd1);caxis([-cmax cmax]);colorbar;
xlabel('$x$','Interpreter','latex','FontSize',18)
ylabel('$y$','Interpreter','latex','FontSize',18)
title('svd mode \#1 of $F$+noise','Interpreter','latex','FontSize',16)
hold off
h=gca;

subplot(312)
contourf(X,Y,Fd2);caxis([-cmax cmax]);colorbar;
xlabel('$x$','Interpreter','latex','FontSize',18)
ylabel('$y$','Interpreter','latex','FontSize',18)
title('svd mode \#2 of $F$+noise','Interpreter','latex','FontSize',16)
hold off
subplot(313)
contourf(X,Y,Fd3,'LineStyle','none');caxis([-cmax cmax]);colorbar;
xlabel('$x$','Interpreter','latex','FontSize',18)
ylabel('$y$','Interpreter','latex','FontSize',18)
title('svd mode \#3 of $F$+noise (almost zero)','Interpreter','latex','FontSize',16)
hold off