%  wavetraj2.m
%  b klinger, 2/20/96
%  Show surface position and some subsurface particles for surface waves.
%  Based on formulation valid for both deep and shallow water waves
%  (see Gill, pp 100-102)
%  based on wavetraj.m
%	Gravitational constant, water depth, amplitude, wavelength, 
%	and (using dispersion relation) frequency.
% Good examples:		deep: H=100, lam=10
%		     kind of shallow: H=2, lam=10
%		        very shallow: H=2, lam=40
% 	MODIFICATIONS
% 2/24/98 matlab 5.0, thickline & marker instead of thinline

g=10;
A=.5;
iwhich=input('Deep, shallow, very shallow (1, 2, 3): ');
if iwhich==1,
  H=100;
  lam=10;
elseif iwhich==2,
  H=2;
  lam=10;
elseif iwhich==3,
  H=2;
  lam=40
end
k=2*pi/lam;
omeg=sqrt(g*k*tanh(k*H));

%	initial particle positions
np=5;
x=[zeros(1,np) 2.5*ones(1,np) -2.5*ones(1,np)];
zp=-2*[0:np-1]/np;
z=[zp zp-2/np zp-2/np];
xx=-5:.25:5;

%	time parameters
ntim=100; ntout=10; dt=.1;

%	do time integration with simple time-stepping.
clf; 
figure(gcf)
for itim=1:ntim
   t=itim*dt;
   h=A*cos(k*xx-omeg*t);
   u=g*(k/omeg)*A*cos(k*x-omeg*t).*cosh(k*(z+H))/cosh(k*H);
   w=g*(k/omeg)*A*sin(k*x-omeg*t).*sinh(k*(z+H))/cosh(k*H);
   x=x+u*dt;
   z=z+w*dt;
   plot(x,z,'.','markersize',24)
   hold on
   plot(xx,h,'b','linewidth',4)
   axis([-5 5 -2 1])
%   thickline(8)
%   marker(6)
   pause(.1)
   hold off
end