Global Index (short | long) | Local contents | Local Index (short | long)
cd /home/disk/hayes2/dvimont/ocean/coup/stoc
| This script calls | |
|---|---|
clear
global RADUS
ctlim = [80 160 -6 6];
cd /home/disk/hayes2/dvimont/ocean/coup/coup
nc = netcdf('sal.cdf', 'nowrite');
lon = nc{'x'}(:);
lat = nc{'y'}(:);
lat = 1000*lat*180/(pi*RADUS);
lon = 1000*lon*180/(pi*RADUS);
[xk, yk] = keep_var(ctlim, lon, lat);
cth = nc{'h',1}(:,:,yk,xk);
cth = squeeze(cth);
nc = close(nc);
cth = squeeze(mean(mean(shiftdim(cth, 1))))-15000;
figure(1); fo(1); clf
subplot(4,1,1);
hh = plot([1:1200]./12, -rave(cth, 3), 'k', ...
[1:1200]./12, -rave(rave(cth, 25), 37), '-k');
set(hh(2), 'linewidth', 2);
set(gca, 'XTick', [1:120:1200]/12, 'XTickLabel', [0:10:100]);
axis([0 100 -5250 5250]);
set(gca, 'YTick', -5000:2500:5000)
grid on;
subplot(4,2,3);
tem = cth(1:600);
nfft = 256; nover = nfft/2; n2 = nfft/2;
[p,f] = spectrum(tem, nfft, nover);
rho = (corr(tem, tem, 1)+sqrt(corr(tem, tem, 2))+...
corr(tem, tem, 3)^(1/3))/3;
rn = (1-rho^2) ./ (1-2*rho*cos((pi/n2)*(0:n2))+rho^2);
rn = rn / mean(rn(3:(n2+1)));
rn = rn * mean(p(3:(n2+1),1));
% rner1 = rn * 1.7; rner5 = rn * 2.5;
rner1 = rn * 2.5; rner5 = rn * 1.92;
dofx = length(tem) / n2;
xt = [128 64 32 16 8 4 2 1 0.5 0.25];
f2 = f.*6;
hh = plot(log(f2), f2.*p(:,1), '-k', ...
log(f2), f2.*rn', '-k', ...
log(f2), f2.*rner5', '--k', ...
log(f2), f2.*rner1', '-.k');
set(hh(1), 'linewidth', 2);
set(gca, 'XTick', log(1./xt), 'XTickLabel', xt);
axis([log(1/20) log(6.5) 1e6*[-0.5 10.5]]);
grid on
title('Years 1-50');
subplot(4,2,5);
tem = cth(601:1200);
nfft = 256; nover = nfft/2; n2 = nfft/2;
[p,f] = spectrum(tem, nfft, nover);
rho = (corr(tem, tem, 1)+sqrt(corr(tem, tem, 2))+...
corr(tem, tem, 3)^(1/3))/3;
rn = (1-rho^2) ./ (1-2*rho*cos((pi/n2)*(0:n2))+rho^2);
rn = rn / mean(rn(3:(n2+1)));
rn = rn * mean(p(3:(n2+1),1));
% rner1 = rn * 1.7; rner5 = rn * 2.5;
rner1 = rn * 2.5; rner5 = rn * 1.92;
dofx = length(tem) / n2;
xt = [128 64 32 16 8 4 2 1 0.5 0.25];
f2 = f.*6;
hh = plot(log(f2), f2.*p(:,1), '-k', ...
log(f2), f2.*rn', '-k', ...
log(f2), f2.*rner5', '--k', ...
log(f2), f2.*rner1', '-.k');
set(hh(1), 'linewidth', 2);
set(gca, 'XTick', log(1./xt), 'XTickLabel', xt);
axis([log(1/20) log(6.5) 1e6*[-0.5 10.5]]);
grid on
title('Years 51-100');
subplot(4,2,7);
tem = cth(1:1200);
nfft = 256; nover = nfft/2; n2 = nfft/2;
[p,f] = spectrum(tem, nfft, nover);
rho = (corr(tem, tem, 1)+sqrt(corr(tem, tem, 2))+...
corr(tem, tem, 3)^(1/3))/3;
rn = (1-rho^2) ./ (1-2*rho*cos((pi/n2)*(0:n2))+rho^2);
rn = rn / mean(rn(3:(n2+1)));
rn = rn * mean(p(3:(n2+1),1));
rner1 = rn * 1.7; rner5 = rn * 2.5;
dofx = length(tem) / n2;
xt = [128 64 32 16 8 4 2 1 0.5 0.25];
f2 = f.*6;
hh = plot(log(f2), f2.*p(:,1), '-k', ...
log(f2), f2.*rn', '-k', ...
log(f2), f2.*rner5', '--k', ...
log(f2), f2.*rner1', '-.k');
set(hh(1), 'linewidth', 2);
set(gca, 'XTick', log(1./xt), 'XTickLabel', xt);
axis([log(1/20) log(6.5) 1e6*[-0.5 10.5]]);
grid on
title('Years 1-100');
% Plot regressions and such
nc = netcdf('sal.cdf', 'nowrite');
h = nc{'h',1}(:,:,:,:);
nc = close(nc);
h = squeeze(h-15000);
h = h(:,:,2:3:180);
h = detrend(h);
[ntim, nlat, nlon] = size(h);
h = reshape(h, ntim, nlat*nlon);
h2 = rave(rave(h, 25), 37);
h = reshape(h, ntim, nlat, nlon);
h2 = reshape(h2, ntim, nlat, nlon);
cth2 = rave(rave(cth, 25), 37);
cd /home/disk/hayes2/dvimont/ocean/coup/coup
%save h_lph_thinned_coup.mat h h2 cth cth2 lat lon
load h_lph_thinned_coup.mat
[reg1, c1] = regress_eof(h, cth, 1);
[reg2, c2] = regress_eof(h-h2, cth-cth2, 1);
[reg3, c3] = regress_eof(h2, cth2, 1);
default_global; FRAME = [0 160 -45 45]; XAX = lon(2:3:180);
subplot(4,2,4); cla;
gcont(reg1(1,:,:), 250, 1);
title('Unfiltered')
ylabel('2.5m Contours')
set(gca, 'XTick', [0:40:160], 'XTickLabel', [-80:40:80]);
subplot(4,2,6); cla;
gcont(reg2(1,:,:), 250, 1);
title('5yr HP')
ylabel('2.5m Contours')
set(gca, 'XTick', [0:40:160], 'XTickLabel', [-80:40:80]);
subplot(4,2,8); cla;
gcont(reg3(1,:,:), 100, 1);
title('5yr LP')
ylabel('1m Contours')
set(gca, 'XTick', [0:40:160], 'XTickLabel', [-80:40:80]);
cd ~/matlab/Ocean/Figs
print -dps2 coup_regs_ct.ps
figure(2); fo(1); clf;
nlag = length(-3:3);
for i = 1:4;
subplot(4,2,(2*i)-1);
gcont(reg2(i,:,:), 200, 1);
subplot(4,2,(2*i));
gcont(reg3(i,:,:), 50, 1);
end
figure(3); fo(1); clf;
nlag = length(-3:3);
for i = 1:3;
subplot(4,2,(2*i)-1);
gcont(reg2(i+4,:,:), 200, 1);
subplot(4,2,(2*i));
gcont(reg3(i+4,:,:), 25, 1);
end
figure(2); fo(1); clf;
[xk1, yk1] = keep_var([-1 161 -20 -15]);
tem = squeeze(mean(h(:,yk1,xk1), 2));
subplot(1,3,1);
surface(lon(2:3:180), [1:1200]./12, tem);
shading flat
axis([0 160 0 100])
set(gca, 'XTick', 0:40:160)
caxis([-3000 3000])
[xk1, yk1] = keep_var([-1 161 -10 -5]);
tem = squeeze(mean(h(:,yk1,xk1), 2));
subplot(1,3,2);
surface(lon(2:3:180), [1:1200]./12, tem);
shading flat
axis([0 160 0 100])
set(gca, 'XTick', 0:40:160)
caxis([-3000 3000])
[xk1, yk1] = keep_var([-1 161 -5 5]);
tem = squeeze(mean(h(:,yk1,xk1), 2));
subplot(1,3,3);
surface(lon(2:3:180), [1:1200]./12, tem);
shading flat
axis([0 160 0 100])
set(gca, 'XTick', 0:40:160)
caxis([-3000 3000])
figure(3); fo(1); clf;
[xk1, yk1] = keep_var([-1 161 -20 -15]);
tem = squeeze(mean(h2(:,yk1,xk1), 2));
subplot(1,3,1);
surface(lon(2:3:180), [1:1200]./12, tem);
shading flat
axis([0 160 0 100])
set(gca, 'XTick', 0:40:160)
caxis([-1000 1000])
[xk1, yk1] = keep_var([-1 161 -10 -5]);
tem = squeeze(mean(h2(:,yk1,xk1), 2));
subplot(1,3,2);
surface(lon(2:3:180), [1:1200]./12, tem);
shading flat
axis([0 160 0 100])
set(gca, 'XTick', 0:40:160)
caxis([-1000 1000])
[xk1, yk1] = keep_var([-1 161 -5 5]);
tem = squeeze(mean(h2(:,yk1,xk1), 2));
subplot(1,3,3);
surface(lon(2:3:180), [1:1200]./12, tem);
shading flat
axis([0 160 0 100])
set(gca, 'XTick', 0:40:160)
caxis([-1000 1000])