Documentation of cospectrum

Global Index (short | long) | Local contents | Local Index (short | long)

Help text

  Note:  a positive phase implies var1 leads var2
         a negative phase implies var2 leads var1
  A one-year lag will be a ~linearly increasing phase difference

Cross-Reference Information

Listing of script cospectrum

clear
cd ~/matlab/CSIRO/Thesis/Data
load slp_eof_npac.mat
tim = 101:1000;
lims = [-0.1 360 -90 90];
var = getflx('psl', lims, tim);
var = detrend(var);
var = getnc('temp', lims, 1, tim);
var = detrend(var);
for i = 1:2
  reg(i,:,:) = regress_eof(var, top, i-2);%rpcs(:,i));
end
[lat, lon, depth, lm] = getll('temp', lims);
default_global; FRAME = [0 360 -90 90];
figure(4); fo;
for i = 1:2
spthes(i)
     gcont(reg(i,:,:), .05);
     dc2(lm);
end
ct = getct;
ct = detrend(ct)./std(detrend(ct));
np = getnc('temp', 
load PDO_detrend_L1_EOF_yr101-1000.mat
load RAW_detrend_L1-7_EOF_yr101-1000.mat
for i = 1:3;
  if i == 1; tim2 = 101:1000; 
  elseif i == 2; tim2 = 101:550;
  elseif i == 3; tim2 = 551:1000; 
  end
tind = find(ismember(tim, tim2));
nfft = 64; noverlap = nfft/2;
[p1, f1] = spectrum(pcs(tind,2), ct(tind), nfft, noverlap);
amp = abs(p1(:,3)); 
phs = atan2(imag(p1(:,3)), real(p1(:,3)));
figure(i); fo(1);
subplot(3,1,1);
  plot(f1./2, p1(:,1), 'r', f1./2, p1(:,2), 'b');
  grid on
subplot(3,1,2);
  plot(f1./2, p1(:,5), 'k');%, f1./2, amp, 'c');
  grid on
  hold on;
    hline(.206 + floor(i/2)*.187, 'k');
  hold off;
  axis([0 0.5 0 1]);
subplot(3,1,3);
  plot(f1./2, 180*phs/pi, '.k');
  grid on
  axis([0 0.5 -180 180]);
  set(gca, 'YTick', -180:90:180)
end

%  look at cospectrum of ML upper level HC, and WSTP lower level HC

clear

tim = 101:1000;
top = getheat([175 205 25 40], 1:3, tim);
bot = getheat([120 165 7.5 22.5], 4:7, tim);
tim = 101:1000;
top = getheat([180 270 -6 6], 1:3, tim);
bot = getheat([120 165 7.5 22.5], 1:3, tim);

top = squeeze(mean2(mean2(shiftdim(top, 1))));
bot = squeeze(mean2(mean2(shiftdim(bot, 1))));

top = detrend(top)./std(detrend(top));
bot = detrend(bot)./std(detrend(bot));

nfft = 128; noverlap = nfft/2;
[p1, f1] = spectrum(top, bot, nfft, noverlap);

amp = abs(p1(:,3)); 
phs = atan2(imag(p1(:,3)), real(p1(:,3)));

figure(2); fo(1);
subplot(3,1,1);
  plot(f1./2, p1(:,1), 'r', f1./2, p1(:,2), 'b');
  grid on
subplot(3,1,2);
  plot(f1./2, p1(:,5), 'k');%, f1./2, amp, 'c');
  grid on
  hold on;
    hline(.206 + floor(i/2)*.187, 'k');
  hold off;
  axis([0 0.5 0 1]);
subplot(3,1,3);
  plot(f1./2, 180*phs/pi, '.k');
  grid on
  axis([0 0.5 -180 180]);
  set(gca, 'YTick', -180:90:180)



[b, a] = butter(9, 2/9);
corr(filtfilt(b, a, top), filtfilt(b, a, bot), 0);