Global Index (short | long) | Local contents | Local Index (short | long)
Load variables for regression
| This script calls | |
|---|---|
clear
cd ~/Papers/mlcsiro/matlab/data
load SVD_curlt_slp_MIX_raw.mat
lims = [89 306 -79 90]
cd /home/disk/hayes2/dvimont/csiro/data
nc = netcdf('spsl_qm1_0036-0365.nc', 'nowrite');
[lat, lon, yk, xk] = get_nclatlon2(lims, nc);
slp = nc{'psl', 1}(31:330, 1:12, yk, xk);
nc = close(nc);
nc = netcdf('stax_qm1_0036-0365.nc', 'nowrite');
[lat, lon, yk, xk] = get_nclatlon2(lims, nc);
tax = nc{'tax', 1}(31:330, 1:12, yk, xk);
nc = close(nc);
nc = netcdf('stsu_qm1_0036-0365.nc', 'nowrite');
[lat, lon, yk, xk] = get_nclatlon2(lims, nc);
tsu = nc{'tsu', 1}(31:330, 1:12, yk, xk);
nc = close(nc);
[nyr, nmo, nlat, nlon] = size(slp);
slp2 = NaN*ones(nyr*nmo, nlat, nlon);
tax2 = NaN*ones(nyr*nmo, nlat, nlon);
tsu2 = NaN*ones(nyr*nmo, nlat, nlon);
for i = 1:nyr;
ind = 12*(i-1)+[1:12];
slp2(ind,:,:) = slp(i,:,:,:);
tax2(ind,:,:) = tax(i,:,:,:);
tsu2(ind,:,:) = tsu(i,:,:,:);
end
clear slp tax tsu;
slp2 = annave(slp2);
tax2 = annave(tax2);
tsu2 = annave(tsu2);
slp2 = detrend(slp2);
tax2 = detrend(tax2);
tsu2 = detrend(tsu2);
[ntim, nlat, nlon] = size(slp2);
[xk, yk] = keep_var(limslp, lon, lat);
lims2 = [lims(1:3) 80];
[xk2, yk2] = keep_var(lims2, lon, lat);
% Reweight patterns by sqrt(cos...)
nlat1 = length(yk); nlon1 = length(xk);
tem = reshape(slpu', 10, nlat1, nlon1);
tem = cosweight(tem, lat(yk));
tem = reshape(tem, 10, nlat1*nlon1)';
use_tim = tem(:,1:10);
% 2-12 mo
[b, a] = butter(9, 2/12);
slp = slp2 - filtfilt(b, a, slp2);
tax = tax2 - filtfilt(b, a, tax2);
tsu = tsu2 - filtfilt(b, a, tsu2);
tim_hp = reshape(slp(:,yk,xk), ntim, nlat1*nlon1) * use_tim;
[pat1aa, c1aa] = regress_eof(slp, tim_hp);
[pat2aa, c2aa] = regress_eof(tax, tim_hp);
[pat3aa, c3aa] = regress_eof(tsu, tim_hp);
% 12-24 mo
[b2, a2] = butter(9, 2/12);
[b, a] = butter(9, 2/24);
slp = filtfilt(b2, a2, slp2) - filtfilt(b, a, slp2);
tax = filtfilt(b2, a2, tax2) - filtfilt(b, a, tax2);
tsu = filtfilt(b2, a2, tsu2) - filtfilt(b, a, tsu2);
tim_bp = reshape(slp(:,yk,xk), ntim, nlat1*nlon1) * use_tim;
[pat1bb, c1bb] = regress_eof(slp, tim_bp);
[pat2bb, c2bb] = regress_eof(tax, tim_bp);
[pat3bb, c3bb] = regress_eof(tsu, tim_bp);
% 24-Inf mo
[b, a] = butter(9, 2/24);
slp = filtfilt(b, a, slp2);
tax = filtfilt(b, a, tax2);
tsu = filtfilt(b, a, tsu2);
tim_lp = reshape(slp(:,yk,xk), ntim, nlat1*nlon1) * use_tim;
[pat1cc, c1cc] = regress_eof(slp, tim_lp);
[pat2cc, c2cc] = regress_eof(tax, tim_lp);
[pat3cc, c3cc] = regress_eof(tsu, tim_lp);
cd ~/Papers/mlcsiro/matlab/data
save MIX_hp_bp_lp_ppcs_regs.mat tim_hp tim_bp tim_lp ...
pat1aa c1aa pat2aa c2aa pat3aa c3aa ...
pat1bb c1bb pat2bb c2bb pat3bb c3bb ...
pat1cc c1cc pat2cc c2cc pat3cc c3cc ...
lims limslp xk yk lat lon
% Plot data
default_global; FRAME = [95 300 -45 75];
[tem, tem, tem, lm] = getll('temp', lims);
global GRDX_SPACING GRDY_SPACING
GRDX_SPACING = 45; GRDY_SPACING = 20;
figure(1); fl(1); clf;
clim = 0.25;
subplot(3,3,1);
gcont(pat1aa, 0.5);
dc2(lm, 0.5, -1000);
color_shade(squeeze(c1aa.^2), clim, 0.7*[1 1 1]);
t1 = title('< ML\_SLP, SLP > 0.5 mb std^-^1');
y1 = ylabel('2-12 Mo.'); set(y1, 'fontsize', 9);
set(gca, 'fontsize', 8); set(t1, 'fontsize', 9);
% h = drawbox(limslp, 'k'); set(h, 'linewidth', 2);
subplot(3,3,2);
gcont(pat2aa, 0.01, 1);
dc2(lm, 0.5, 1000);
color_shade(squeeze(c2aa.^2), clim, 0.7*[1 1 1]);
t1 = title('< ML\_SLP, \tau x > 0.1 dyne m^-^2 std^-^1');
set(gca, 'fontsize', 8); set(t1, 'fontsize', 9);
subplot(3,3,3);
gcont(pat3aa, 0.01);
dc2(lm, 0.5, 1000);
color_shade(squeeze(c3aa.^2), clim, 0.7*[1 1 1]);
t1 = title('< ML\_SLP, SST > 0.01 K std^-^1');
set(gca, 'fontsize', 8); set(t1, 'fontsize', 9);
subplot(3,3,4);
gcont(pat1bb, 0.25);
dc2(lm, 0.5, -1000);
color_shade(squeeze(c1bb.^2), clim, 0.7*[1 1 1]);
t1 = title('< ML\_SLP, SLP > 0.25 mb std^-^1');
y1 = ylabel('12-24 Mo.'); set(y1, 'fontsize', 9);
set(gca, 'fontsize', 8); set(t1, 'fontsize', 9);
subplot(3,3,5);
gcont(pat2bb, 0.0025, 1);
dc2(lm, 0.5, 1000);
color_shade(squeeze(c2bb.^2), clim, 0.7*[1 1 1]);
t1 = title('< ML\_SLP, \tau x > 0.025 dyne m^-^2 std^-^1');
set(gca, 'fontsize', 8); set(t1, 'fontsize', 9);
subplot(3,3,6);
gcont(pat3bb, 0.05);
dc2(lm, 0.5, 1000);
color_shade(squeeze(c3bb.^2), clim, 0.7*[1 1 1]);
t1 = title('< ML\_SLP, SST > 0.05 K std^-^1');
set(gca, 'fontsize', 8); set(t1, 'fontsize', 9);
subplot(3,3,7);
gcont(pat1cc, 0.25);
dc2(lm, 0.5, -1000);
color_shade(squeeze(c1cc.^2), clim, 0.7*[1 1 1]);
t1 = title('< ML\_SLP, SLP > 0.25 mb std^-^1');
y1 = ylabel('24-Inf Mo.'); set(y1, 'fontsize', 9);
set(gca, 'fontsize', 8); set(t1, 'fontsize', 9);
subplot(3,3,8);
gcont(pat2cc, 0.0025, 1);
dc2(lm, 0.5, 1000);
color_shade(squeeze(c2cc.^2), clim, 0.7*[1 1 1]);
t1 = title('< ML\_SLP, \tau x > 0.025 dyne m^-^2 std^-^1');
set(gca, 'fontsize', 8); set(t1, 'fontsize', 9);
subplot(3,3,9);
gcont(pat3cc, 0.1);
dc2(lm, 0.5, 1000);
color_shade(squeeze(c3cc.^2), clim, 0.7*[1 1 1]);
t1 = title('< ML\_SLP, SST > 0.1 K std^-^1');
set(gca, 'fontsize', 8); set(t1, 'fontsize', 9);
cd ~/Papers/mlcsiro/matlab/model/Figs
print -dps2 HP_BP_LP_SLP_TX_SLP_pslpx.ps
nfft = 64*12; nover = 7*nfft/8;
[p, f] = spectrum(slp_tim, nfft, nover);
f2 = f*6;
n2 = nfft/2;
rho = (corr(slp_tim, slp_tim, 1)+sqrt(corr(slp_tim,slp_tim,2)))/2;
rn = (1-rho^2) ./ (1-2*rho*cos((pi/n2)*(0:n2))+rho^2);
rn = rn / mean(rn);
rn = rn * mean(p(:,1));
rner1 = rn * 2.05; rner5 = rn * 1.68;
ind = [1:n2]+1;
figure(4); fl(1);
subplot(2,1,1);
plot(log(f2(ind)), f2(ind).*p(ind,1), '-k', ...
log(f2(ind)), f2(ind).*rn(ind)', '--k');
plot(f2, p(:,1), 'k', f2, rn, '--k');
% Determine equal variances
n2 = 4096;
rho = (corr(slp_tim, slp_tim, 1)+sqrt(corr(slp_tim,slp_tim,2)))/2;
rn = (1-rho^2) ./ (1-2*rho*cos((pi/n2)*(0:n2))+rho^2);
rn = rn / mean(rn);
rn = rn * mean(p(:,1));
ftem = 6*[0:(1/n2):1];