Global Index (short | long) | Local contents | Local Index (short | long)
Load PCS:
This script calls | |
---|---|
clear biff2 = 1; cd /home/disk/hayes2/dvimont/csiro/matlab_data/Heat_Content if biff2 == 1; load LP10_detrend_L1-7_yr550-1000.mat; tit = 'Lowpass Filtered Data ( > 10 Years )'; % load BP10-60yr_L1-7_CEOF.mat; tit = 'Lowpass Filtered Data ( > 10 Years )'; elseif biff2 == 2; load HP8_L1-7_CEOF.mat; tit = 'Highpass Filtered Data ( < 8 Years )'; elseif biff2 == 3; load RAW_L1-7_CEOF.mat; tit = 'Unfiltered Data'; end for levind = [3]; cd ~/matlab/CSIRO/Heat/Old_routines nfrm = [-15:15]; %tim = 101:550; tim = 551: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 [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); ubtp_wstp = -1 * A * reshape(utem, ntem, nlat*nlon)' ./ sum(A); vbtp_wstp = -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); uptb_wstp = -1 * A * reshape(utem, ntem, nlat*nlon)' ./ sum(A); vptb_wstp = -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); wbtp_wstp = -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); wptb_wstp = -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(lev, tim2, lims); [ntim, nlat, nlon] = size(heat); [heat, tem] = remove_mean(heat); heat = reshape(heat, ntim, nlat*nlon); % if biff == 1; % [b, a] = butter(6, 2/8); heat = heat - filtfilt(b, a, heat); % elseif biff == 2; % [b, a] = butter(6, 2/10); heat = filtfilt(b, a, heat); % end heat2 = (heat(3:ntim,:)-heat(1:(ntim-2),:)) / (2*24*3600*365); heat2 = regress_eof(heat2, pcs(1:ntim-2, :), nfrm); heat = regress_eof(heat(tim-min(tim)+1, :), 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) * 2*3600*24*365); 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_eof(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, [ubtp_wstp3 -1*ubtp_wstp3 ubtp_wstp3(1)], 'x-k', ... ph, [vbtp_wstp3 -1*vbtp_wstp3 vbtp_wstp3(1)], 'o-k', ... ph, [wbtp_wstp3 -1*wbtp_wstp3 wbtp_wstp3(1)], 'v-k', ...; ph, [uptb_wstp3 -1*uptb_wstp3 uptb_wstp3(1)], '*--k', ... ph, [vptb_wstp3 -1*vptb_wstp3 vptb_wstp3(1)], 's--k', ... ph, [wptb_wstp3 -1*wptb_wstp3 wptb_wstp3(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