%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 
%   Find the circulation during/before the onset in the 3 regions
%   1. Mean Circulation during 0+5; -1-5 -6-10 and -11-15 before onset 
%   2. Pearson and Spearman Correlation with significvant test using
%   permutation
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear;
close all;
load start20mmonly.txt;
year=start20mmonly(:,1);
nord=start20mmonly(:,2);
cent=start20mmonly(:,3);
sud=start20mmonly(:,4);
anord=nord(:)-mean(nord(:));
acent=cent(:)-mean(cent(:));
asud=sud(:)-mean(sud(:));
% actual rainfall data file 
load indicenorth.txt;
load indicecent.txt;
load indicesud.txt;
loaddods('http://iridl.ldeo.columbia.edu/expert/SOURCES/.NOAA/.NCEP-NCAR/.CDAS-1/.DAILY/.Intrinsic/.PressureLevel/.u/P/925/VALUE%5BP%5Daverage/T/%281%20Jan%201968%29%2831%20Dec%202004%29RANGE/T/%28May-Oct%29VALUES/Y/%2850N%29%2840S%29RANGEEDGES/X/%28180W%29%28180E%29RANGEEDGES/dods');
loaddods('http://iridl.ldeo.columbia.edu/expert/SOURCES/.NOAA/.NCEP-NCAR/.CDAS-1/.DAILY/.Intrinsic/.PressureLevel/.v/P/925/VALUE%5BP%5Daverage/T/%281%20Jan%201968%29%2831%20Dec%202004%29RANGE/T/%28May-Oct%29VALUES/Y/%2850N%29%2840S%29RANGEEDGES/X/%28180W%29%28180E%29RANGEEDGES/dods');

nx=size(X,1);
ny=size(Y,1);
[x,y]=meshgrid(X,Y);
% creating rainfall daily anomalies
            r1=indicenorth;
            r2=indicecent;
            r3=indicesud;
            mdays=mean(r1,2);  % daily mean
            rmean=mdays(:,ones(1,37)); % repaeting
            r1a=r1-rmean; % verified OK
            mdays=mean(r2,2);  % daily mean
            rmean=mdays(:,ones(1,37)); % repaeting
            r2a=r2-rmean; % verified OK
            mdays=mean(r3,2);  % daily mean
            rmean=mdays(:,ones(1,37)); % repaeting
            r3a=r3-rmean; % verified OK   
 
% Daily anomalies of U and V by removing anual cycle
    for ix=1:nx
        for iy=1:ny
            u1(1,:)=u(iy,ix,:);
            mat=reshape(u1,184,37);
            mdays=mean(mat,2);  % verified OK
            zmean=mdays(:,ones(1,37));
            meanserie=reshape(zmean,37*184,1);
            ua(iy,ix,:)=u1-meanserie'; % verified OK
            
            v1(1,:)=v(iy,ix,:);
            mat=reshape(v1,184,37);
            mdays=mean(mat,2);  % verified OK
            zmean=mdays(:,ones(1,37));
            meanserie=reshape(zmean,37*184,1);
            va(iy,ix,:)=v1-meanserie'; % verified OK            
        end
    end
clear X Y T mdays rmean mat zmean meanserie start20mmonly indicenorth indicecent indicesud;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%   Finding the Precip on dates of onsets t=00 to +05
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 
tot1=0.;
tot2=0.;
tot3=0.;
tota1=0.;
tota2=0.;
tota3=0.;
for i=1:37
    tota1=tota1+sum(r1a(nord(i):nord(i)+5,i));
    tot1=tot1+sum(r1(nord(i):nord(i)+5,i));  
    tota2=tota2+sum(r2a(cent(i):cent(i)+5,i));    
    tot2=tot2+sum(r2(cent(i):cent(i)+5,i));  
    tota3=tota3+sum(r3a(sud(i):sud(i)+5,i));        
    tot3=tot3+sum(r3(sud(i):sud(i)+5,i));  
end

%%%%%%%%%% Find the circulation dates 
n1=32;
for i=1:37
    for ix=1:nx
        for iy=1:ny
            indan=(i-1)*184;
            n1=nord(i)-121; % 121= 01Jan-30Apr
            ijour=indan+n1;  % exact day in the matrix [37*184]==6808
            uanord(ix,iy,i)=mean(ua(iy,ix,ijour:ijour+5));
            vanord(ix,iy,i)=mean(va(iy,ix,ijour:ijour+5));
            unord(ix,iy,i)=mean(u(iy,ix,ijour:ijour+5));
            vnord(ix,iy,i)=mean(v(iy,ix,ijour:ijour+5));            

            n1=cent(i)-121; 
            ijour=indan+n1;
            uacent(ix,iy,i)=mean(ua(iy,ix,ijour:ijour+5));
            vacent(ix,iy,i)=mean(va(iy,ix,ijour:ijour+5));
            ucent(ix,iy,i)=mean(u(iy,ix,ijour:ijour+5));
            vcent(ix,iy,i)=mean(v(iy,ix,ijour:ijour+5));            

            n1=sud(i)-121;   
            ijour=indan+n1;
            uasud(ix,iy,i)=mean(ua(iy,ix,ijour:ijour+5));
            vasud(ix,iy,i)=mean(va(iy,ix,ijour:ijour+5));
            usud(ix,iy,i)=mean(u(iy,ix,ijour:ijour+5));
            vsud(ix,iy,i)=mean(v(iy,ix,ijour:ijour+5));            
        end
    end
end
 
% % averaging the circulations during those dates  
unordm=mean(unord,3);
vnordm=mean(vnord,3);
ucentm=mean(ucent,3);
vcentm=mean(vcent,3);
usudm=mean(usud,3);
vsudm=mean(vsud,3);

uanordm=mean(uanord,3);
vanordm=mean(vanord,3);
uacentm=mean(uacent,3);
vacentm=mean(vacent,3);
uasudm=mean(uasud,3);
vasudm=mean(vasud,3);

% load coast;
% figure(1);
% clf;
% subplot(3,2,1)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,unordm',vnordm');
% set(h,'color',[0,0,0]);
% title('Nord -1-5 dys tot')
% ylabel(tot1)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,2)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,uanordm',vanordm');
% set(h,'color',[0,0,0]);
% title('Nord -1-5 dys ano')
% ylabel(tota1)
% axis([-180 180 -40 50]);
% 
% 
% subplot(3,2,3)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ucentm',vcentm');
% set(h,'color',[0,0,0]);
% title('Cent -1-5 dys tot')
% ylabel(tot2)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,4)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,uacentm',vacentm');
% set(h,'color',[0,0,0]);
% title('Cent -1-5 dys ano')
% ylabel(tota2)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,5)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,usudm',vsudm');
% set(h,'color',[0,0,0]);
% title('Sud -1-5 dys tot')
% ylabel(tot3)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,6)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,uasudm',vasudm');
% set(h,'color',[0,0,0]);
% title('Sud -1-5 dys ano')
% ylabel(tota3)
% axis([-30 30 0 30])
% suptitle('Mean circulation at 925hPa during 00 to +05 days after onset')
% orient tall;
% eval(sprintf('print circ92start00-05.ps -dpsc'))
% 
% %pause;
fid=fopen('circ92start00-05.bin','wb');
a=fwrite(fid,unordm,'float32');
a=fwrite(fid,vnordm,'float32');
a=fwrite(fid,uanordm,'float32');
a=fwrite(fid,vanordm,'float32');
a=fwrite(fid,ucentm,'float32');
a=fwrite(fid,vcentm,'float32');
a=fwrite(fid,uacentm,'float32');
a=fwrite(fid,vacentm,'float32');
a=fwrite(fid,usudm,'float32');
a=fwrite(fid,vsudm,'float32');
a=fwrite(fid,uasudm,'float32');
a=fwrite(fid,vasudm,'float32');
fclose(fid);
%%%%%%%%%%%%%

% correlating using Spearman rank correlation  
for i=1:nx
    for j=1:ny
       s=uanord(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(anord,serie);
        if(sig>=90) 
            ruan(i,j)=r;
        else
            ruan(i,j)=-999;            
        end 
       s=unord(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(nord,serie);
        if(sig>=90) 
            run(i,j)=r;
        else
            run(i,j)=-999;            
        end 
        
       s=uacent(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(acent,serie);
        if(sig>=90) 
            ruac(i,j)=r;
        else
            ruac(i,j)=-999;            
        end 
       s=ucent(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(cent,serie);
        if(sig>=90) 
            ruc(i,j)=r;
        else
            ruc(i,j)=-999;            
        end 
        
       s=uasud(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(asud,serie);
        if(sig>=90) 
            ruas(i,j)=r;
        else
            ruas(i,j)=-999;            
        end 
       s=usud(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(sud,serie);
        if(sig>=90) 
            rus(i,j)=r;
        else
            rus(i,j)=-999;            
        end 

%%% V component
       s=vanord(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(anord,serie);
        if(sig>=90) 
            rvan(i,j)=r;
        else
            rvan(i,j)=-999;            
        end 
       s=vnord(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(nord,serie);
        if(sig>=90) 
            rvn(i,j)=r;
        else
            rvn(i,j)=-999;            
        end 
        
       s=vacent(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(acent,serie);
        if(sig>=90) 
            rvac(i,j)=r;
        else
            rvac(i,j)=-999;            
        end 
       s=vcent(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(cent,serie);
        if(sig>=90) 
            rvc(i,j)=r;
        else
            rvc(i,j)=-999;            
        end 
        
       s=vasud(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(asud,serie);
        if(sig>=90) 
            rvas(i,j)=r;
        else
            rvas(i,j)=-999;            
        end 
       s=vsud(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(sud,serie);
        if(sig>=90) 
            rvs(i,j)=r;
        else
            rvs(i,j)=-999;            
        end 
        
    end
end

% load coast;
% figure(1);
% clf;
% subplot(3,2,1)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,run',rvn');
% set(h,'color',[0,0,0]);
% title('Nord -1-5 dys tot')
% ylabel(tot1)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,2)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ruan',rvan');
% set(h,'color',[0,0,0]);
% title('Nord -1-5 dys ano')
% ylabel(tota1)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,3)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ruc',rvc');
% set(h,'color',[0,0,0]);
% title('Cent -1-5 dys tot')
% ylabel(tot2)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,4)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ruac',rvac');
% set(h,'color',[0,0,0]);
% title('Cent -1-5 dys ano')
% ylabel(tota2)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,5)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,rus',rvs');
% set(h,'color',[0,0,0]);
% title('Sud -1-5 dys tot')
% ylabel(tot3)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,6)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ruas',rvas');
% set(h,'color',[0,0,0]);
% title('Sud -1-5 dys ano')
% ylabel(tota3)
% axis([-30 30 0 30])
% suptitle('Spearman Rank Corr onset date vs circulation at 925hPa during 00 to +05 days after onset')
% orient tall;
% eval(sprintf('print ruv92start00-05spearm.ps -dpsc'))

fid = fopen('ruv92start00-05spearm.bin','wb');
a=fwrite(fid,run,'float32');
a=fwrite(fid,rvn,'float32');
a=fwrite(fid,ruan,'float32');
a=fwrite(fid,rvan,'float32');
a=fwrite(fid,ruc,'float32');
a=fwrite(fid,rvc,'float32');
a=fwrite(fid,ruac,'float32');
a=fwrite(fid,rvac,'float32');
a=fwrite(fid,rus,'float32');
a=fwrite(fid,rvs,'float32');
a=fwrite(fid,ruas,'float32');
a=fwrite(fid,rvas,'float32');
fclose(fid);

for i=1:nx
    for j=1:ny
        [r,sig]=signpermut(anord,uanord(i,j,:),100);
        if(sig>=90) 
            ruan(i,j)=r;
        else
            ruan(i,j)=-999.;
        end 
        [r,sig]=signpermut(nord,unord(i,j,:),100);
        if(sig>=90) 
            run(i,j)=r;
        else
            run(i,j)=-999.;
        end 
        
        [r,sig]=signpermut(acent,uacent(i,j,:),100);
        if(sig>=90) 
            ruac(i,j)=r;
        else
            ruac(i,j)=-999.;
        end 
        [r,sig]=signpermut(cent,ucent(i,j,:),100);
        if(sig>=90) 
            ruc(i,j)=r;
        else
            ruc(i,j)=-999.;
        end 
        
        [r,sig]=signpermut(asud,uasud(i,j,:),100);
        if(sig>=90) 
            ruas(i,j)=r;
        else
            ruas(i,j)=-999.;
        end 
        [r,sig]=signpermut(sud,usud(i,j,:),100);
        if(sig>=90) 
            rus(i,j)=r;
        else
            rus(i,j)=-999.;
        end

%%%%%%%%%%%%%
        [r,sig]=signpermut(anord,vanord(i,j,:),100);
        if(sig>=90) 
            rvan(i,j)=r;
        else
            rvan(i,j)=-999.;
        end 
        [r,sig]=signpermut(nord,vnord(i,j,:),100);
        if(sig>=90) 
            rvn(i,j)=r;
        else
            rvn(i,j)=-999.;
        end 
        
        [r,sig]=signpermut(acent,vacent(i,j,:),100);
        if(sig>=90) 
            rvac(i,j)=r;
        else
            rvac(i,j)=-999.;
        end 
        [r,sig]=signpermut(cent,vcent(i,j,:),100);
        if(sig>=90) 
            rvc(i,j)=r;
        else
            rvc(i,j)=-999.;
        end 
        
        [r,sig]=signpermut(asud,vasud(i,j,:),100);
        if(sig>=90) 
            rvas(i,j)=r;
        else
            rvas(i,j)=-999.;
        end 
        [r,sig]=signpermut(sud,vsud(i,j,:),100);
        if(sig>=90) 
            rvs(i,j)=r;
        else
            rvs(i,j)=-999.;
        end 
    
    end
end

% load coast;
% figure(1);
% clf;
% subplot(3,2,1)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,run',rvn');
% set(h,'color',[0,0,0]);
% title('Nord -1-5 dys tot')
% ylabel(tot1)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,2)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ruan',rvan');
% set(h,'color',[0,0,0]);
% title('Nord -1-5 dys ano')
% ylabel(tota1)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,3)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ruc',rvc');
% set(h,'color',[0,0,0]);
% title('Cent -1-5 dys tot')
% ylabel(tot2)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,4)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ruac',rvac');
% set(h,'color',[0,0,0]);
% title('Cent -1-5 dys ano')
% ylabel(tota2)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,5)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,rus',rvs');
% set(h,'color',[0,0,0]);
% title('Sud -1-5 dys tot')
% ylabel(tot3)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,6)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ruas',rvas');
% set(h,'color',[0,0,0]);
% title('Sud -1-5 dys ano')
% ylabel(tota3)
% axis([-30 30 0 30])
% suptitle('Corr onset date vs circulation at 925hPa 00 to +05 days after onset')
% orient tall;
% eval(sprintf('print ruv92start00-05corr.ps -dpsc'))

fid = fopen('ruv92start00-05corr.bin','wb');
a=fwrite(fid,run,'float32');
a=fwrite(fid,rvn,'float32');
a=fwrite(fid,ruan,'float32');
a=fwrite(fid,rvan,'float32');
a=fwrite(fid,ruc,'float32');
a=fwrite(fid,rvc,'float32');
a=fwrite(fid,ruac,'float32');
a=fwrite(fid,rvac,'float32');
a=fwrite(fid,rus,'float32');
a=fwrite(fid,rvs,'float32');
a=fwrite(fid,ruas,'float32');
a=fwrite(fid,rvas,'float32');
fclose(fid);


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%   Finding the Precip on dates of onsets t=-5 to -1
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 
tot1=0.;
tot2=0.;
tot3=0.;
tota1=0.;
tota2=0.;
tota3=0.;
for i=1:37
    tota1=tota1+sum(r1a(nord(i)-5:nord(i)-1,i));
    tot1=tot1+sum(r1(nord(i)-5:nord(i)-1,i));  
    tota2=tota2+sum(r2a(cent(i)-5:cent(i)-1,i));    
    tot2=tot2+sum(r2(cent(i)-5:cent(i)-1,i));  
    tota3=tota3+sum(r3a(sud(i)-5:sud(i)-1,i));        
    tot3=tot3+sum(r3(sud(i)-5:sud(i)-1,i));  
end

%%%%%%%%%% Find the circulation dates 
n1=32;
for i=1:37
    for ix=1:nx
        for iy=1:ny
            indan=(i-1)*184;
            n1=nord(i)-121; % 121= 01Jan-30Apr
            ijour=indan+n1;  % exact day in the matrix [37*184]==6808
            uanord(ix,iy,i)=mean(ua(iy,ix,ijour-5:ijour-1));
            vanord(ix,iy,i)=mean(va(iy,ix,ijour-5:ijour-1));
            unord(ix,iy,i)=mean(u(iy,ix,ijour-5:ijour-1));
            vnord(ix,iy,i)=mean(v(iy,ix,ijour-5:ijour-1));            

            n1=cent(i)-121; 
            ijour=indan+n1;
            uacent(ix,iy,i)=mean(ua(iy,ix,ijour-5:ijour-1));
            vacent(ix,iy,i)=mean(va(iy,ix,ijour-5:ijour-1));
            ucent(ix,iy,i)=mean(u(iy,ix,ijour-5:ijour-1));
            vcent(ix,iy,i)=mean(v(iy,ix,ijour-5:ijour-1));            

            n1=sud(i)-121;   
            ijour=indan+n1;
            uasud(ix,iy,i)=mean(ua(iy,ix,ijour-5:ijour-1));
            vasud(ix,iy,i)=mean(va(iy,ix,ijour-5:ijour-1));
            usud(ix,iy,i)=mean(u(iy,ix,ijour-5:ijour-1));
            vsud(ix,iy,i)=mean(v(iy,ix,ijour-5:ijour-1));            
        end
    end
end
 
% % averaging the circulations during those dates  
unordm=mean(unord,3);
vnordm=mean(vnord,3);
ucentm=mean(ucent,3);
vcentm=mean(vcent,3);
usudm=mean(usud,3);
vsudm=mean(vsud,3);

uanordm=mean(uanord,3);
vanordm=mean(vanord,3);
uacentm=mean(uacent,3);
vacentm=mean(vacent,3);
uasudm=mean(uasud,3);
vasudm=mean(vasud,3);

% load coast;
% figure(1);
% clf;
% subplot(3,2,1)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,unordm',vnordm');
% set(h,'color',[0,0,0]);
% title('Nord -1-5 dys tot')
% ylabel(tot1)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,2)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,uanordm',vanordm');
% set(h,'color',[0,0,0]);
% title('Nord -1-5 dys ano')
% ylabel(tota1)
% axis([-180 180 -40 50]);
% 
% 
% subplot(3,2,3)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ucentm',vcentm');
% set(h,'color',[0,0,0]);
% title('Cent -1-5 dys tot')
% ylabel(tot2)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,4)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,uacentm',vacentm');
% set(h,'color',[0,0,0]);
% title('Cent -1-5 dys ano')
% ylabel(tota2)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,5)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,usudm',vsudm');
% set(h,'color',[0,0,0]);
% title('Sud -1-5 dys tot')
% ylabel(tot3)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,6)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,uasudm',vasudm');
% set(h,'color',[0,0,0]);
% title('Sud -1-5 dys ano')
% ylabel(tota3)
% axis([-30 30 0 30])
% suptitle('Mean circulation at 925hPa during -1-5 days before onset')
% orient tall;
% eval(sprintf('print circ92start-1.ps -dpsc'))
% 
% %pause;
fid=fopen('circ92start-01-05.bin','wb');
a=fwrite(fid,unordm,'float32');
a=fwrite(fid,vnordm,'float32');
a=fwrite(fid,uanordm,'float32');
a=fwrite(fid,vanordm,'float32');
a=fwrite(fid,ucentm,'float32');
a=fwrite(fid,vcentm,'float32');
a=fwrite(fid,uacentm,'float32');
a=fwrite(fid,vacentm,'float32');
a=fwrite(fid,usudm,'float32');
a=fwrite(fid,vsudm,'float32');
a=fwrite(fid,uasudm,'float32');
a=fwrite(fid,vasudm,'float32');
fclose(fid);
%%%%%%%%%%%%%

% correlating using Spearman rank correlation  
for i=1:nx
    for j=1:ny
       s=uanord(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(anord,serie);
        if(sig>=90) 
            ruan(i,j)=r;
        else
            ruan(i,j)=-999;            
        end 
       s=unord(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(nord,serie);
        if(sig>=90) 
            run(i,j)=r;
        else
            run(i,j)=-999;            
        end 
        
       s=uacent(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(acent,serie);
        if(sig>=90) 
            ruac(i,j)=r;
        else
            ruac(i,j)=-999;            
        end 
       s=ucent(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(cent,serie);
        if(sig>=90) 
            ruc(i,j)=r;
        else
            ruc(i,j)=-999;            
        end 
        
       s=uasud(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(asud,serie);
        if(sig>=90) 
            ruas(i,j)=r;
        else
            ruas(i,j)=-999;            
        end 
       s=usud(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(sud,serie);
        if(sig>=90) 
            rus(i,j)=r;
        else
            rus(i,j)=-999;            
        end 

%%% V component
       s=vanord(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(anord,serie);
        if(sig>=90) 
            rvan(i,j)=r;
        else
            rvan(i,j)=-999;            
        end 
       s=vnord(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(nord,serie);
        if(sig>=90) 
            rvn(i,j)=r;
        else
            rvn(i,j)=-999;            
        end 
        
       s=vacent(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(acent,serie);
        if(sig>=90) 
            rvac(i,j)=r;
        else
            rvac(i,j)=-999;            
        end 
       s=vcent(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(cent,serie);
        if(sig>=90) 
            rvc(i,j)=r;
        else
            rvc(i,j)=-999;            
        end 
        
       s=vasud(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(asud,serie);
        if(sig>=90) 
            rvas(i,j)=r;
        else
            rvas(i,j)=-999;            
        end 
       s=vsud(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(sud,serie);
        if(sig>=90) 
            rvs(i,j)=r;
        else
            rvs(i,j)=-999;            
        end 
        
    end
end

% load coast;
% figure(1);
% clf;
% subplot(3,2,1)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,run',rvn');
% set(h,'color',[0,0,0]);
% title('Nord -1-5 dys tot')
% ylabel(tot1)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,2)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ruan',rvan');
% set(h,'color',[0,0,0]);
% title('Nord -1-5 dys ano')
% ylabel(tota1)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,3)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ruc',rvc');
% set(h,'color',[0,0,0]);
% title('Cent -1-5 dys tot')
% ylabel(tot2)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,4)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ruac',rvac');
% set(h,'color',[0,0,0]);
% title('Cent -1-5 dys ano')
% ylabel(tota2)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,5)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,rus',rvs');
% set(h,'color',[0,0,0]);
% title('Sud -1-5 dys tot')
% ylabel(tot3)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,6)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ruas',rvas');
% set(h,'color',[0,0,0]);
% title('Sud -1-5 dys ano')
% ylabel(tota3)
% axis([-30 30 0 30])
% suptitle('Spearman Rank Corr onset date vs circulation at 925hPa during -1-5 days before')
% orient tall;
% eval(sprintf('print ruv92start-1spear.ps -dpsc'))

fid = fopen('ruv92start-01-05spearm.bin','wb');
a=fwrite(fid,run,'float32');
a=fwrite(fid,rvn,'float32');
a=fwrite(fid,ruan,'float32');
a=fwrite(fid,rvan,'float32');
a=fwrite(fid,ruc,'float32');
a=fwrite(fid,rvc,'float32');
a=fwrite(fid,ruac,'float32');
a=fwrite(fid,rvac,'float32');
a=fwrite(fid,rus,'float32');
a=fwrite(fid,rvs,'float32');
a=fwrite(fid,ruas,'float32');
a=fwrite(fid,rvas,'float32');
fclose(fid);

for i=1:nx
    for j=1:ny
        [r,sig]=signpermut(anord,uanord(i,j,:),100);
        if(sig>=90) 
            ruan(i,j)=r;
        else
            ruan(i,j)=-999.;
        end 
        [r,sig]=signpermut(nord,unord(i,j,:),100);
        if(sig>=90) 
            run(i,j)=r;
        else
            run(i,j)=-999.;
        end 
        
        [r,sig]=signpermut(acent,uacent(i,j,:),100);
        if(sig>=90) 
            ruac(i,j)=r;
        else
            ruac(i,j)=-999.;
        end 
        [r,sig]=signpermut(cent,ucent(i,j,:),100);
        if(sig>=90) 
            ruc(i,j)=r;
        else
            ruc(i,j)=-999.;
        end 
        
        [r,sig]=signpermut(asud,uasud(i,j,:),100);
        if(sig>=90) 
            ruas(i,j)=r;
        else
            ruas(i,j)=-999.;
        end 
        [r,sig]=signpermut(sud,usud(i,j,:),100);
        if(sig>=90) 
            rus(i,j)=r;
        else
            rus(i,j)=-999.;
        end

%%%%%%%%%%%%%
        [r,sig]=signpermut(anord,vanord(i,j,:),100);
        if(sig>=90) 
            rvan(i,j)=r;
        else
            rvan(i,j)=-999.;
        end 
        [r,sig]=signpermut(nord,vnord(i,j,:),100);
        if(sig>=90) 
            rvn(i,j)=r;
        else
            rvn(i,j)=-999.;
        end 
        
        [r,sig]=signpermut(acent,vacent(i,j,:),100);
        if(sig>=90) 
            rvac(i,j)=r;
        else
            rvac(i,j)=-999.;
        end 
        [r,sig]=signpermut(cent,vcent(i,j,:),100);
        if(sig>=90) 
            rvc(i,j)=r;
        else
            rvc(i,j)=-999.;
        end 
        
        [r,sig]=signpermut(asud,vasud(i,j,:),100);
        if(sig>=90) 
            rvas(i,j)=r;
        else
            rvas(i,j)=-999.;
        end 
        [r,sig]=signpermut(sud,vsud(i,j,:),100);
        if(sig>=90) 
            rvs(i,j)=r;
        else
            rvs(i,j)=-999.;
        end 
    
    end
end

% load coast;
% figure(1);
% clf;
% subplot(3,2,1)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,run',rvn');
% set(h,'color',[0,0,0]);
% title('Nord -1-5 dys tot')
% ylabel(tot1)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,2)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ruan',rvan');
% set(h,'color',[0,0,0]);
% title('Nord -1-5 dys ano')
% ylabel(tota1)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,3)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ruc',rvc');
% set(h,'color',[0,0,0]);
% title('Cent -1-5 dys tot')
% ylabel(tot2)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,4)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ruac',rvac');
% set(h,'color',[0,0,0]);
% title('Cent -1-5 dys ano')
% ylabel(tota2)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,5)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,rus',rvs');
% set(h,'color',[0,0,0]);
% title('Sud -1-5 dys tot')
% ylabel(tot3)
% axis([-180 180 -40 50]);
% 
% subplot(3,2,6)
% h=plot(long,lat,'LineWidth',1);
% set(h,'color',[.5,.5,.5]);
% hold on
% h=quiver(x,y,ruas',rvas');
% set(h,'color',[0,0,0]);
% title('Sud -1-5 dys ano')
% ylabel(tota3)
% axis([-30 30 0 30])
% suptitle('Corr onset date vs circulation at 925hPa during -1-5 days before')
% orient tall;
% eval(sprintf('print ruv92start-1corr.ps -dpsc'))

fid = fopen('ruv92start-01-05corr.bin','wb');
a=fwrite(fid,run,'float32');
a=fwrite(fid,rvn,'float32');
a=fwrite(fid,ruan,'float32');
a=fwrite(fid,rvan,'float32');
a=fwrite(fid,ruc,'float32');
a=fwrite(fid,rvc,'float32');
a=fwrite(fid,ruac,'float32');
a=fwrite(fid,rvac,'float32');
a=fwrite(fid,rus,'float32');
a=fwrite(fid,rvs,'float32');
a=fwrite(fid,ruas,'float32');
a=fwrite(fid,rvas,'float32');
fclose(fid);


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%   Finding the Precip on dates of onsets t= -6 to -10
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 
tot1=0.;
tot2=0.;
tot3=0.;
tota1=0.;
tota2=0.;
tota3=0.;
for i=1:37
    tota1=tota1+sum(r1a(nord(i)-10:nord(i)-6,i));
    tot1=tot1+sum(r1(nord(i)-10:nord(i)-6,i));  
    tota2=tota2+sum(r2a(cent(i)-10:cent(i)-6,i));    
    tot2=tot2+sum(r2(cent(i)-10:cent(i)-6,i));  
    tota3=tota3+sum(r3a(sud(i)-10:sud(i)-6,i));        
    tot3=tot3+sum(r3(sud(i)-10:sud(i)-6,i));  
end

%%%%%%%%%% Find the circulation dates 
n1=32;
for i=1:37
    for ix=1:nx
        for iy=1:ny
            indan=(i-1)*184;
            n1=nord(i)-121; % 121= 01Jan-30Apr
            ijour=indan+n1;  % exact day in the matrix [37*184]==6808
            uanord(ix,iy,i)=mean(ua(iy,ix,ijour-10:ijour-6));
            vanord(ix,iy,i)=mean(va(iy,ix,ijour-10:ijour-6));
            unord(ix,iy,i)=mean(u(iy,ix,ijour-10:ijour-6));
            vnord(ix,iy,i)=mean(v(iy,ix,ijour-10:ijour-6));            

            n1=cent(i)-121; 
            ijour=indan+n1;
            uacent(ix,iy,i)=mean(ua(iy,ix,ijour-10:ijour-6));
            vacent(ix,iy,i)=mean(va(iy,ix,ijour-10:ijour-6));
            ucent(ix,iy,i)=mean(u(iy,ix,ijour-10:ijour-6));
            vcent(ix,iy,i)=mean(v(iy,ix,ijour-10:ijour-6));            

            n1=sud(i)-121;   
            ijour=indan+n1;
            uasud(ix,iy,i)=mean(ua(iy,ix,ijour-10:ijour-6));
            vasud(ix,iy,i)=mean(va(iy,ix,ijour-10:ijour-6));
            usud(ix,iy,i)=mean(u(iy,ix,ijour-10:ijour-6));
            vsud(ix,iy,i)=mean(v(iy,ix,ijour-10:ijour-6));            
        end
    end
end
 
% % averaging the circulations during those dates  
unordm=mean(unord,3);
vnordm=mean(vnord,3);
ucentm=mean(ucent,3);
vcentm=mean(vcent,3);
usudm=mean(usud,3);
vsudm=mean(vsud,3);

uanordm=mean(uanord,3);
vanordm=mean(vanord,3);
uacentm=mean(uacent,3);
vacentm=mean(vacent,3);
uasudm=mean(uasud,3);
vasudm=mean(vasud,3);

fid=fopen('circ92start-06-10.bin','wb');
a=fwrite(fid,unordm,'float32');
a=fwrite(fid,vnordm,'float32');
a=fwrite(fid,uanordm,'float32');
a=fwrite(fid,vanordm,'float32');
a=fwrite(fid,ucentm,'float32');
a=fwrite(fid,vcentm,'float32');
a=fwrite(fid,uacentm,'float32');
a=fwrite(fid,vacentm,'float32');
a=fwrite(fid,usudm,'float32');
a=fwrite(fid,vsudm,'float32');
a=fwrite(fid,uasudm,'float32');
a=fwrite(fid,vasudm,'float32');
fclose(fid);
%%%%%%%%%%%%%

% correlating using Spearman rank correlation  
for i=1:nx
    for j=1:ny
       s=uanord(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(anord,serie);
        if(sig>=90) 
            ruan(i,j)=r;
        else
            ruan(i,j)=-999;            
        end 
       s=unord(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(nord,serie);
        if(sig>=90) 
            run(i,j)=r;
        else
            run(i,j)=-999;            
        end 
        
       s=uacent(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(acent,serie);
        if(sig>=90) 
            ruac(i,j)=r;
        else
            ruac(i,j)=-999;            
        end 
       s=ucent(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(cent,serie);
        if(sig>=90) 
            ruc(i,j)=r;
        else
            ruc(i,j)=-999;            
        end 
        
       s=uasud(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(asud,serie);
        if(sig>=90) 
            ruas(i,j)=r;
        else
            ruas(i,j)=-999;            
        end 
       s=usud(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(sud,serie);
        if(sig>=90) 
            rus(i,j)=r;
        else
            rus(i,j)=-999;            
        end 

%%% V component
       s=vanord(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(anord,serie);
        if(sig>=90) 
            rvan(i,j)=r;
        else
            rvan(i,j)=-999;            
        end 
       s=vnord(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(nord,serie);
        if(sig>=90) 
            rvn(i,j)=r;
        else
            rvn(i,j)=-999;            
        end 
        
       s=vacent(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(acent,serie);
        if(sig>=90) 
            rvac(i,j)=r;
        else
            rvac(i,j)=-999;            
        end 
       s=vcent(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(cent,serie);
        if(sig>=90) 
            rvc(i,j)=r;
        else
            rvc(i,j)=-999;            
        end 
        
       s=vasud(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(asud,serie);
        if(sig>=90) 
            rvas(i,j)=r;
        else
            rvas(i,j)=-999;            
        end 
       s=vsud(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(sud,serie);
        if(sig>=90) 
            rvs(i,j)=r;
        else
            rvs(i,j)=-999;            
        end 
        
    end
end


fid = fopen('ruv92start-06-10spearm.bin','wb');
a=fwrite(fid,run,'float32');
a=fwrite(fid,rvn,'float32');
a=fwrite(fid,ruan,'float32');
a=fwrite(fid,rvan,'float32');
a=fwrite(fid,ruc,'float32');
a=fwrite(fid,rvc,'float32');
a=fwrite(fid,ruac,'float32');
a=fwrite(fid,rvac,'float32');
a=fwrite(fid,rus,'float32');
a=fwrite(fid,rvs,'float32');
a=fwrite(fid,ruas,'float32');
a=fwrite(fid,rvas,'float32');
fclose(fid);

for i=1:nx
    for j=1:ny
        [r,sig]=signpermut(anord,uanord(i,j,:),100);
        if(sig>=90) 
            ruan(i,j)=r;
        else
            ruan(i,j)=-999.;
        end 
        [r,sig]=signpermut(nord,unord(i,j,:),100);
        if(sig>=90) 
            run(i,j)=r;
        else
            run(i,j)=-999.;
        end 
        
        [r,sig]=signpermut(acent,uacent(i,j,:),100);
        if(sig>=90) 
            ruac(i,j)=r;
        else
            ruac(i,j)=-999.;
        end 
        [r,sig]=signpermut(cent,ucent(i,j,:),100);
        if(sig>=90) 
            ruc(i,j)=r;
        else
            ruc(i,j)=-999.;
        end 
        
        [r,sig]=signpermut(asud,uasud(i,j,:),100);
        if(sig>=90) 
            ruas(i,j)=r;
        else
            ruas(i,j)=-999.;
        end 
        [r,sig]=signpermut(sud,usud(i,j,:),100);
        if(sig>=90) 
            rus(i,j)=r;
        else
            rus(i,j)=-999.;
        end

%%%%%%%%%%%%%
        [r,sig]=signpermut(anord,vanord(i,j,:),100);
        if(sig>=90) 
            rvan(i,j)=r;
        else
            rvan(i,j)=-999.;
        end 
        [r,sig]=signpermut(nord,vnord(i,j,:),100);
        if(sig>=90) 
            rvn(i,j)=r;
        else
            rvn(i,j)=-999.;
        end 
        
        [r,sig]=signpermut(acent,vacent(i,j,:),100);
        if(sig>=90) 
            rvac(i,j)=r;
        else
            rvac(i,j)=-999.;
        end 
        [r,sig]=signpermut(cent,vcent(i,j,:),100);
        if(sig>=90) 
            rvc(i,j)=r;
        else
            rvc(i,j)=-999.;
        end 
        
        [r,sig]=signpermut(asud,vasud(i,j,:),100);
        if(sig>=90) 
            rvas(i,j)=r;
        else
            rvas(i,j)=-999.;
        end 
        [r,sig]=signpermut(sud,vsud(i,j,:),100);
        if(sig>=90) 
            rvs(i,j)=r;
        else
            rvs(i,j)=-999.;
        end 
    
    end
end


fid = fopen('ruv92start-06-10corr.bin','wb');
a=fwrite(fid,run,'float32');
a=fwrite(fid,rvn,'float32');
a=fwrite(fid,ruan,'float32');
a=fwrite(fid,rvan,'float32');
a=fwrite(fid,ruc,'float32');
a=fwrite(fid,rvc,'float32');
a=fwrite(fid,ruac,'float32');
a=fwrite(fid,rvac,'float32');
a=fwrite(fid,rus,'float32');
a=fwrite(fid,rvs,'float32');
a=fwrite(fid,ruas,'float32');
a=fwrite(fid,rvas,'float32');
fclose(fid);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%   Finding the Precip on dates of onsets t= -11 to -15
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 
tot1=0.;
tot2=0.;
tot3=0.;
tota1=0.;
tota2=0.;
tota3=0.;
for i=1:37
    tota1=tota1+sum(r1a(nord(i)-15:nord(i)-11,i));
    tot1=tot1+sum(r1(nord(i)-15:nord(i)-11,i));  
    tota2=tota2+sum(r2a(cent(i)-15:cent(i)-11,i));    
    tot2=tot2+sum(r2(cent(i)-15:cent(i)-11,i));  
    tota3=tota3+sum(r3a(sud(i)-15:sud(i)-11,i));        
    tot3=tot3+sum(r3(sud(i)-15:sud(i)-11,i));  
end

%%%%%%%%%% Find the circulation dates 
n1=32;
for i=1:37
    for ix=1:nx
        for iy=1:ny
            indan=(i-1)*184;
            n1=nord(i)-121; % 121= 01Jan-30Apr
            ijour=indan+n1;  % exact day in the matrix [37*184]==6808
            uanord(ix,iy,i)=mean(ua(iy,ix,ijour-15:ijour-11));
            vanord(ix,iy,i)=mean(va(iy,ix,ijour-15:ijour-11));
            unord(ix,iy,i)=mean(u(iy,ix,ijour-15:ijour-11));
            vnord(ix,iy,i)=mean(v(iy,ix,ijour-15:ijour-11));            

            n1=cent(i)-121; 
            ijour=indan+n1;
            uacent(ix,iy,i)=mean(ua(iy,ix,ijour-15:ijour-11));
            vacent(ix,iy,i)=mean(va(iy,ix,ijour-15:ijour-11));
            ucent(ix,iy,i)=mean(u(iy,ix,ijour-15:ijour-11));
            vcent(ix,iy,i)=mean(v(iy,ix,ijour-15:ijour-11));            

            n1=sud(i)-121;   
            ijour=indan+n1;
            uasud(ix,iy,i)=mean(ua(iy,ix,ijour-15:ijour-11));
            vasud(ix,iy,i)=mean(va(iy,ix,ijour-15:ijour-11));
            usud(ix,iy,i)=mean(u(iy,ix,ijour-15:ijour-11));
            vsud(ix,iy,i)=mean(v(iy,ix,ijour-15:ijour-11));            
        end
    end
end
 
% % averaging the circulations during those dates  
unordm=mean(unord,3);
vnordm=mean(vnord,3);
ucentm=mean(ucent,3);
vcentm=mean(vcent,3);
usudm=mean(usud,3);
vsudm=mean(vsud,3);

uanordm=mean(uanord,3);
vanordm=mean(vanord,3);
uacentm=mean(uacent,3);
vacentm=mean(vacent,3);
uasudm=mean(uasud,3);
vasudm=mean(vasud,3);

fid=fopen('circ92start-11-15.bin','wb');
a=fwrite(fid,unordm,'float32');
a=fwrite(fid,vnordm,'float32');
a=fwrite(fid,uanordm,'float32');
a=fwrite(fid,vanordm,'float32');
a=fwrite(fid,ucentm,'float32');
a=fwrite(fid,vcentm,'float32');
a=fwrite(fid,uacentm,'float32');
a=fwrite(fid,vacentm,'float32');
a=fwrite(fid,usudm,'float32');
a=fwrite(fid,vsudm,'float32');
a=fwrite(fid,uasudm,'float32');
a=fwrite(fid,vasudm,'float32');
fclose(fid);
%%%%%%%%%%%%%

% correlating using Spearman rank correlation  
for i=1:nx
    for j=1:ny
       s=uanord(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(anord,serie);
        if(sig>=90) 
            ruan(i,j)=r;
        else
            ruan(i,j)=-999;            
        end 
       s=unord(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(nord,serie);
        if(sig>=90) 
            run(i,j)=r;
        else
            run(i,j)=-999;            
        end 
        
       s=uacent(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(acent,serie);
        if(sig>=90) 
            ruac(i,j)=r;
        else
            ruac(i,j)=-999;            
        end 
       s=ucent(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(cent,serie);
        if(sig>=90) 
            ruc(i,j)=r;
        else
            ruc(i,j)=-999;            
        end 
        
       s=uasud(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(asud,serie);
        if(sig>=90) 
            ruas(i,j)=r;
        else
            ruas(i,j)=-999;            
        end 
       s=usud(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(sud,serie);
        if(sig>=90) 
            rus(i,j)=r;
        else
            rus(i,j)=-999;            
        end 

%%% V component
       s=vanord(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(anord,serie);
        if(sig>=90) 
            rvan(i,j)=r;
        else
            rvan(i,j)=-999;            
        end 
       s=vnord(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(nord,serie);
        if(sig>=90) 
            rvn(i,j)=r;
        else
            rvn(i,j)=-999;            
        end 
        
       s=vacent(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(acent,serie);
        if(sig>=90) 
            rvac(i,j)=r;
        else
            rvac(i,j)=-999;            
        end 
       s=vcent(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(cent,serie);
        if(sig>=90) 
            rvc(i,j)=r;
        else
            rvc(i,j)=-999;            
        end 
        
       s=vasud(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(asud,serie);
        if(sig>=90) 
            rvas(i,j)=r;
        else
            rvas(i,j)=-999;            
        end 
       s=vsud(i,j,:);serie=reshape(s,37,1);
       [r,sig]=spear2(sud,serie);
        if(sig>=90) 
            rvs(i,j)=r;
        else
            rvs(i,j)=-999;            
        end 
        
    end
end


fid = fopen('ruv92start-11-15spearm.bin','wb');
a=fwrite(fid,run,'float32');
a=fwrite(fid,rvn,'float32');
a=fwrite(fid,ruan,'float32');
a=fwrite(fid,rvan,'float32');
a=fwrite(fid,ruc,'float32');
a=fwrite(fid,rvc,'float32');
a=fwrite(fid,ruac,'float32');
a=fwrite(fid,rvac,'float32');
a=fwrite(fid,rus,'float32');
a=fwrite(fid,rvs,'float32');
a=fwrite(fid,ruas,'float32');
a=fwrite(fid,rvas,'float32');
fclose(fid);

for i=1:nx
    for j=1:ny
        [r,sig]=signpermut(anord,uanord(i,j,:),100);
        if(sig>=90) 
            ruan(i,j)=r;
        else
            ruan(i,j)=-999.;
        end 
        [r,sig]=signpermut(nord,unord(i,j,:),100);
        if(sig>=90) 
            run(i,j)=r;
        else
            run(i,j)=-999.;
        end 
        
        [r,sig]=signpermut(acent,uacent(i,j,:),100);
        if(sig>=90) 
            ruac(i,j)=r;
        else
            ruac(i,j)=-999.;
        end 
        [r,sig]=signpermut(cent,ucent(i,j,:),100);
        if(sig>=90) 
            ruc(i,j)=r;
        else
            ruc(i,j)=-999.;
        end 
        
        [r,sig]=signpermut(asud,uasud(i,j,:),100);
        if(sig>=90) 
            ruas(i,j)=r;
        else
            ruas(i,j)=-999.;
        end 
        [r,sig]=signpermut(sud,usud(i,j,:),100);
        if(sig>=90) 
            rus(i,j)=r;
        else
            rus(i,j)=-999.;
        end

%%%%%%%%%%%%%
        [r,sig]=signpermut(anord,vanord(i,j,:),100);
        if(sig>=90) 
            rvan(i,j)=r;
        else
            rvan(i,j)=-999.;
        end 
        [r,sig]=signpermut(nord,vnord(i,j,:),100);
        if(sig>=90) 
            rvn(i,j)=r;
        else
            rvn(i,j)=-999.;
        end 
        
        [r,sig]=signpermut(acent,vacent(i,j,:),100);
        if(sig>=90) 
            rvac(i,j)=r;
        else
            rvac(i,j)=-999.;
        end 
        [r,sig]=signpermut(cent,vcent(i,j,:),100);
        if(sig>=90) 
            rvc(i,j)=r;
        else
            rvc(i,j)=-999.;
        end 
        
        [r,sig]=signpermut(asud,vasud(i,j,:),100);
        if(sig>=90) 
            rvas(i,j)=r;
        else
            rvas(i,j)=-999.;
        end 
        [r,sig]=signpermut(sud,vsud(i,j,:),100);
        if(sig>=90) 
            rvs(i,j)=r;
        else
            rvs(i,j)=-999.;
        end 
    
    end
end


fid = fopen('ruv92start-11-15corr.bin','wb');
a=fwrite(fid,run,'float32');
a=fwrite(fid,rvn,'float32');
a=fwrite(fid,ruan,'float32');
a=fwrite(fid,rvan,'float32');
a=fwrite(fid,ruc,'float32');
a=fwrite(fid,rvc,'float32');
a=fwrite(fid,ruac,'float32');
a=fwrite(fid,rvac,'float32');
a=fwrite(fid,rus,'float32');
a=fwrite(fid,rvs,'float32');
a=fwrite(fid,ruas,'float32');
a=fwrite(fid,rvas,'float32');
fclose(fid);


quit;
