function [ex,ey,pcx,pcy,expvar,Lambda]= calSVD(xdata,ydata,N)

[xlon, xlat, xtime]=size(xdata);
xdata=reshape(xdata, [xlon*xlat xtime]); % [space, time]
[ylon, ylat, ytime]=size(ydata);
ydata=reshape(ydata, [ylon*ylat ytime]); % [space, time]

innx= find(~isnan(xdata(:,10)));
xdata=xdata(innx,:);
xdata=permute(xdata, [2 1]);
xdata=detrend(xdata,'constant');
xdata(isnan(xdata))=0;
inny= find(~isnan(ydata(:,10)));
ydata=ydata(inny,:);
ydata=permute(ydata, [2 1]);
ydata=detrend(ydata,'constant');
ydata=detrend(ydata);
xdata(isnan(xdata))=0;

%%
S=xdata;
P=ydata;
% Form the covariance matrix:
C = S'*P;

% Find eigenvectors and singular values
[U,Lambda,V] = svds(C,N);

% PC
a = S*U;
b = P*V;

% Make them clear for output
for iN=1:N
    e1(iN,:) = squeeze( U(:,iN) )';
   pcx(iN,:) = squeeze( a(:,iN) )';
    e2(iN,:) = squeeze( V(:,iN) )';
   pcy(iN,:) = squeeze( b(:,iN) )';
end

% Amount of variance explained a 0.1 pres et en %
L2=Lambda.^2;
dsum=diag(L2)/trace(L2);

for iN=1:N
   expvar(iN)=fix( ( dsum(iN)*100/sum(dsum) )*10 ) /10;
end

%% reshape zforp to 3-D physical space
ex=repmat(NaN,[xlon*xlat N]);
for i=1:N
ex(innx,i)=e1(i,:);
end
ex=squeeze(reshape(ex, [xlon xlat N]));

ey=repmat(NaN,[ylon*ylat N]);
for i=1:N
ey(inny,i)=e2(i,:);
end
ey=squeeze(reshape(ey, [ylon ylat N]));
 



转载本文请联系原作者获取授权，同时请注明本文来自王国杰科学网博客。
链接地址：https://blog.sciencenet.cn/blog-569118-628306.html

上一篇：estimating sensitivity of rainfall frequencies to vegetation
下一篇：fast runmean smoothing