% function [z,t,A] = gausspulse
% Plots the amplitude and instantaneous frequency of a gaussian
% pulse propagating through a dispersive medium (Eq. (10.96) in 
% the textbook by Rao), with dispersion  Dnu = tau0/vg 
% tau0 is the initial 1/e-with of the pulse amplitude, 
% and vg is group velocity
%
%  Output:
%  z: sampling points in propagation direction
%  t: sampling points in time
%  A: matrix of complex amplitude values

z  =  0:.1:10;  % 101 position samples, z is in units of z0
t  = -1:.1:10;  % 111 samples for initial Gaussian
                % time is scaled to yield unity group velocity
tau0  =  1/pi;  % initial pulse witdth
Dnu   = -tau0;  % dispersion chosen 
[Z,T] =  meshgrid(z,t); 
q  =  Z + j;    % note that z0 is unity
A  =  sqrt(j./q) .* exp(j.*pi.*(T - Z).^2./(Dnu.*q));
f  =  real ((T - Z)./(Dnu.*q)); % instantaneous frequency
A0 = find(abs(A) < 1.e-16);
f(A0) = 0.;
mesh (z, t, abs(A), f);
pause
% make movie
for t = 1:12
   plot (z, abs(A(t*10-9,:)));
   title  ('Gaussian pulse in dispersive medium')
   xlabel ('Propagation distance')
   ylabel ('Amplitude')
   axis ([0 10 0 1])
   M(t) = getframe;
end
movie(M)

   



