function [y,n_firsty]=TND_conv(x, k_firstx, h, k_firsth, plots,e_zeros);
%
% FUNCTION [y,n_firsty]=TND_conv(x,k_firstx,h,k_firsth,plots,e_zeros);
%
% A function to compute convolution using the graphical method,
% with the option of incremental graphical output
%
% VERSION:  0.9, September 2001.
%
% INPUTS:
%
% x:        a vector containing the input sequence
% k_firstx: time index corresponding to first element of x
% h:        a vector containing the impulse response 
% k_firsth: time index corresponding to first element of h
%
% plots:    boolean variable. Set to 1 to produce plots; 0 to avoid them
% e_zeros:  Number of extra zeros padded to x to make the plots look nice
%           I suggest you set this to be between 0 and 5, inclusive.
%
% OUTPUTS:
%
% y:        a vector containing the output sequence 
% n_firsty: the time index corresponding to the first element of y.      
%
% The vector y contains the output sequence only for time indices
% n from k_firstx+k_firsth to k_firstx+length(x)-1 + k_firsth+length(h)-1.
% If the input and the impulse response are zero outside the elements
% given above, then the output is zero outside the range of time indices
% covered by the vector y.
%   

%%%%%%%%%%%%%%%%% Legal Statements %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Copyright (C) September 2001.
% Timothy N. Davidson,
% Dept. Elec. & Comp. Engineering,
% McMaster University,
% Hamilton, Ontario, Canada.
% davidson@mcmaster.ca

% This program is free software; you can redistribute it and/or
% modify it under the terms of the GNU General Public License
% as published by the Free Software Foundation; either version 2
% of the License, or (at your option) any later version.

% This program is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
% GNU General Public License for more details.

% A copy of the GNU General Public License is available on the web at:
% http://www.gnu.org/copyleft/gpl.html; or by writing to the Free Software
% Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.

 

incr_plots=plots;
extra_zeros=e_zeros;

Lx=length(x);
Lh=length(h);

Ly= Lx + Lh - 1;

n_firsty = k_firstx + k_firsth;

x_padded = [zeros(1,extra_zeros),x, zeros(1,extra_zeros),zeros(1,Lh-1)];

k_firstxpadded=k_firstx-extra_zeros;

Ly_padded=Lx +Lh-1 +2*extra_zeros;

n_firstypadded= n_firsty - extra_zeros;

h_reversed=fliplr(h);


for nn=1:Ly_padded,

 n=n_firstypadded+nn-1;

 no_elements_h=min(Lh,n -k_firsth-k_firstxpadded+1);
 no_zeros_padded_left=max(0,n- k_firsth-Lh-k_firstxpadded+1);
 no_zeros_padded_right=Ly_padded-no_elements_h-no_zeros_padded_left;

 h_rev_overlapped=[zeros(1, no_zeros_padded_left),h_reversed(Lh-no_elements_h+1:Lh),zeros(1,no_zeros_padded_right)];
  
 w= x_padded.*h_rev_overlapped;

 y_padded(nn)=sum(w);


 if incr_plots,

  kleft=min(k_firstxpadded -(k_firsth+Lh-1), k_firstxpadded);
  kright=max(k_firstxpadded+Ly_padded-1,k_firstxpadded +length(x_padded)-1);

  subplot(4,1,1);
  
  stem(k_firstx:k_firstx+Lx-1,x,'r');

  axis([kleft,kright,min([min(x)-.5,0,1.25*min(x)]),max([max(x)+.5,1,1.25*max(x)])]); 
  ylabel('x[k]');
  xlabel('k');

  subplot(4,1,2);  
  stem(n - k_firsth-Lh+1:n-k_firsth,h_reversed,'b');
  axis([kleft,kright,min([min(h)-.5,0,1.25*min(h)]),max([max(h)+.5,1,1.25*max(h)])]); 
  ylabel_string=['h[n-k], n = ',int2str(n)];
  ylabel(ylabel_string);
  xlabel('k');

  subplot(4,1,3);
  stem(k_firstxpadded:k_firstxpadded+Ly_padded-1,w,'m');
  axis([kleft,kright,min([min(w)-.5,0,1.25*min(w)]),max([max(w)+.5,1,1.25*max(w)])]); 
  ylabel_string=['w_n[k], n = ',int2str(n)];
  ylabel(ylabel_string);
  xlabel('k');

  subplot(4,1,4);
  stem(n_firstypadded:n_firstypadded+length(y_padded)-1,y_padded,'k');
  axis([n_firstypadded,n_firstypadded+Ly_padded-1,min([min(y_padded)-.5,0,1.25*min(y_padded)]),max([max(y_padded)+.5,1,1.25*max(y_padded)])]);
  ylabel('y[n]');
  xlabel('n');

  if nn< Ly_padded,
    disp('Press return to continue');
    pause;
  end;

 end; %if incr_plots

end; % nn

y=y_padded(extra_zeros+1:extra_zeros+Ly);
