Global Index (short | long) | Local contents | Local Index (short | long)
[udhdx, vdhdy, lat_out, lon_out, depth_out] = ...
[up_dhbardx, vp_dhbardy, lat_out, lon_out, depth_out] = ...
get_uprime_gradtbar(pcs, lims, lags, tim, lev);
| This function calls | |
|---|---|
function [udhdx, vdhdy, lat_out, lon_out, depth_out] = ...
get_uprime_gradtbar(pcs, lims, lags, tim, lev);
nfrm = length(lags);
if nargin == 1;
lims = [110 300 -60 60]; nfrm = 6; tim = 101:550; lev = 1:7;
elseif nargin == 2;
nfrm = 6; tim = 101:550; lev = 1:7;
elseif nargin == 3;
tim = 101:550; lev = 1:7;
elseif nargin == 4;
lev = 1:7;
end
cdtem = ['cd ' eval('pwd')];
cd /home/disk/hayes2/dvimont/csiro/data
% First get tbar
filin = 'temp_L1-10.nc';
nc = netcdf(filin, 'nowrite');
depth = nc{'depth'}(:);
latt = nc{'latitude'}(:);
lont = nc{'longitude'}(:);
[xk, yk] = keep_var(lims, lont, latt);
temp = nc{'temp'}(tim, lev, yk, xk);
mv = nc{'temp'}.missing_value(:);
nc = close(nc);
temp(temp == mv) = NaN;
latt = latt(yk); lont = lont(xk);
temp = squeeze(mean2(temp));
% Now get u' and v'
filin = 'u_L1-10.nc';
nc = netcdf(filin, 'nowrite');
latu = nc{'latitude'}(:);
lonu = nc{'longitude'}(:);
[xk, yk] = keep_var(lims, lonu, latu);
u = nc{'u'}(tim,lev,yk,xk);
mv = nc{'u'}.missing_value(:);
nc = close(nc);
u(u == mv) = NaN;
latu = latu(yk); lonu = lonu(xk);
filin = 'ue_L1-10.nc';
nc = netcdf(filin, 'nowrite');
yk2 = [yk; max(yk)+1];
ue = nc{'ue'}(tim,lev,yk2,xk);
mv = nc{'ue'}.missing_value(:);
nc = close(nc);
ue(ue == mv) = NaN;
[ntim, nlev, nlat, nlon] = size(ue);
ue = (ue(:,:,1:(nlat-1),:) + ue(:,:,2:nlat,:))/2;
u = u + ue;
clear ue;
[u, climu] = remove_mean(u/100);
filin = 'v_L1-10.nc';
nc = netcdf(filin, 'nowrite');
v = nc{'v'}(tim,lev,yk,xk);
mv = nc{'v'}.missing_value(:);
nc = close(nc);
v(v == mv) = NaN;
filin = 've_L1-10.nc';
nc = netcdf(filin, 'nowrite');
xk2 = [xk; max(xk)+1];
ve = nc{'ve'}(tim,lev,yk,xk2);
mv = nc{'ve'}.missing_value(:);
nc = close(nc);
ve(ve == mv) = NaN;
[ntim, nlev, nlat, nlon] = size(ve);
ve = (ve(:,:,:,1:(nlon-1)) + ve(:,:,:,2:nlon))/2;
v = v + ve;
clear ve;
[v, climv] = remove_mean(v/100);
u = squeeze(u); v = squeeze(v);
[ntim, nlev, nlat, nlon] = size(u);
depth = depth/100;
% PCS should be sent in the command line
%cd /home/disk/hayes2/dvimont/csiro/matlab_data/Heat_Content
%load LP10_L1-7_CEOF.mat; tit = 'Lowpass Filtered Data ( > 10 Years )';
%load HP8_L1-7_CEOF.mat; tit = 'Highpass Filtered Data ( < 8 Years )';
%load RAW_L1-7_CEOF.mat; tit = 'Unfiltered Data';
% Get regressions
ureg = regress_eof(u, pcs, lags);
vreg = regress_eof(v, pcs, lags);
%lag = 0; lg = lag*pi/180; lg2 = 1;
%num = 1; lind = 1;
% Store all regressions under one variable
%[nlev, nlat, nlon] = size(squeeze(climu));
%j = sqrt(-1);
%timeseries = sqrt(2)*pcs(:,num)./std(pcs(:,num));
%u = reshape(u, ntim, nlev*nlat*nlon);
%v = reshape(v, ntim, nlev*nlat*nlon);
%clear temtim ureg vreg
%for i = 1:nfrm
% wgt = conj(exp(j * ((i-1) * pi/(lg2*nfrm) + lg) ));
% temtim(:,i) = squeeze(real(wgt .* timeseries));
% ureg(i, :) = temtim(:,i)' * u ./ ntim;
% vreg(i, :) = temtim(:,i)' * v ./ ntim;
%end
%ureg = reshape(ureg, nfrm, nlev, nlat, nlon);
%vreg = reshape(vreg, nfrm, nlev, nlat, nlon);
% Get dtdx and dtdy
[nlev, nlat, nlon] = size(temp);
global DEGREE RADIAN RADUS
clear ty tx dtdy dtdx j
for lind = 1:nlev;
for j = 1:nlat;
ty(j,:) = interp1(lont, squeeze(temp(lind,j,:)), lonu)';
end
for j = 1:nlon;
tx(:,j) = interp1(latt, squeeze(temp(lind,:,j)), latu);
end
[dtdy(lind,:,:) temlat] = sph_grady1(ty, RADIAN*latt, RADIAN*lonu, 0);
[dtdx(lind,:,:) temlon] = sph_gradx1(tx, RADIAN*latu, RADIAN*lont, 0);
end
% Get thickness of layer, to convert t' to Heat'
[nlev, nlat, nlon] = size(dtdx);
d1 = depth(1);
depthw = [];
for i = 1:nlev;
depthw = [depthw; depth(i)+d1];
d1 = depth(i+1) - depth(i) - d1;
end
dz = diff([0; depthw]);
dz = reshape((dz * ones(1, nlat*nlon)), nlev, nlat, nlon);
dtdx = dtdx.*dz;
dtdy = dtdy.*dz;
% Get ubar*dtdx and vbar*dtdy, and sum vertically
[nlev, nlat, nlon] = size(dtdx);
ureg = reshape(ureg, nfrm, nlev*nlat*nlon);
dtdx = reshape(dtdx, 1, nlev*nlat*nlon);
udhdx = ureg .* (ones(nfrm, 1) * dtdx);
udhdx = shiftdim(reshape(udhdx, nfrm, nlev, nlat, nlon), 1);
udhdx = shiftdim(squeeze(sum(udhdx)), 2);
[nlev, nlat, nlon] = size(dtdy);
vreg = reshape(vreg, nfrm, nlev*nlat*nlon);
dtdy = reshape(dtdy, 1, nlev*nlat*nlon);
vdhdy = vreg .* (ones(nfrm, 1) * dtdy);
vdhdy = shiftdim(reshape(vdhdy, nfrm, nlev, nlat, nlon), 1);
vdhdy = shiftdim(squeeze(sum(vdhdy)), 2);
% Convert units to heat
cwat = 4.218e3; % heat capacity of liquid water, J/(kg K)
rhowat = 1e3; % density of liquid water, kg/m^3
vdhdy = cwat * rhowat * vdhdy;
udhdx = cwat * rhowat * udhdx;
lat_out = latu;
lon_out = lonu;
depth_out = depth(lev);
eval(cdtem);