Global Index (short | long) | Local contents | Local Index (short | long)
ctlim = [115 285 -30 30];
| This script calls | |
|---|---|
clear
filin = 'temp_A_L1-10.nc';
cd /home/disk/hayes2/dvimont/csiro/data
tim = 101:550;
nc = netcdf(filin, 'nowrite');
lat = nc{'latitude'}(:);
lon = nc{'longitude'}(:);
ctlim = [110 300 -60 60];
[xk, yk] = keep_var(ctlim, lon, lat);
temp = nc{'temp'}(tim,1,yk,xk);
mv = nc{'temp'}.missing_value(:);
nc = close(nc);
temp = squeeze(temp);
[ntim, nlat, nlon] = size(temp);
temp(find(temp == mv)) = NaN * ones(size(find(temp == mv)));
lat = lat(yk); lon = lon(xk);
if 1;
[ntim, nlat, nlon] = size(temp);
temp = shiftdim(temp, 1);
for i = 1:nlat/2
indy = 2*(i-1)+[1:2];
for j = 1:nlon/2
indx = 2*(j-1)+[1:2];
temp2(:,i,j) = squeeze(mean(mean(temp(indy,indx,:))));
end
end
lat = mean([lat(1:2:nlat)'; lat(2:2:nlat)']);
lon = mean([lon(1:2:nlon)'; lon(2:2:nlon)']);
temp = temp2; clear temp2;
[ntim, nlat, nlon] = size(temp);
end;
get_global;
XAX = lon; YAX = lat; FRAME = ctlim;
temp = reshape(temp, ntim, nlat*nlon);
[temp, clim] = remove_mean(temp);
temp = detrend(temp);
temp = temp ./ (ones(ntim,1)*std(temp));
temp = cosweight(temp, lat);
kp = find(~isnan(clim));
temp = reshape(temp, ntim, nlat*nlon);
temp = temp(:, kp);
% Low Pass Filter the data:
if 1;
[b,a]=butter(6,2/4.5);
temp = filtfilt(b, a, temp);
%timekp = 11:1:(ntim-10);
%temp = temp(timekp, :);
[ntim, nkp] = size(temp);
temp = detrend(temp);
end;
[lam, lds, pcs, per] = complex_eof(temp);
cd /home/disk/tao/dvimont/matlab/CSIRO/Data
% save LPceof_T1_yr101-550.mat lam lds pcs per kp lat lon
if 0;
clear
cd /home/disk/hayes2/dvimont/csiro/matlab_data
load LPregmap_3deof_-67.5to60pac.mat
ctlim = FRAME;
get_global; default_global; FRAME = ctlim;
ntim = size(pcs, 1);
[nlev, nlat, nlon] = size(clim);
tim = [101:550];
elseif 0;
clear
cd /home/disk/tao/dvimont/matlab/CSIRO/Data
load LPceof_T1_yr101-550.mat
ctlim = [110 300 -60 60];
get_global; default_global; FRAME = ctlim;
ntim = size(pcs, 1);
tim = [101:550];
else
clear
cd /home/disk/hayes2/dvimont/csiro/matlab_data
load 8yrLPtemp_3dceof_-67.5to60pac.mat
ctlim = [110 300 -67.5 60];
get_global; default_global; FRAME = ctlim;
ntim = size(pcs, 1);
tim = [101:550];
end
cd /home/disk/hayes2/dvimont/csiro/data
tim = [101:550];
filin = ['temp_A_L1-10.nc'];
nc = netcdf(filin, 'nowrite');
lind = 1:10;
depth = nc{'depth'}(lind);
nc = close(nc);
% Get regressions of fluxes:
cd /home/disk/hayes2/dvimont/csiro/data
nc = netcdf('zonal_stress_A_1000_years.nc', 'nowrite');
lat1 = nc{'latitude'}(:);
lon1 = nc{'longitude'}(:);
[xk, yk] = keep_var(FRAME, lon1, lat1);
taux = nc{'smfzon'}(tim,yk,xk);
mv = nc{'smfzon'}.missing_value(:);
nc = close(nc);
taux(find(taux == mv)) = NaN * ones(size(find(taux == mv)));
nc = netcdf('merid_stress_A_1000_years.nc', 'nowrite');
tauy = nc{'smfmer'}(tim,yk,xk);
mv = nc{'smfmer'}.missing_value(:);
nc = close(nc);
tauy(find(tauy == mv)) = NaN * ones(size(find(tauy == mv)));
lat1 = lat1(yk); lon1 = lon1(xk);
nc = netcdf('heat_flux_A_1000_years.nc', 'nowrite');
lat2 = nc{'latitude'}(:);
lon2 = nc{'longitude'}(:);
[xk, yk] = keep_var(FRAME, lon2, lat2);
hflx = nc{'stfht'}(tim,yk,xk);
mv = nc{'stfht'}.missing_value(:);
nc = close(nc);
hflx(find(hflx == mv)) = NaN * ones(size(find(hflx == mv)));
lat2 = lat2(yk); lon2 = lon2(xk);
nc = netcdf('temp_A_L1-10.nc', 'nowrite');
lat3 = nc{'latitude'}(:);
lon3 = nc{'longitude'}(:);
[xk, yk] = keep_var(FRAME, lon3, lat3);
clear temp
for i = 1:4;
temp(:,i,:,:) = nc{'temp'}(tim,(3*(i-1)+1),yk,xk);
end
% temp = nc{'temp'}(tim,:,yk,xk);
mv = nc{'temp'}.missing_value(:);
nc = close(nc);
temp(find(temp == mv)) = NaN * ones(size(find(temp == mv)));
lat3 = lat3(yk); lon3 = lon3(xk);
depth = depth(1:3:10);
[ntim, nlev, nlat, nlon] = size(temp);
hflx2 = hflx + (5.67e-8 * (squeeze(temp(:,1,:,:)+273.15).^4));
temp(:,2,:,:) = hflx2;
f = (2*7.292e-5)*sin(lat1*pi/180);
for i = 1:ntim;
temx = squeeze(taux(i,:,:)) ./ (f * ones(1, length(lon1)));
temy = squeeze(tauy(i,:,:)) ./ (f * ones(1, length(lon1)));
curlt(i,:,:) = sph_curl(temx, temy, lon1, lat1);
end
taux = curlt;
yk = keep_var3([-8 8], lat1);
taux(:,yk,:) = NaN;
hflx = reshape(hflx, ntim, size(hflx,2)*size(hflx,3));
taux = reshape(taux, ntim, size(taux,2)*size(taux,3));
tauy = reshape(tauy, ntim, size(tauy,2)*size(tauy,3));
temp = reshape(temp, ntim, nlev, size(temp,3)*size(temp,4));
[hflx, clim] = remove_mean(hflx);
[taux, clim] = remove_mean(taux); kpx = find(~isnan(clim));
[tauy, clim] = remove_mean(tauy); kpy = find(~isnan(clim));
[temp, clim] = remove_mean(temp);
if 0;
[b, a] = butter(6, 2/8);
temp = temp - filtfilt(b, a, temp);
hflx = hflx - filtfilt(b, a, hflx);
taux = taux - filtfilt(b, a, taux);
tauy = tauy - filtfilt(b, a, tauy);
else
[b, a] = butter(6, 2/8);
temp = filtfilt(b, a, temp);
hflx = filtfilt(b, a, hflx);
taux = filtfilt(b, a, taux);
tauy = filtfilt(b, a, tauy);
end
% Plot the data
timeseries = sqrt(2)*pcs(:,1)./std(pcs(:,1));
cint1 = [0.05 0.05 0.05 0.025 1.2 .0005];
for i = 6:6;
fig1 = floor(i/2 - 1/2);
if i < 5;
tem1 = squeeze(temp(:,i,:));
tem2 = std(tem1);
cint = cint1(i); clev = [-1:cint:1];
xlab = [num2str(cint) ' K (std)^-^1'];
XAX = lon3; YAX = lat3;
nlat = length(YAX); nlon = length(XAX);
tit = ['Depth: ' num2str(depth(i)/100) ' m; '];
elseif i == 5;
tem1 = taux;
tem2 = tauy;
cint = cint1(i);
XAX = lon1; YAX = lat1;
nlat = length(YAX); nlon = length(XAX);
tit = ['Wind Stress, '];
elseif i == 6;
tem1 = taux;
yk = keep_var3([-10 10], lat1);
tem1(:,yk,:) = NaN;
tem2 = std(tem1);
cint = cint1(i); clev = [-.02:cint:.02];
xlab = [num2str(cint) ' K (std)^-^1'];
XAX = lon1; YAX = lat1;
nlat = length(YAX); nlon = length(XAX);
tit = ['Curl( TAU / f ), '];
% elseif i == 6;
% tem1 = 1*hflx;
% tem2 = std(hflx);
% cint = cint1(i); clev = [-18:cint:18];
% xlab = [num2str(cint) ' W m^-^2 (std)^-^1'];
% XAX = lon2; YAX = lat2;
% nlat = length(YAX); nlon = length(XAX);
% tit = ['TOT Heat Flux (Pos. Up), '];
end;
% dofy = get_dof(tem1) - 3;
if i ~= 5; % Perform regular t-test on correllation coefficient
lag = 0; lg = lag*pi/180;
num = 1; nfrm = 6;
for ind = 1:nfrm;
fig2 = floor(ind/2 - 0.5);
wgt = conj(exp(j * ((ind-1) * pi/nfrm + lg) ));
temtim = real(wgt .* timeseries);
% dofx = get_dof(temtim) - 3;
tem = temtim'*tem1./ntim;
% ccoef = tem./tem2;
% tstat = ccoef ./ (sqrt((1-ccoef.^2)./dofy));
% score = tscore(dofx, 2.5);
[ccoef, tstat] = corr_sig(temtim, tem1);
tem = reshape(tem, nlat, nlon);
tstat = reshape(tstat, nlat, nlon);
% figure(3*fig1 + fig2+1); figure_orient(1);
% subplot(2,2,(rem(ind-1,2)+1+2*rem(i-1,2)));
subplot(3,2,ind);
cla;
gshade2(abs(tstat), 0.9); hold on;
gcont(tem, clev); hold off;
dc;
title([tit 'Phase = ' num2str((ind-1)*180/nfrm + lag)]);
xlabel(['Cont. Int.: ' xlab]);
end
% subplot(3,2,3);
% ylabel([tit1 ': < ' tit ' CPC1 of 20 deg. T'' >; Years ' ttit]);
else
lag = 0; lg = lag*pi/180;
num = 1; nfrm = 6;
for ind = 1:nfrm;
fig2 = floor(ind/2 - 0.5);
wgt = conj(exp(j * ((ind-1) * pi/nfrm + lg) ));
temtim = real(wgt .* timeseries);
tem3 = temtim' * taux ./ ntim; tem3 = reshape(tem3, nlat, nlon);
tem4 = temtim' * tauy ./ ntim; tem4 = reshape(tem4, nlat, nlon);
% subplot(3,2,ind);
figure(3*fig1 + fig2+1); figure_orient(1);
subplot(2,2,(rem(ind-1,2)+1+2*rem(i-1,2)));
% figure(2*fig1 + fig2+1); figure_orient;
% subplot(2,2,(rem(ind-1,2)+1+2*rem(i-1,2)));
XAX = thin(lon3, 2); YAX = thin(lat3, 2);
tem3 = thin(tem3, 2); tem4 = thin(tem4, 2);
[h, a, c] = gquiv(tem3, tem4, 1.2); dc;
title([tit 'Phase = ' num2str((ind-1)*180/nfrm + lag)]);
xlabel(['Max = ' num2str(round(c(2)*100)/100) ' dynes cm^-^1'])
end
% subplot(3,2,3);
% ylabel([tit1 ': < ' tit ' CPC1 of 20 deg. T'' >; Years ' ttit]);
end
end
%for i = 9:-1:1; figure(i); figure_orient(1); end;
cd /home/disk/tao/dvimont/matlab/CSIRO/Plots4/yr101-550
if 0;
for i = 1:9;
figure(i); figure_orient(1);
pcommand = ['print -dps2 fig_' num2str(i) '_8yr_HP_3dCEOF.ps'];
eval(pcommand);
end
end