% routine to simulate tide given freq, amp, and phase
% and, in phase 2, test for coh&ph calculations of 2 series

% define amplitude
ua=40;
va=70;
period=12.42;

% generate 29 days (696 hours) of a sine function
k=0;
hrph=1;
ang=0;
for i=hrph:695+hrph;
  k=k+1;
  ph(i)=2*pi*(1/period)*i;
  hr(k)=i/24;
  u(k)=ua*cos(ph(i))+rand*20;
  v(k)=va*sin(ph(i))+rand*10;

% add another freq to this
  ph2(i)=2*pi*(1/(period*2))*i;
  u(k)=u(k)+80*cos(ph2(i))+rand*20;
  v(k)=v(k)+70*sin(ph2(i))+rand*10;


% rot axis
  [ur(k),vr(k)]=rot(ang,u(k),v(k));
end;


plot(hr,vr);
axis([0 5 -100 100])

% save a file
yy=[ur; vr]';
fploopy(yy,'uv.dat','o')

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  PHASE TWO %%%%%%%%%%%%%
%  Testing coh&phase calculation

% generate another series

% define amplitude
ua=20;
va=70;
period=24.84;


% generate 29 days (696 hours) of a sine function
k=0;
hrph=19;
ang=0;
for i=hrph:695+hrph;
  k=k+1;
  ph(i)=2*pi*(1/period)*i;
  hr(k)=i/24;
  u(k)=ua*cos(ph(i))+rand*20;
  v(k)=va*sin(ph(i))+rand*10;

% add another freq to this
%  u(k)=u(k)+ua/8*cos(ph(i)/2)+rand*20;
%  v(k)=v(k)+va/8*sin(ph(i)/2)+rand*10;

% rot axis
  [ur2(k),vr2(k)]=rot(ang,u(k),v(k));
end;

y=[yy ur2' vr2'];
t=hr;

plot(hr,ur,hr,ur2);
axis([0 5 -100 100])

foi