%% convolution theorem
% show 1D pure wave forms and their convolution with a mask 
% illustrating the eigenvalue property

% S. Allen Broughton
 
N = 64;
FD = 0:(N-1);            % non-shifted frequency domain
T = FD/N;                % discrete time points
h1 = [2,-1];             % convolution mask causal part
h2 = [-1];                % convolution mask anti-causal part
h = [h1,zeros(1,N-length(h1)-length(h2)),h2]; % create full length mask
h = h/sum(abs(h(:)));    % normalize
dt = 0.25;               % time between pictures


% iteration parameters
waitforkey = [0,1,2,3,4,5]

% display results 
figure(1) 
clf
for k = 0:(N-1)
   E = exp(2*pi*i*k*T);
   hE = circconv(E,h); 
   rE = real(E);
   iE = imag(E);
   rhE = real(hE);
   ihE = imag(hE);
   fh  = fft(h);
   afh = abs(fh);
   tfh = angle(fh); 
   
   
   subplot(2,3,1);
    plot(T,rE,'r');
    axis([-0.1,1.1,-1.1,1.1]);
    title('original - real part');
    drawnow
    
   subplot(2,3,2);
    plot(T,iE,'g');
    axis([-0.1,1.1,-1.1,1.1]);
    title('original - imaginary part')
   
   subplot(2,3,4);
    plot(T,rhE,'r');
    axis([-0.1,1.1,-1.1,1.1]);
    title('filtered - real part');
    drawnow
    
   subplot(2,3,5);
    plot(T,ihE,'g');
    axis([-0.1,1.1,-1.1,1.1]);
    title('filtered - imaginary part')
    
   subplot(2,3,3);
    hold on 
    plot(FD,afh,'r.');
    plot(k,afh(k+1),'k.'); 
    axis([0,N-1,-1.1,1.1]);
    title('DFT(h) - amplitude');
    hold off
    drawnow
    
   subplot(2,3,6);
    hold on 
    plot(FD,tfh,'g.');
    plot(k,tfh(k+1),'k.'); 
    axis([0,N-1,-3.15,3.15]);
    title('DFT(h) - phase')
    hold off
    drawnow
    
   pauseforkey = sum(find(waitforkey==k)); 
   if pauseforkey > 0
      disp('Press a key to proceed'), pause,
      else pause(dt)
      end
   end