Global Index (short | long) | Local contents | Local Index (short | long)
Start with the CT index
| This script calls | This script is called by |
|---|---|
clear
% Load the CT index
cd /home/disk/tao/dvimont/Data/Indices
load -ascii CT.asc
% The first column is the year, so retain only 1948--1998
year = CT(:,1);
keepyear = find(ismember(year, 1948:1998));
year = year(keepyear);
ct = CT(keepyear,2:13)/100;
% Reshape the CT index (now nyr X nmo) to a vector
ct = repmat(NaN, [length(year)*12, 1]);
for i = 1:length(year);
ind = 12*(i-1)+[1:12];
ct(ind) = CT(keepyear(i),2:13);
end
ct = ct/100;
% Load SST data for regressions. We'll use NCEP Reanalysis
cd /home/disk/tao/data/nmc.reanalysis/monthly
filename = 'air.mon.mean.sfc.48mar01.nc';
varname = 'air';
xylims = [0 360 -90 90]; % xylims = [minlon maxlon minlat maxlat];
level = 1; % This is arbitrary for a 3-D (but not for a 4-D) field
timelims = [1:612]; % Monthly data, 1948:1998
[sst, lat, lon] = getnc2(filename, varname, xylims, level, timelims);
% Start by removing the annual cycle from sst and ct index
[sst, clim] = annave(sst);
[ct, ctclim] = annave(ct);
% Next, reshape the data, so we can generate covariance maps
% Note: we can reshape the data back to 3-D at any time by
% the command: sst = reshape(sst, ntim, nlat, nlon);
[ntim, nlat, nlon] = size(sst);
sst = reshape(sst, ntim, nlat*nlon);
% Get regression map of sst onto ct index
regmap = ct'*sst/(var(ct)*ntim-1);
regmap = reshape(regmap, nlat, nlon);
% Plot the data:
% First, set up map window and such
figure(1); fl(1); clf;
get_global; global_figs(9, 6, 0, 0, 1.5);
dg(lat, lon);
sp(1,1); cla
maxes('mollweid', [0 180]);
% Next, contour data using only solid, black contours
[c, h] = contorm(YAX2, XAX2, regmap, [-2:.1:2], 'k');
% Block out land:
dcmfill(1, 0.7);
gridm on;
axis_limits
% Label the max and mins manually
cl = clabelm(c, 'manual');
sp(1,1); cla
maxes('mollweid', [0 180]);
% Next, contour data using only solid for positive, dashed for negative
[c, h] = contorm(YAX2, XAX2, regmap, [.1:.1:2], '-k');
hold on;
[c2, h2] = contorm(YAX2, XAX2, regmap, [-2:.1:-.1], '--k');
hold off;
% Block out land:
dcmfill(1, 0.7);
gridm on;
axis_limits
% Mess around with shading
color_trans; colormap(cmap);
sp(1,1); cla;
maxes('mollweid', [0 200]);
surfacem(YAX2, XAX2, regmap);
caxis([-1 1]);
dcmfill(1);
gridm on;
axis_limits;
regmap = ct'*sst/(std(ct)*(ntim-1));
regmap = reshape(regmap, nlat, nlon);
sp(1,1); cla;
maxes('mollweid', [0 200]);
surfacem(YAX2, XAX2, regmap);
caxis([-1 1]);
dcmfill(1);
gridm on;
axis_limits;
regmap = standardize(ct)'*standardize(sst)./(ntim-1);
regmap = reshape(regmap, nlat, nlon);
FRAME = [0 360 -90 90];
sp(1,1); cla;
maxes('mercator', [0 200]);
surfacem(YAX3, XAX3, [regmap regmap(:,1)]);
caxis([-1 1]);
dcmfill(1);
gridm on;
axis_limits;
for i = 1:3; cmap2(:,i) = interp(cmap(:,i), 5); end;
cmap2(cmap2>1) = 1;