%Error estimation: 

imax=max(size(Mdreq)); p=min(size(Mdreq));  

%Calculating the average: 

cmean = mean(Mdreq')';ru=ones( size(Mdreq(1,:)) ); 
Mmean=cmean*ru; Mdif=Mdreq-Mmean;  
figure(1),mesh(Mdif),view(0,90),colormap('mastcolors'),colorbar,
title('Dev. from mean')

%Performing SVD: every row must be a sample. 

[U,S,V]=svd(Mdif',0); 
%[U,S,V]=svd(Mdreq',0); %Then Mdreq' = U*S*V'!  


%Let's take approximations: 


for j=1:p 
s0=num2str(j);st=['Ratio to dev. from mean, ' s0 ' EOFs'];
s1=['U(:,1:' s0 ')'];s2=['S(1:' s0 ',1:' s0 ')'];s3=['V(:,1:' s0 ')'];
eval(['Ma=' s1 ';' 'Mb=' s2 ';' 'Mc=' s3 ';']);Mprod=(Ma*Mb*Mc')';
eval(['M' s0 '=Mprod;']); 
for i=1:imax
a=norm(Mprod(:,i))^2;b=norm(Mdif(:,i))^2;percent(j,i)=a/b; 
end 
pmean(j)=100*mean(percent(j,:));
%eval(['M' s0 '=Mprod./Mdif;']);  
%eval(['figure(' s0 '+1),mesh(M' s0 '),']);
%view(0,90),shading interp,colormap('mastcolors'),colorbar
%title(st) 
end
plot(pmean),
xlabel('Number of EOFs');ylabel('Percent');title('Percent as function of EOFs') 
