%SSTa for jan
lon=[-350 -10]; lat=[-65 65];
time=[datenum(1950,1,15) datenum(2008,12, 15)];
forcd=rnc_Extract_NOAA_SST(lon, lat, time, 0);
%pcolor(forcd.lon, forcd.lat, forcd.mask); shading flat; colorbar
[year,month, day]=dates_datenum(forcd.datenum);
index=find(year == 1997 & month == 12);
SSTmean=mean(forcd.SST,3);
[SSTa, SSTm] =RemoveSeas(forcd.SST, month, 3);
SSTseas=SSTm- repmat(SSTmean, [1 1 12]);
%mfig;
%clf;
pcolor(forcd.lon, forcd.lat, SSTseas(:,:,1).*forcd.mask); shading flat; colorbar;
hold on;
world_coast
title 'SST mean for Jan'
caxis([-5 5]);
colorbar;
shading interp
%print -dpng  dec.png


SST=forcd.SST(:,:,index);
SSTmean=mean(forcd.SST,3);
% Removes just the annual mean but retain the seasonal component
SSTa=SST-SSTmean;
mfig;
clf;
pcolor(forcd.lon, forcd.lat, SSTa.*forcd.mask); shading flat; colorbar;
hold on;
world_coast


% Remove mean and seasonal components
[SSTa, SSTm] =RemoveSeas(forcd.SST, month, 3);
SSTseas=SSTm- repmat(SSTmean, [1 1 12]);
mfig;
clf;
pcolor(forcd.lon, forcd.lat, SSTa(:,:,index).*forcd.mask); shading flat; colorbar;
hold on;
world_coast
title 'SSTa anomaly for Dec 1997'
mfig;
clf;
pcolor(forcd.lon, forcd.lat, SSTseas(:,:,1).*forcd.mask); shading flat; colorbar;
hold on;
world_coast
title 'SST mean for Jan'


% Make 12 figures of the average seasonal cycle (SSTseas)
% Look at the temproal evolution of the anomalies during the period 
% Feb 1997 to Feb of 1998 to Feb 1999 (El Nino to La Nina) (SSTa)
% Make a correlation map between the SSTa and the NPGO index (SSTa and
% NPGO)
[npgo_time, npgo] = rout_npgo;
pdo=rout_pdo;



% Feb 1999 (La Nina) (SSTa)
lon=[-350 -10]; lat=[-65 65];
time=[datenum(1950,1,15) datenum(2008,12, 15)];
forcd=rnc_Extract_NOAA_SST(lon, lat, time, 0);

pcolor(forcd.lon, forcd.lat, forcd.mask); shading flat; colorbar

[year,month, day]=dates_datenum(forcd.datenum);
index=find(year == 1999 & month == 3);

SST=forcd.SST(:,:,index);
SSTmean=mean(forcd.SST,3);

SSTa=SST-SSTmean;

pcolor(forcd.lon, forcd.lat, SSTa.*forcd.mask); shading flat; colorbar;
hold on;
world_coast
% Remove mean and seasonal components
[SSTa, SSTm] =RemoveSeas(forcd.SST, month, 3);
SSTseas=SSTm- repmat(SSTmean, [1 1 12]);
mfig;
pcolor(forcd.lon, forcd.lat, SSTa(:,:,index).*forcd.mask); shading flat; colorbar;
hold on;
world_coast
title 'SSTa anomaly for Feb 1999 (mean and seasonal removed)'


% Make a correlation map between the SSTa and the NPGO index (SSTa and
% NPGO)
lon=[-350 -10]; lat=[-65 65];
time=[datenum(1950,1,15) datenum(2008,12, 15)];
forcd=rnc_Extract_NOAA_SST(lon, lat, time, 0);
[year,month, day]=dates_datenum(forcd.datenum);
SSTmean=mean(forcd.SST,3);
[SSTa, SSTm] =RemoveSeas(forcd.SST, month, 3);
SSTseas=SSTm- repmat(SSTmean, [1 1 12]);
in=find(year == 1993 & month == 1);
in2=find(year == 2004 & month == 12);
SSTa=SSTa(:,:,in:in2);

[npgo_time, npgo] = rout_npgo;
[year_npgo,month_npgo, day_npgo]=dates_datenum(npgo_time);
year_npgo=year_npgo(1:708);
month_npgo=month_npgo(1:708);
day_npgo=day_npgo(1:708);
index=find(year_npgo == 1993 & month_npgo == 1);
index2=find(year_npgo == 2004 & month_npgo == 12);
npgo=npgo(index:index2);
npgo_time=npgo_time(index:index2);

SSTa_mean=zeros(173,69);
i=1;
s=[];
for i=1:173
    for j=1:69
        for k=1:144
            s(k)=SSTa(i,j,k);
            k=k+1;
        end
        q = corrcoef(s,npgo);
        SSTa_mean(i,j) = q(1,2);
        j=j+1;
        
    end
    i=i+1;
end
SSTa_mean = imrotate(SSTa_mean,90);
pcolor(SSTa_mean);
shading flat;
colorbar;
hold on;
world_coast