Global Index (short | long) | Local contents | Local Index (short | long)
Get climatology for October through March
| This script calls | |
|---|---|
cd /home/disk/hayes/dvimont/ccm3.6/data
filin = 'T42M5079.nc'
nc = netcdf(filin, 'nowrite');
sstcyc = nc{'SST'}(:);
lat = nc{'lat'}(:);
lon = nc{'lon'}(:);
nc = close(nc)
sstcyc = sstcyc([9:12 1:4],:,:);
[ntim, nlat, nlon] = size(sstcyc);
% Get CT and GR patterns
cd /home/disk/hayes2/dvimont/ccm/ccm3.6/data
load ct_gr.mat
whos
% Normalize ctstar
ctstar = (ctstar-mean(ctstar))/std(ctstar);
% Try some plots...
subplot(2,1,1)
plot(ctstar)
set(gca, 'XTick',1:12:409,'XTickLabel',[1:12:409])
axis([0 408 -3 3.5])
grid
subplot(2,1,2)
plot(ctstar)
set(gca, 'XTick',1:12:409,'XTickLabel',[61:1:100])
axis([0 408 -3 3.5])
grid
% Get average cycle of CTstar and such...
cyr = [49 85 121 157 181 277 289 337 ];
wyr = [37 61 97 109 145 193 265 325 373]
cave = [];
wave = [];
for i = 1:ntim
cave = [cave ctstar(cyr - 5 + i)];
wave = [wave ctstar(wyr - 5 + i)];
end
% Compare warm to cold years...
figure(1)
subplot(1,2,1)
plot(1:ntim,mean(wave),'-ok',1:ntim,wave,'-k')
axis([0 ntim+1 0 3.5])
grid
set(gca, 'XTick', 1:2:7, 'XTickLabel',['SEP';'NOV';'JAN';'MAR']);
subplot(1,2,2)
plot(1:ntim,mean(cave),'-xk',1:ntim,cave,'-k')
axis([0 ntim+1 -3.5 0])
grid
set(gca, 'XTick', 1:2:7, 'XTickLabel',['SEP';'NOV';'JAN';'MAR']);
% Get average cycle of warm and cold years
moave = mean([mean(wave); -1*mean(cave)]);
% Compare warm to cold years:
subplot(1,2,1)
plot(1:ntim,mean(-1*cave),'-ok',1:ntim,mean(wave),'-xk')
axis([0 ntim+1 0 2.5])
set(gca, 'XTick', 1:2:7, 'XTickLabel',['SEP';'NOV';'JAN';'MAR']);
grid
subplot(1,2,2)
plot(1:ntim,moave,'-xk')
axis([0 ntim+1 0 2.5])
set(gca, 'XTick', 1:2:7, 'XTickLabel',['SEP';'NOV';'JAN';'MAR']);
grid
% Add CTpat to climatology to get idealized pattern to force CCM with:
sstcyc = reshape(sstcyc, ntim, nlat*nlon);
ctpat = reshape(ctpat, 1, nlat*nlon);
iceflag = min(min(min(sstcyc)));
sstnew = sstcyc;
for i = 1:ntim
noice = find(sstcyc(i,:) ~= iceflag);
sstnew(i,noice) = sstcyc(i,noice) + moave(i) * ctpat(1, noice);
end
sstnew = reshape(sstnew, ntim, nlat, nlon);
%sstcyc = reshape(sstcyc, ntim, nlat, nlon);
% Try some more plots...
lims = [0 360 -90 90]
num = 3;
figure(1);
gcont(squeeze(sstnew(num,:,:)), lims, [-2:2:38])
dc
figure(2);
gcont(squeeze(sstnew(num,:,:)-sstcyc(num,:,:)), lims, [-5:.25:5])
dc
% Write away (ct_warm_sst.nc)
%
% Note: The following files were created using the following command:
%
% ncrcat -d time,190,197 sst_nov1981-jul1998.nc sst_temp.nc
%
% First, get time axis
cd /home/disk/hayes/dvimont/ccm3.6/data
filin = 'ct_warm_sst.nc'
nc = netcdf(filin, 'nowrite');
timeunits = nc{'time'}.units(:);
nc = close(nc);
timeunits = ['days since 1981-11-14 00:00:00']
nc = netcdf(filin, 'write');
nc{'SST'}(:) = sstnew;
nc = close(nc)
ccmid = ncmex('OPEN', filin, 'NC_WRITE');
[datatype, len, status] = ncmex('ATTINQ', ccmid, 'time', 'units');
status = ncmex('ATTPUT', ccmid, 'time', 'units', datatype, len, timeunits)
ccmid = ncmex('CLOSE', ccmid);
% Write ct_cold_sst.nc
sstnew = sstcyc;
for i = 1:ntim
noice = find(sstcyc(i,:) ~= iceflag);
sstnew(i,noice) = sstcyc(i,noice) - moave(i) * ctpat(1, noice);
end
sstnew = reshape(sstnew, ntim, nlat, nlon);
filin = 'ct_cold_sst.nc'
nc = netcdf(filin, 'write');
nc{'SST'}(:) = sstnew;
nc = close(nc)
ccmid = ncmex('OPEN', filin, 'NC_WRITE');
[datatype, len, status] = ncmex('ATTINQ', ccmid, 'time', 'units');
status = ncmex('ATTPUT', ccmid, 'time', 'units', datatype, len, timeunits)
ccmid = ncmex('CLOSE', ccmid);
% Do the same for the GR patterns:
% Plot GR...
figure(1)
subplot(2,1,1)
plot(myrunning_ave(myrunning_ave(-1*gr,9),5));
set(gca, 'XTick',1:12:409,'XTickLabel',[1:12:409],'YTick',-3:3)
axis([0 408 -3 3])
grid
subplot(2,1,2)
plot(myrunning_ave(myrunning_ave(-1*gr,9),5));
set(gca, 'XTick',1:12:409,'XTickLabel',[61:1:100])
axis([0 408 -3 3])
grid
% Assume no seasonal variability in magnitude of the GR pattern. So, just
% add 1 STD of the GR pattern to the climatology.
grpat = reshape(grpat, 1, nlat*nlon);
sstnew = sstcyc;
for i = 1:ntim
noice = find(sstcyc(i,:) ~= iceflag);
sstnew(i,noice) = sstcyc(i,noice) - grpat(1, noice);
end
sstnew = reshape(sstnew, ntim, nlat, nlon);
sstnew(find(sstnew <= iceflag))=iceflag*ones(size(find(sstnew <= iceflag)));
filin = 'gr_warm_sst.nc'
nc = netcdf(filin, 'write');
nc{'SST'}(:) = sstnew;
nc = close(nc)
ccmid = ncmex('OPEN', filin, 'NC_WRITE');
[datatype, len, status] = ncmex('ATTINQ', ccmid, 'time', 'units');
status = ncmex('ATTPUT', ccmid, 'time', 'units', datatype, len, timeunits)
ccmid = ncmex('CLOSE', ccmid);
% Do the same for a cold phase:
sstnew = sstcyc;
for i = 1:ntim
noice = find(sstcyc(i,:) ~= iceflag);
sstnew(i,noice) = sstcyc(i,noice) + grpat(1, noice);
end
sstnew = reshape(sstnew, ntim, nlat, nlon);
sstnew(find(sstnew <= iceflag))=iceflag*ones(size(find(sstnew <= iceflag)));
filin = 'gr_cold_sst.nc'
nc = netcdf(filin, 'write');
nc{'SST'}(:) = sstnew;
nc = close(nc)
ccmid = ncmex('OPEN', filin, 'NC_WRITE');
[datatype, len, status] = ncmex('ATTINQ', ccmid, 'time', 'units');
status = ncmex('ATTPUT', ccmid, 'time', 'units', datatype, len, timeunits)
ccmid = ncmex('CLOSE', ccmid);