Global Index (short | long) | Local contents | Local Index (short | long)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
The following is for previously averaged data
23 July, 1998
| This script calls | |
|---|---|
clear
cd /home/disk/hayes2/dvimont/csm/300yr/srf
load sst1-150.mat
load oro1-150.mat
load lims1-150.mat
whos
lims = [0 360 -90 90];
sst = subset(sst, nxa, nya, lims);
oro = subset(oro, nxa, nya, lims);
sst = sst(1:1800,:,:);
oro = oro(1:1800,:,:);
nxkeep = nxa(find(nxa >= lims(1) & nxa <= lims(2)));
nykeep = nya(find(nya >= lims(3) & nya <= lims(4)));
[ntmstp, nypts, nxpts] = size(sst);
whos
sst = reshape(sst, ntmstp, nypts*nxpts);
oro = reshape(oro, ntmstp, nypts*nxpts);
oro = mean(oro);
dokeep = find(oro == 0);
sst = sst(:,dokeep);
%
% Alter sst for ice points. Notice, the annual cycle feature was
% not included for figures 1 and 2 in ~dvimont/pacific
%
%noelim = [];
%for i = 1:length(dokeep);
% noelimtmp = [];
% for j = 1:12
% icept = (find(sst(j:12:ntmstp,i) <= 273.14)-1).*12 + j;
% noice = (find(sst(j:12:ntmstp,i) >= 273.14)-1).*12 + j;
% if (~isempty(icept) & ~isempty(noice));
% sst(icept,i) = mean(sst(noice,i)) * ones(length(icept),1);
% end
% if ~isempty(noice);
% noelimtmp = [noelimtmp j];
% else
% sst(icept,i) = mean(sst(icept,i)) * ones(length(icept),1);
% end
% end
% if ~isempty(noelimtmp);
% noelim = [noelim i];
% end
%end
%sst = sst(:,noelim);
%dokeep = dokeep(noelim);
% Remove any point that ever had ice
noelim = [];
for i = 1:length(dokeep);
icept = find(sst(1:ntmstp,i) <= 273.14);
if isempty(icept);
noelim = [noelim i];
end
end
sst = sst(:,noelim);
dokeep = dokeep(noelim);
whos
% Remove climatology
[sst, clim] = annave(sst);
sstemp = zeros(ntmstp, nxpts*nypts);
sstemp(:,dokeep) = sst;
sst = sstemp;
clear sstemp
sst = cosweight(sst, nya);
sst = sst(:, dokeep);
% Running Mean
sst3yr = myrunning_ave(sst, 36);
% This might make things faster
sst3yr = sst3yr(1:3:ntmstp,:);
% Take EOF's. (Might alter the dimensions if need be...)
c = sst3yr * sst3yr' / (ntmstp);
[lam, eof, per] = eof(c);
pc1 = eof(:,1:10);
eof1 = pc1'*sst3yr;
eof1 = eof1' * diag(1 ./ sqrt(diag(eof1*eof1'/length(dokeep))));
eof1 = eof1';
rescale = std(pc1);
pc1 = pc1 * diag(1./rescale);
eof1 = diag(rescale) * eof1;
temp = zeros(10, nypts*nxpts);
temp(1:10,dokeep) = eof1(1:10,:);
%
% The following command can be used to get a general idea
% of what the loadings and pc's look like. For graphics,
% refer to the file in graphics/pac_plot.m
%
%contour(squeeze(reshape(temp,1,m,n)));
%figure(1)
%subplot(2,1,1)
%gcont(squeeze(reshape(temp(1,:),1,nypts,nxpts)),lims,10)
%dc
for i = 1:length(dokeep);
a(i) = 0.5*(corr(sst3yr(:,i),sst3yr(:,i),1) + sqrt(abs(corr(sst3yr(:,i),sst3yr(:,i),2))));
end
a = mean(a);
dof = -1 * 600 * log(a) / 2;
dlam = lam * sqrt(2/dof);