Global Index (short | long) | Local contents | Local Index (short | long)
rho = (corr(pcs3(:,1), pcs3(:,1), 1));% + sqrt(corr(ctann, ctann, 2))) / 2 ;% + ...
| This script calls | |
|---|---|
clear
ct = getnc('temp', [180 270 -6 6], 1, 101:1000);
ct = squeeze(mean(mean(ct, 2), 3));
ct = detrend(ct); ct = ct./std(ct);
pcs1 = ct(1:450);
pcs2 = ct(451:900);
pcs3 = ct;
cd /home/disk/tao/dvimont/matlab/CSIRO/Thesis/Data
load LP9_detrend_L1-7_EOF_yr101-1000.mat; pcs1 = pcs;
load HP10_detrend_L1-7_EOF_yr101-1000.mat; pcs2 = pcs;
load RAW_detrend_L1-7_EOF_yr101-1000.mat; pcs3 = pcs;
num = 1;
for i = 1:3;
if i == 2;
ctann = pcs3(1:450,num);
tit = ['Years 101 to 550'];
elseif i == 3;
ctann = pcs3(451:900,num);
tit = ['Years 551 to 1000'];
elseif i == 1;
ctann = pcs3(:,num);
tit = ['Years 101 to 1000'];
end
nx = length(ctann);
nfft = 128;
noverlap = 0.5*nfft;
[p, f] = spectrum(ctann, nfft, noverlap);
n2 = nfft/2;
f2 = 0.5 * (0:n2)/n2;
rho = (corr(ctann, ctann, 1));% + sqrt(corr(ctann, ctann, 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));
if nfft == 128;
if i ~= 1;
rner1 = rn * 2.64; rner5 = rn * 2.01;
else
rner1 = rn * 2.05; rner5 = rn * 1.68;
end
else nfft == 256;
if i ~= 1;
rner1 = rn * 3.55; rner5 = rn * 2.50;
else
rner1 = rn * 2.64; rner5 = rn * 2.01;
end
end
dofx = nx / n2;
figure(1);% figure_orient;
sptalk(4,1,i)
% pos = get(gca, 'Position'); set(gca, 'Position', [0.25 pos(2) 0.5350 pos(4)]);
semilogy(f2, p(:,1), 'b-', f2, rn, 'b-', ...
f2, rner5, 'b--', f2, rner1, 'b-.')
% set(gca, 'YTick', [.1 1 10 100], 'XTick', [0:.05:.5]);
% set(gca, 'YTickLabel', [0.1 1 10])
% axis([f2(3) 0.5 .05 15])
grid
ylabel(tit);
if i == 1;
title(['Power Spectra for PC' num2str(num) ' of 0-270m HC, ' tit]);
end
% xlabel(['Frequency: yr^-^1; Degrees of Freedom: ' ...
% num2str(round(dofx))]);%' ...
% ' Bandwidth: 7.8 x 10^-^3 yr^-^1'])
% l=legend('PC Spectrum', 'Red Noise', '95% Confidence', '99% Confidence');
end
% xlabel(['Frequency: yr^-^1']);%; Degrees of Freedom: ' ...
% num2str(round(dofx))]);%' ...
cd /home/disk/tao/dvimont/matlab/CSIRO/Thesis/Chap1/Plots
for i = 1:2;
yr = 8+i;
[b,a]=butter(9, 2/yr)
[h,w]=freqz(b,a,64);
r=abs(h);
rlow=r.^2;
if i == 1; rlow1 = rlow;
elseif i == 2; rlow2 = 1-rlow;
end
end
f2=w/(2*pi);
figure(2); fo(1);
sptalk(4,1,1)
% pos = get(gca, 'Position'); set(gca, 'Position', [0.25 pos(2) 0.5350 pos(4)]);
h2 = plot(f2, rlow2, '--k', f2, rlow1, '-k');
axis([f2(3) 0.5 -0.1 1.1])
set(gca, 'XTick', [0:.05:.5], 'YTick', 0:.2:1);
grid on
title('Butterworth Filter Response')
legend(h2, '10 Year High Pass', '9 Year Low Pass')
xlabel(['Frequency: yr^-^1']);
cd /home/disk/tao/dvimont/matlab/CSIRO/Thesis/Chap1/Plots
cd /home/disk/tao/dvimont/Thesis/Chap2