load Orbit91.txt
figure

time = Orbit91(:,1);
ecc = Orbit91(:,2);
omega = Orbit91(:,3);
obl=Orbit91(:,4);
prec = Orbit91(:,5);

[f,P,Pred_norm,Pwhite_norm,P95,P95chi2]=myspectrum(time,ecc);
subplot(2,2,1);
loglog(f,P,'k'); hold on;
loglog(f,Pred_norm,'r');
loglog(f,Pwhite_norm,'b');
title('Eccentricity')

[f,P,Pred_norm,Pwhite_norm,P95,P95chi2]=myspectrum(time,omega);
subplot(2,2,2);
loglog(f,P,'k'); hold on;
loglog(f,Pred_norm,'r');
loglog(f,Pwhite_norm,'b');
title('Omega');

[f,P,Pred_norm,Pwhite_norm,P95,P95chi2]=myspectrum(time,prec);
subplot(2,2,3);
loglog(f,P,'k'); hold on;
loglog(f,Pred_norm,'r');
loglog(f,Pwhite_norm,'b');
title('Precession');

[f,P,Pred_norm,Pwhite_norm,P95,P95chi2]=myspectrum(time,obl);
subplot(2,2,4);
loglog(f,P,'k'); hold on;
loglog(f,Pred_norm,'r');
loglog(f,Pwhite_norm,'b');
title('Obliquity');




