%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                     GV2
%          A Gaussian wave-packet.
%
% SECOND PROGRAM TO STUDY GROUP VELOCITY.
%
% Calculate the envelope of the solution
% using the Hilbert transform technique.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

XLEN = 40E6; XNUM = 1001;
x=linspace(0,XLEN,XNUM);
TDAYS = 40; TDAYS = 20; TDAYS = 12; TDAYS = 20;
TLEN=TDAYS*24*3600; TNUM = 501;
t=linspace(0,TLEN,TNUM);

ubar = 2.5; ubar = 5.0;   
ubar = 12.5;
       %  ubar = 0; %  <<<<<<<<<<<<<<<<<<<<<<<

m0 = 24;  m0 = 12;    m0=20;   m0=12;
k0 = 2*pi*m0/XLEN;

beta=1.6E-11;
omega0 = k0*ubar - beta/k0;
c0 = omega0/k0;
cg = ubar + beta/k0^2;
sigma = 2.5E+6;    %   5*(2*pi/k0);  sigma = 0.5E+6 ;
x0 = XLEN/5;

w0 = sin(k0*x);
w = sin(k0*x).*exp(-(x-x0).^2/(2*sigma^2));
XU = ubar*t(1);
axis([0,XLEN,-1,+1]);
axis('off')
subplot(2,1,1);
plot(x,w0,'r','LineWidth',1); 
title('SINGLE WAVE');
hold on; plot(XU,0,'--rs','LineWidth',1,...
                'MarkerEdgeColor','k',...
                'MarkerFaceColor','r',...
                'MarkerSize',10);  hold off
subplot(2,1,2); 
plot(x,w,'b','LineWidth',1); 
title('WAVE PACKET');
hold on; plot(XU,0,'--rs','LineWidth',1,...
                'MarkerEdgeColor','k',...
                'MarkerFaceColor','r',...
                'MarkerSize',10);  % hold off
pause

% Get the Envelope
XH = (XNUM-1)/2;
mask = zeros(size(w)); mask(1:XH) = 1;
%
WT = fft(w); WT = WT.*mask;
wenv = 2*abs(ifft(WT));
plot(x,wenv,':k'); plot(x,-wenv,':k'); 
axis([0,XLEN,-1,+1]); pause; hold off

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

for nt=1:TNUM
   tt = t(nt);
   w0 = sin(k0*x-omega0*tt);
   w = sin(k0*x-omega0*tt).*exp(-((x-x0)-cg*tt).^2/(2*sigma^2));
   XU = ubar*tt; XU = rem(XU,XLEN);
   subplot(2,1,1); 
   plot(x,w0,'r','LineWidth',1); 
   title('SINGLE WAVE');
     hold on; plot(XU,0,'--rs','LineWidth',2,...
                'MarkerEdgeColor','k',...
                'MarkerFaceColor','r',...
                'MarkerSize',10);  hold off
   subplot(2,1,2); 
   plot(x,w,'b','LineWidth',1); 
   title('WAVE PACKET');
     hold on; plot(XU,0,'--rs','LineWidth',2,...
                'MarkerEdgeColor','k',...
                'MarkerFaceColor','r',...
                'MarkerSize',10);  % hold off
     %
     WT = fft(w); WT = WT.*mask;
     wenv = 2*abs(ifft(WT));
     plot(x,wenv,':k'); plot(x,-wenv,':k'); 
     axis([0,XLEN,-1,+1]); hold off
     drawnow
end

return   
