
%  waves6
%  Some elementary wave exercises to demonstrate the
%  effects of superposition, interferance, etc., between
%  two sinusoidal waves by way of a Matlab movie. 
%
%  Jim Price, January, 2000

clear
clear memory
close all

%  Set some default graphics things.
set(0,'DefaultTextFontSize',14)
set(0,'DefaultAxesFontSize',14)
set(0,'DefaultAxesLineWidth',1.6)
set(0,'DefaultLineLineWidth',1.4)


figure(1)
clf reset
orient tall

x = -15:0.05:15;   %  the distance domain
dt = 0.1;          %  the time step
nstep = 200;       %  number of time steps to take
M = moviein(nstep);

%  define the wavelengths and wavenumbers of two pure waves
%  the base wave (the second pair is selected with the mouse)
k1 = 5.;
omega1 = 5.;

%  plot the dispersion diagram
subplot(2,1,1)
plot(k1, omega1, 'r*')
hold on; plot([0 10], [0 10], ':b')
axis([4 6 4 6])
axis('square')
grid
xlabel('wavenumber, k')
ylabel('angular freq., \omega')
set(gca,'fontsize',14)
set(gca,'xminortick','on')
set(gca,'yminortick','on')
set(gca,'xtick',[4:0.5:6])
set(gca,'ytick',[4:0.5:6])


grid on

title('Do you believe in group velocity?')

aaa = '  Use the mouse to select the (k, \omega) of the second wave.'

[k2, omega2] = ginput(1)
hold on
plot(k2, omega2, 'g*')


c1 = omega1/k1;
c2 = omega2/k2;
cgc = (omega2 + omega1)/(k2 + k1);
cg = (omega2 - omega1)/(k2 - k1);
deltg = 10; if cg <= 0; deltg = -10.; end


for j = 1:nstep
   
   t = j*dt;
   arg1 = k1*x - omega1*t;
   arg2 = k2*x - omega2*t;
   
   s1 = cos(arg1);
   s2 = cos(arg2);   
   sum = s1 + s2;

sg = 2*cos((k2-k1)*x/2 - (omega2 - omega1)*t/2);

   subplot(2,1,2)  
%plot(x, sum, ':b', x, s1, 'r', x, s2, 'g', x, sg, 'b', x, -sg, 'b'); 
plot(x, s1, 'r'); hold on; grid on
plot(x, s2, 'g')
plot(x, sg, 'b--')
plot(x, -sg, 'b--'); 
plot(x, sum, 'b','linewidth',3); 

hold on
plot(t*c1-10, 0, '*r', t*c2-10, 0, '*g', t*cgc-10, 0, '*b')
hp = text(t*cg-deltg, 0, '\bf g');
hold off

 xlabel('x distance (m)')
 ylabel('amplitude') 
 axis([-15 15 -2 2]);
 text(10,2.5,['Time: ' num2str(t) ' s'])
   
 M(:,j) = getframe;
   
end

%  rerun the movie
%movie(M,2,6);