Global Index (short | long) | Local contents | Local Index (short | long)
Load PCS:
| This script calls | |
|---|---|
clear
cd /home/disk/tao/dvimont/matlab/CSIRO/Thesis/Data
load HP10_detrend_L1-7_CEOF_yr101-1000.mat
for levind = [3];
cd ~/matlab/CSIRO/Heat/Old_routines
%nfrm = [-15:15];
nfrm = 6;
%tim = 101:550;
tim = 101:1000;
if levind == 1;
lev = 1:3;
elseif levind == 2;
lev = 1:4;
elseif levind == 3;
lev = 1:7;
end
%for biff = 2:2;
biff = 2;
% First, get WSTP values
global RADUS DEGREE RADIAN
if biff == 1;
lims = [117.5 169 7.5 24];
% lims = [155 210 -15 3];
% lims = [178 210 29 47];
elseif biff == 2;
% lims = [116 148 4.5 25];
% lims = [178 272 -5 5];
end
% Get ubar_dHCprimedx and vbar_dHCprimedy
cd ~/matlab/CSIRO/Heat/Old_routines
[utem, vtem, lat_out, lon_out, depth_out] = ...
get_ubar_gradtprime(pcs, lims, nfrm, tim, lev);
[ntem, nlat, nlon] = size(utem);
dx = mean(diff(lon_out)); dy = mean(diff(lat_out));
lon_out = [lon_out; max(lon_out)+dx] - dx/2;
lat_out = [lat_out; max(lat_out)+dy] - dy/2;
A = RADUS * ( diff(sin(RADIAN * lat_out)) ) * ( RADIAN * diff(lon_out) )';
A = reshape(A, 1, nlat*nlon);
ubtpw13 = -1 * A * reshape(utem, ntem, nlat*nlon)' ./ sum(A);
vbtpw13 = -1 * A * reshape(vtem, ntem, nlat*nlon)' ./ sum(A);
% Get uprime_dHCbardx and vprime_dHCbardy
[utem, vtem, lat_out, lon_out, depth_out] = ...
get_uprime_gradtbar(pcs, lims, nfrm, tim, lev);
[ntem, nlat, nlon] = size(utem);
dx = mean(diff(lon_out)); dy = mean(diff(lat_out));
lon_out = [lon_out; max(lon_out)+dx] - dx/2;
lat_out = [lat_out; max(lat_out)+dy] - dy/2;
A = RADUS * ( diff(sin(RADIAN * lat_out)) ) * ( RADIAN * diff(lon_out) )';
A = reshape(A, 1, nlat*nlon);
uptbw13 = -1 * A * reshape(utem, ntem, nlat*nlon)' ./ sum(A);
vptbw13 = -1 * A * reshape(vtem, ntem, nlat*nlon)' ./ sum(A);
% Get wbar_dHCprimedz
[wtem, lat_out, lon_out] = ...
get_wbar_dtprimedz(pcs, lims, nfrm, tim, lev);
[ntem, nlat, nlon] = size(wtem);
dx = mean(diff(lon_out)); dy = mean(diff(lat_out));
lon_out = [lon_out; max(lon_out)+dx] - dx/2;
lat_out = [lat_out; max(lat_out)+dy] - dy/2;
A = RADUS * ( diff(sin(RADIAN * lat_out)) ) * ( RADIAN * diff(lon_out) )';
A = reshape(A, 1, nlat*nlon);
wbtpw13 = -1 * A * reshape(wtem, ntem, nlat*nlon)' ./ sum(A);
% Get wprime_dHCbardz
[wtem, lat_out, lon_out] = ...
get_wprime_dtbardz(pcs, lims, nfrm, tim, lev);
[ntem, nlat, nlon] = size(wtem);
dx = mean(diff(lon_out)); dy = mean(diff(lat_out));
lon_out = [lon_out; max(lon_out)+dx] - dx/2;
lat_out = [lat_out; max(lat_out)+dy] - dy/2;
A = RADUS * ( diff(sin(RADIAN * lat_out)) ) * ( RADIAN * diff(lon_out) )';
A = reshape(A, 1, nlat*nlon);
wptbw13 = -1 * A * reshape(wtem, ntem, nlat*nlon)' ./ sum(A);
% Get dHC/dt
tim2 = [min(tim)-1 tim];
[heat, lat_out, lon_out, depth, middepth] = getheat(lims, lev, tim2);
[ntim, nlat, nlon] = size(heat);
[heat, tem] = remove_mean(heat);
heat = reshape(heat, ntim, nlat*nlon);
heat2 = (heat(3:ntim,:)-heat(1:(ntim-2),:)) / (2*24*3600*365);
heat2 = regress_ceof(heat2, pcs(1:ntim-2, :), nfrm);
heat = regress_ceof(heat(tim-min(tim)+2, :), pcs, nfrm);
dx = mean(diff(lon_out)); dy = mean(diff(lat_out));
lon_out = [lon_out; max(lon_out)+dx] - dx/2;
lat_out = [lat_out; max(lat_out)+dy] - dy/2;
A = RADUS * ( diff(sin(RADIAN * lat_out)) ) * ( RADIAN * diff(lon_out) )';
A = reshape(A, 1, nlat*nlon);
heat = A * heat' ./ (sum(A));
heat2 = A * heat2' ./ sum(A);
% Get HFLX
hflx = getnc('heat', lims, 1, tim);
[hflx, tem] = remove_mean(hflx);
[lat_out, lon_out] = getll('heat', lims);
hflx = regress_ceof(hflx, pcs, nfrm);
[ntem, nlat, nlon] = size(hflx);
[tem, nlat, nlon] = size(hflx);
dx = mean(diff(lon_out)); dy = mean(diff(lat_out));
lon_out = [lon_out; max(lon_out)+dx] - dx/2;
lat_out = [lat_out; max(lat_out)+dy] - dy/2;
A = RADUS * ( diff(sin(RADIAN * lat_out)) ) * ( RADIAN * diff(lon_out) )';
A = reshape(A, 1, nlat*nlon);
hflx = 1 * A * reshape(hflx, ntem, nlat*nlon)' ./ sum(A);
% Plot the data
if 0;
ph = 180/ntem * [0:(2*nfrm)];
%figure(1); figure_orient(1);
sp([]);
h1 = plot(...
ph, [ubtpw133 -1*ubtpw133 ubtpw133(1)], 'x-k', ...
ph, [vbtpw133 -1*vbtpw133 vbtpw133(1)], 'o-k', ...
ph, [wbtpw133 -1*wbtpw133 wbtpw133(1)], 'v-k', ...;
ph, [uptbw133 -1*uptbw133 uptbw133(1)], '*--k', ...
ph, [vptbw133 -1*vptbw133 vptbw133(1)], 's--k', ...
ph, [wptbw133 -1*wptbw133 wptbw133(1)], '^--k', ...
ph, [heat23 -1*heat23 heat23(1) ], 'd-.k', ...
ph, [hflx3 -1*hflx3 hflx3(1) ], '>-.k');
hold on;
h3 = plot(ph, 0.5*[heat, -1*heat, heat(1)], '.-k');
set(h3, 'linewidth', 2);
hold off;
grid on
if biff == 1;
axis([0 360 -.75 .75]);
% title([tit ': WSTP [117.5E 170E 7.5N 22.5N] 80-270m Heat Budget Terms']);
% title([tit ': Western EQ PAC [155E 150W 5S 5N] Heat Budget Terms']);
title([tit ': ML anomaly [180 150W 30N 45N] Heat Budget Terms']);
elseif biff == 2;
axis([0 360 -1.5 1.5]);
% title([tit ': CT [180 270E 5S 5N] Level 80:270m Heat Budget Terms']);
title([tit ': Eastern EQ PAC [150W 90W 5S 5N] Heat Budget Terms']);
end
xlabel('Phase');
ylabel('W m^-^2');
set(gca, 'XTick', [0:30:360]);
h2 = legend(h1, 'Ubar dHC''dx', 'Vbar dHC''dy', 'Wbar dHC''dz', 'U'' dHCbardx',...
'V'' dHCbardy', 'W'' dHCbardz', 'dHCdt', 'Net HFLX');
end;
if 0;
cd /home/disk/tao/dvimont/matlab/CSIRO/Heat/Old_routines/Plot_HCadv
if biff2 == 1;
print -dps2 LP10_Tendency_Terms_WSTP_EQ.ps
elseif biff2 == 2;
print -dps2 HP8_Tendency_Terms_WSTP_EQ.ps
elseif biff2 == 3;
print -dps2 RAW_Tendency_Terms_WSTP_EQ.ps
end
end
% Swap variables: var3 means lev 1-3
if levind == 1;
ubtp_wstp3 = ubtp_wstp;
vbtp_wstp3 = vbtp_wstp ;
wbtp_wstp3 = wbtp_wstp;
uptb_wstp3 = uptb_wstp;
vptb_wstp3 = vptb_wstp ;
wptb_wstp3 = wptb_wstp ;
heat23 = heat2 ;
hflx3 = hflx ;
heat3 = heat;
% Swap variables: var4 means lev 1-4
elseif levind == 2;
ubtp_wstp4 = ubtp_wstp;
vbtp_wstp4 = vbtp_wstp ;
wbtp_wstp4 = wbtp_wstp;
uptb_wstp4 = uptb_wstp;
vptb_wstp4 = vptb_wstp ;
wptb_wstp4 = wptb_wstp ;
heat24 = heat2 ;
hflx4 = hflx ;
heat4 = heat;
% Get lower level: var7 means lev 1-7
elseif levind == 3;
ubtp_wstp7 = ubtp_wstp;
vbtp_wstp7 = vbtp_wstp ;
wbtp_wstp7 = wbtp_wstp;
uptb_wstp7 = uptb_wstp;
vptb_wstp7 = vptb_wstp ;
wptb_wstp7 = wptb_wstp ;
heat27 = heat2 ;
hflx7 = hflx ;
heat7 = heat;
end
end
% Get lower level: var37 means lev 3-7
ubtp_wstp37 = ubtp_wstp7 - ubtp_wstp3;
vbtp_wstp37 = vbtp_wstp7 - vbtp_wstp3;
wbtp_wstp37 = wbtp_wstp7 - wbtp_wstp3;
uptb_wstp37 = uptb_wstp7 - uptb_wstp3;
vptb_wstp37 = vptb_wstp7 - vptb_wstp3;
wptb_wstp37 = wptb_wstp7 - wptb_wstp3;
heat37 = heat7 - heat3 ;
heat237 = heat27 - heat23;
hflx37 = hflx7 - hflx3 ;
tot37 = ubtp_wstp37 + vbtp_wstp37 + wbtp_wstp37 + ...
uptb_wstp37 + vptb_wstp37 + wptb_wstp37 + hflx37;
% Normalize if necessary
tot3 = ubtp_wstp3 + vbtp_wstp3 + wbtp_wstp3 + ...
uptb_wstp3 + vptb_wstp3 + wptb_wstp3 + hflx3;
tot7 = ubtp_wstp7 + vbtp_wstp7 + wbtp_wstp7 + ...
uptb_wstp7 + vptb_wstp7 + wptb_wstp7 + hflx7;
tot4 = ubtp_wstp4 + vbtp_wstp4 + wbtp_wstp4 + ...
uptb_wstp4 + vptb_wstp4 + wptb_wstp4 + hflx4;
% For lagged regressions:
% Plot the data
ph = nfrm;
figure(1); figure_orient(1);
sp(1);
h1 = plot(...
ph, -1*ubtp_wstp3, 'x-k', ...
ph, -1*vbtp_wstp3, 'o-k', ...
ph, -1*wbtp_wstp3, 'v-k', ...;
ph, -1*uptb_wstp3, '*--k', ...
ph, -1*vptb_wstp3, 's--k', ...
ph, -1*wptb_wstp3, '^--k', ...
ph, -1*heat23, 'd-.k', ...
ph, -1*hflx3, '>-.k');
hold on;
h3 = plot(ph, -1*0.5*heat3, '.-k');
set(h3, 'linewidth', 2);
hold off;
grid on
if biff == 1;
axis([-15 15 -.4 .4]);
% title([tit ': WSTP [115E 150E 5N 25N] 0-80m Heat Budget Terms']);
% title([tit ': Western EQ PAC [155E 150W 5S 5N] Heat Budget Terms']);
title([tit ': ML anomaly [180 150W 30N 45N] Heat Budget Terms']);
elseif biff == 2;
axis([-15 15 -1.2 1.2]);
% title([tit ': CT [180 270E 5S 5N] level 80-270m Heat Budget Terms']);
title([tit ': Eastern EQ PAC [150W 90W 5S 5N] Heat Budget Terms']);
end
xlabel('Lag');
ylabel('W m^-^2');
set(gca, 'XTick', [-20:2:20]);
h2 = legend(h1, 'Ubar dHC''dx', 'Vbar dHC''dy', 'Wbar dHC''dz', 'U'' dHCbardx',...
'V'' dHCbardy', 'W'' dHCbardz', 'dHCdt', 'Net HFLX');
end;
ph = nfrm;
figure(1); figure_orient(1);
sp(2);
h1 = plot(...
ph, ubtp_wstp3./heat23, 'x-k', ...
ph, vbtp_wstp3./heat23, 'o-k', ...
ph, wbtp_wstp3./heat23, 'v-k', ...;
ph, uptb_wstp3./heat23, '*--k', ...
ph, vptb_wstp3./heat23, 's--k', ...
ph, wptb_wstp3./heat23, '^--k', ...
ph, heat23./heat23, 'd-.k', ...
ph, hflx3./heat23, '>-.k');
cd /home/disk/tao/dvimont/matlab/CSIRO/Heat/Old_routines/Plot_HCadv/LP10_ttend
print -dps2 CT_ttend_lag_reg_top_bot_yr551-1000.ps
dhcw17dt = heat2;
hcw17 = heat;