function tau_ac = autocorrelation_time(A)
% Estimate normalized fluctuation autocorrelation time for a timeseries.
%
% ARGUMENTS
%   A - timeseries to compute normalized fluctuation autocorrelation function for
%
% RETURNS
%   tau_ac - estimate of autocorrelation time
%
% REFERENCES
%
% [1] Chodera JD et al. JCTC 2007

% Promote timeseries to double.
A = double(A);

% Determine length of timeseries.
T = length(A);

% Compute sample mean.
EA = mean(A);

% Compute fluctuation timeseries.
dA = A - EA;

% Compute t = 0 term.
CAA0 = mean(dA.*dA);

% Initialize tau.
tau_ac = 0.0;

% Set "safe" tau_max.
tau_max = round(T/2);

% Estimate normalized fluctuation correlation function.
step = 1;
tau = 1;
while (tau <= tau_max)
  % Compute normalized fluctuation correlation function.
  CAA = mean(dA(1:T-tau) .* dA(1+tau:T)) / CAA0;

  % Stop accumulating when we hit a negative value.
  if (CAA <= 0.0)
    break;
  end
  
  % Accumulate contribution to g.
  tau_ac = tau_ac + step * CAA * (1.0 - tau/T);
  
  % Increment step size.
  tau = tau + step;
  step = step + 1;
end

% Ensure tau_ac > 0.
if (tau_ac < 0.0)
  tau_ac = 0.0;
end

return;

