%
% function ACFunction=acf(xSeries,nLags,graphFlag);
%
% Calculates ACF for x, includes Bartlett bounds
% updated 10/4/02 Jirka

function [ACFunction, ACVFunction, bartlettBound]=acf(xSeries,nLags,graphFlag);
if nargin==2; graphFlag=1; end; if nargin==1; graphFlag=1; nLags=20; end;
T=length(xSeries);
for i=1:nLags
    xCurrent=xSeries(1:T-i); xFuture=xSeries(1+i:T);
    CorrMatrix=corrcoef(xCurrent,xFuture)*(T-i)/T; CovMatrix=cov(xCurrent,xFuture)*(T-i)/T; 
    ACFunction(i,1)=CorrMatrix(2,1); ACVFunction(i,1)=CovMatrix(2,1);
end;
bartlettBound=2*ones(nLags,1)/sqrt(T);

if graphFlag==1;
   figure;
    liteGB=[171,219,224]/255;
    h=bar(ACFunction); set(h,'FaceColor',liteGB);
    titleString='Autocorrelation Function'; if ~isempty(inputname(1)); titleString=[titleString ' of ' inputname(1)]; end;
    title(titleString); hold on;
    t=1:nLags; plot(t,bartlettBound,':k'); hold on; plot(t,-bartlettBound,':k'); 
    axis tight; axis 'auto y'; hold off;
end;
