function [Y,t] = gensig(T,N,B,tstep)
% [Y,t] = gensig(T,N,B)
%  generates N low-passed zero-mean unit-variance Gaussian noise waveforms,
%  each T milliseconds long.
%
%  T     = duration of input stimulus I  (in msec)
%  N     = no. of trials
%  B     = bandwidth of I                (in Hz)
%  tstep = sampling interval             (in sec)
%
%  Y     = generated waveform
%  t     = corresponding time vector (in msec)
%
% If the user does not specify a tstep, the program computes the tstep
% automatically

order = 5;                                 % order of the filter
t_s = sqrt(2^(1/(order+1)) - 1)/(2*pi*B);   % signal correlation time

if (nargin < 4)
   tstep = t_s/10;                             % filter sampling interval
end

tf =  (0:tstep:15*t_s);                   % samples for smoothing filter
h_s = tf.^order .* exp(-tf/t_s);          % filter gives sharp cut-off


t = [0:tstep*1000:T]';                    % Output time vector (in msec)
L = length(t);                            % Length of the output vector
rand('seed',sum(100*clock));

I   = randn(L,N);                         % start w/ white noise
Y   = filterC(h_s,I);                     % convolve w/ kernel
Y   = Y - ones(L,1)*mean(Y);              % set to zero mean
Y   = Y./(ones(L,1)*std(Y));              % set to variance 1