clear
% Extract SLP and SST in Pacific tropical region
trange=[datenum(1966, 1,1) datenum(1999, 11, 15)];
pac.xlim=[-250  -100]; pac.ylim=[20 62];
% Load SLP
tmp=rnc_Extract_NCEP_SLP(pac.xlim,pac.ylim,trange,0);
SLP = rnc_forcdStats(tmp,'SLP');
% Loap SST
tmp=rnc_Extract_NOAA_SST(pac.xlim,pac.ylim,trange,0);
SST = rnc_forcdStats(tmp,'SST');

% Load Zooplankton
zoo=koe_zooplantlon

% Correlation with NPGO
[tt,npgo,pdo]=rout_npgo;
[year,month]=dates_datenum(tt);
in=find(month==3);
red_signi(zoo.year, zoo.sbww, year(in)+2.5, npgo(in), 0.4, 1500, 1,'r',[30 134 69]/255,0);
red_signi(year(in)+2.5, -npgo(in), zoo.year, zoo.sbww, 0.4, 1500, 1,'k','r',0);
subplot(2,1,1); set(gca, 'ylim',[-3 3]);

% try integration model
s=AR1_model(tt, npgo, 1/12);
zz.time=datenum(zoo.year, 3, 15);
mfig

for it=0.1:0.2:2.2
red_signi(zz.time, zoo.sbww, s.time+it*360, s.sig, 30, 500, 1,'r',[30 134 69]/255,0);
it
pause
end

% KOE modes Taguchi
xlim=[143 180]; ylim=[30 45]
[ktime, kssha, ktrans, ktrans_lon, ksshap, koe2]=OFES_KOE_ssha(xlim, ylim);
ssha2=mean(koe2.ssha,1);

%ssha=dtrend2d(ssha);
[eof,pc,vexp]=rnt_doEof(ssha2, ssha2(:,:,1)*0+1);


figure(1);
p=[2.3472      2.6673      3.7876      5.6547];
set(gcf, 'paperposition',p);
red_signi(year(in)+2.5, -npgo(in), zoo.year, zoo.sbww, 0.4, 1500, 1,'k','r',0);
subplot(2,1,1); set(gca, 'ylim',[-3 4]);


clf
plot(zoo.year, detrend(zoo.sbww), 'linewidth',1.5, 'color', 'r');
hold on
plot(year(in)+2.5, -npgo(in), 'linewidth',1.5, 'color', 'b');


zoo=koe_zooplantlon
p=interp1(datenum(zoo.year, 7,15), zoo.sbww, SLP.datenum);
p(isnan(p))=0;

o=rnc_corr3(SLP.datenum, SLP.ano, SLP.mask, SLP.datenum, p);

o2=rnc_corr3(SST.datenum, SST.ano, SST.mask, SLP.datenum, p);

clf; rnc_map(o.corr, SLP);

clf; rnc_map(o2.corr, SST);


for imon=0:2:36
    o2=rnc_corr3(SST.datenum+30*imon, SST.ano, SST.mask, SLP.datenum, p);
    clf; rnc_map(o2.corr, SST); caxis([-.3 .3]); colorbar 'h'
    imon
    pause(1)
    
end

imon=12
o2=rnc_corr3(SST.datenum+30*imon, SST.ano, SST.mask, SLP.datenum, p);
clf; rnc_map(o2.corr, SST); caxis([-.3 .3]); colorbar 'h'
imon

for imon=0:2:36
    o=rnc_corr3(SLP.datenum+30*imon, SLP.ano, SLP.mask, SLP.datenum, p);
    clf; rnc_map(o.corr, SLP); caxis([-.3 .3]); colorbar 'h'
    imon
    pause(1)
    
end





for imon=1:12
    in=find(month==imon);
    red_signi(zoo.year, zoo.sbww, year(in)+2.5, npgo(in), 0.5, 1500, 1,'r','b',1);
    imon
    pause
end


in=find(month==3);
for lag=0:0.2:3
    red_signi(zoo.year, zoo.sbww, year(in)+2.5, npgo(in), 0.4, 1500, 1,'r','b',0);
    lag
    pause
end

p=interp1(datenum(zoo.year, 7,15), zoo.sbww, SLP.datenum);
p(isnan(p))=0;

p2=interp1(datenum(zoo.year, 7,15), zoo.sbww, SLP.datenum);
p2(isnan(p))=0;

%% OFES SSHa analysis
%/ocean/drive/edl/NEPD/process_data/SalinityMODE/KOE
addpath /drive/edl/NEPD/OFES
xlim=[143 180]; ylim=[32 34];
[ktime, kssha32]=OFES_KOE_ssha(xlim, ylim);
[year,month]=dates_datenum(ktime);
k=0;
for yr=1966:1999
    in=find(month>=3 & month <=6 & year==yr);
    k=k+1;
    koei(k)=mean(kssha32(in));
end
red_signi(zoo.year, zoo.sbww, 1966:1999, koei, 0.2, 1500, 1,'r','b',0);
red_signi(year(in)+2.5, -npgo(in), 1966:1999, koei, 0.2, 1500, 1,'r','b',0);

figure(1);
p=[2.3472      2.6673      3.7876      5.6547];
set(gcf, 'paperposition',p);
red_signi(1966:1999, -koei, zoo.year, zoo.sbww, 0.2, 1500, 1,'b','r',0);
subplot(2,1,1); set(gca, 'ylim',[-3 4]);


clf
plot(ktime, detrend(kssha32),'color','k','linewidth',2);
hold on; rnt_font('Helvetica',16);
datetick; set(gca,'xlim',[datenum(1950,1,1) datenum(2008,1,1)]);
[ntime, npgo]=rout_npgo;
plot(ntime+365*2.5, npgo, 'color' ,'b', 'linewidth',2);

zoo=koe_zooplantlon
p=interp1(datenum(zoo.year, 7,15), zoo.sbww, SLP.datenum);
p(isnan(p))=0;
red_signi(SLP.datenum, p, ktime, kssha32, 30, 1500, 1,'r','b',1);

[year,month]=dates_datenum(ktime);
for imon=1:12
    in=find(month==imon);
    red_signi(zoo.year, zoo.sbww, year(in), pc(in,2), 0.2, 1500, 1,'r','b',0);
    imon
    pause
end




% KOE modes Taguchi
xlim=[143 180]; ylim=[30 45]
[ktime, kssha, ktrans, ktrans_lon, ksshap, koe]=OFES_KOE_ssha(xlim, ylim);
ssha2=mean(koe.ssha,1);

%ssha=dtrend2d(ssha);
[eof,pc,vexp]=rnt_doEof(ssha2, ssha2(:,:,1)*0+1);
red_signi(ktime, pc(:,2), ktime, kssha32, 30, 1500, 1,'r','b',1);