EOF程序问题

学者

请教各位大神，帮我看看这三个图分别是程序的哪部分（程序是我的，图是从别人文献里粘过来的），我也想做这三个部分，运行完后的数据没有看懂。谢谢，不胜感激








a）是载荷向量图，（b）是对应的时间系数图

sst1 = load('SST.txt')
    sst = sst1'
    [m,it] = size(sst)% m = 50 年 ；it = 24 基站数
    ix =1;iy =15;  %3*7 = 21
       %time = 1644 + 0.5/12:1911+11.5/12;
    time = 1971:2015
   
   
    %   for i=1:m
   %   sst(i,:)=deseason(sst(i,:));
    %   end

% (a) get the mean for each month:
        for i=1:m % m is the number of spatial grids (m=25*18=450)
        for k=1:12  % average for each month
        aa(i,k)=mean(sst(i,k:12:it)); % n is the length of sample
        end
        end
%  (b) remove seasonal cycle
        for i=1:m
        for j=1:it
        k=mod(j-1,12)+1;
        sst(i,j)=sst(i,j)-aa(i,k);
        end
        end
        
% (4) calculate covariance matrix
     %s=cov(sst',1);  
     %sst = bsxfun(@minus,sst,sum(sst,1)/n);
   
option='eof'
     
     %option='svd'
     
     if(option == 'eof' )  
     
      
     
     C=(sst*sst')/it; %remove mean by column (i.e. time dimension, =1776)
    %[a,b] = size(C)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%   
[E,g]=eig(C);   %特征值对角阵g 特征向量矩阵E
  E=fliplr(E);  %E左右翻转
  g=rot90(g,2); %逆时针旋转180度   
  gd=diag(g);   %取主对角线元素
  PCA=E'*sst;
  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  e1=E(:,1);    f1=reshape(e1,ix,iy)
  e2=E(:,2);    f2=reshape(e2,ix,iy);
  e3=E(:,3);    f3=reshape(e3,ix,iy);
  e4=E(:,4);    f4=reshape(e4,ix,iy)
  e5=E(:,5);    f5=reshape(e5,ix,iy);
  e6=E(:,6);    f6=reshape(e6,ix,iy);
   PC1=PCA(1,:);PC2=PCA(2,:);PC3=PCA(3,:);
   
   gd(1:3)
   total=sum(gd);
   first_eof=gd(1)/total*100; %
   second_eof=gd(2)/total*100; %
   third_eof=gd(3)/total*100; %
   forth_eof=gd(4)/total*100; %
   fifth_eof=gd(5)/total*100; %
   sixth_eof=gd(6)/total*100; %
   seventh_eof=gd(7)/total*100; %
   eighth_eof=gd(8)/total*100; %
   ninth_eof=gd(9)/total*100; %
   tenth_eof=gd(10)/total*100; %
   
  
     %%%%%%%%%%%%%%%%%      
    figure(1);clf;
    subplot(3,2,1)
    contour(f1')
    title (['The First mode of EOF  '  num2str( round(first_eof*10)/10 )   '%'])
   
    subplot(3,2,3)
    contour(f2')
    title (['The second mode of EOF  '   num2str( round(second_eof*10)/10 )   '%'])
    subplot(3,2,5)
    contour(f3')
    title (['The third mode of EOF  '   num2str( round(third_eof*10)/10 )   '%'])
   
   
    subplot(3,2,2)
    plot(time,PC1)
    grid
    title ('PC1 of EOF')
    xlim([1971 2015])
    subplot(3,2,4)
    plot(time,PC2)
    grid
    title ('PC2 of EOF')
    xlim([1971 2015])
   
    subplot(3,2,6)
    plot(time,PC3)
    grid
    title ('PC3 of EOF')
    xlim([1971 2015])
     
     else
         
         
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  SVD

[u1,v,v1]=svd(sst);
   whos u1 v1
   gd_svd=diag(v);    sigma=gd_svd.^2/it; total_svd=sum(sigma);  
   sigma(1:3)
   first_svd=sigma(1)/total_svd*100;
   second_svd=sigma(2)/total_svd*100;
   third_svd=sigma(3)/total_svd*100;
   
    left1=u1(:,1);left2=u1(:,2);left3=u1(:,3);
    svd1=reshape(left1,ix,iy);
    svd2=(1)*reshape(left2,ix,iy);
    svd3=(1)*reshape(left3,ix,iy);
    %PCA_svd=sst'*u1;
    PCA_svd=v*v1';
    PC1_svd=PCA_svd(1,:);
    PC2_svd=PCA_svd(2,:);
    PC3_svd=PCA_svd(3,:);
%%%%%%%%%%%%%%%%%
    figure(2);clf;
    subplot(3,2,1)
    contour(svd1')
    title (['The First mode of SVD  '  num2str( round(  first_svd*10)/10 )   '%'])
   
    subplot(3,2,2)
    plot(time,PC1_svd)
    grid
    title ('PC1 of SVD')
    xlim([1971 2015])
   
    subplot(3,2,3)
    contour(svd2')
    title (['The second mode of SVD  '   num2str( round(  second_svd*10)/10 )   '%'])
   
   
    subplot(3,2,4)
    plot(time,PC2_svd)
    grid
    title ('PC2 of SVD')
    xlabel('year')
    xlim([1971 2015])
   
    subplot(3,2,5)
    contour(svd3')
    title (['The third mode of SVD  '   num2str( round(  third_svd*10)/10 )   '%'])
   
   
    subplot(3,2,6)
    plot(time,PC3_svd)
    grid
    title ('PC3 of SVD')
    xlabel('year')
    xlim([1971 2015])
     end