function autoco = autocorrelation(s_train, dt, lags, zero_mean_flag)

% Constructs the autocorrelation histogram of a spike train
%   If zero_mean_flag is set to 1, the output is the "auto-covariogram"     

%   NEUR 1680
%   Last updated 02/08/2012

s_times = find(s_train);            % spike times in bin units
n = length(s_times);                % number of spikes in record
N = length(s_train);                % number of bins in record
T = N*dt;                           % length of record in sec
autoco = zeros(1,length(lags));

m = 0;
for lag = lags
    if lag ~= 0
        m = m + 1;
        lst = s_times + lag;                    % lagged spike times
        lst = lst( lst >= 1 & lst <= N );       % discard times outside record
        N_m = sum(s_train(lst));                % number of spike pairs that fall into bin m of histogram
        autoco(m) = N_m;
    end
end

if zero_mean_flag
    autoco = autoco - n^2*dt/T;                 % subtract the mean to get a covariance rather than a correlation
end

autoco = autoco/T;                  % divide by length of record to normalize

