; $Id: a_correlate.pro,v 1.13 2001/01/15 22:27:59 scottm Exp $
; Copyright (c) 1995-2001, Research Systems, Inc.  All rights reserved.
;       Unauthorized reproduction prohibited.
;+
; NAME:
;       A_CORRELATE
;
; PURPOSE:
;       This function computes the autocorrelation Px(L) or autocovariance
;       Rx(L) of a sample population X as a function of the lag (L).
;
; CATEGORY:
;       Statistics.
;
; CALLING SEQUENCE:
;       Result = A_correlate(X, Lag)
;
; INPUTS:
;       X:    An n-element vector of type integer, float or double.
;
;     LAG:    A scalar or n-element vector, in the interval [-(n-2), (n-2)],
;             of type integer that specifies the absolute distance(s) between 
;             indexed elements of X.
;
; KEYWORD PARAMETERS:
;       COVARIANCE:    If set to a non-zero value, the sample autocovariance
;                      is computed.
;
;       DOUBLE:        If set to a non-zero value, computations are done in
;                      double precision arithmetic.
;
; EXAMPLE
;       Define an n-element sample population.
;         x = [3.73, 3.67, 3.77, 3.83, 4.67, 5.87, 6.70, 6.97, 6.40, 5.57]
;
;       Compute the autocorrelation of X for LAG = -3, 0, 1, 3, 4, 8
;         lag = [-3, 0, 1, 3, 4, 8]
;         result = a_correlate(x, lag)
;
;       The result should be:
;         [0.0146185, 1.00000, 0.810879, 0.0146185, -0.325279, -0.151684]
;
; PROCEDURE:
;       See computational formula published in IDL manual.
;
; REFERENCE:
;       INTRODUCTION TO STATISTICAL TIME SERIES
;       Wayne A. Fuller
;       ISBN 0-471-28715-6
;
; MODIFICATION HISTORY:
;       Written by:  GGS, RSI, October 1994
;       Modified:    GGS, RSI, August 1995
;                    Corrected a condition which excluded the last term of the
;                    time-series.
;       Modified:    GGS, RSI, April 1996
;                    Simplified AUTO_COV function. Added DOUBLE keyword.
;                    Modified keyword checking and use of double precision.
;-

FUNCTION Auto_Cov, X, M, nX, Double = Double, XMEAN=Xmean
  ;Sample autocovariance function.

if n_elements(Xmean) eq 0 then $ ;Compute mean if undefined.
  Xmean = TOTAL(X, Double = Double) / nX
RETURN, TOTAL((X[0:nX - M - 1L] - Xmean) * (X[M:nX - 1L] - Xmean), $
              Double = Double)
END

FUNCTION A_Correlate, X, Lag, Covariance = Covariance, Double = Double

  ;Compute the sample-autocorrelation or autocovariance of (Xt, Xt+l)
  ;as a function of the lag (l).

  ON_ERROR, 2

  TypeX = SIZE(X)
  nX = TypeX[TypeX[0]+2]

  ;Check length.
  if nX lt 2 then $
    MESSAGE, "X array must contain 2 or more elements."
 
  ;If the DOUBLE keyword is not set then the internal precision and
  ;result are identical to the type of input.
  if N_ELEMENTS(Double) eq 0 then $
    Double = (TypeX[TypeX[0]+1] eq 5)

  nLag = N_ELEMENTS(Lag)

  if nLag eq 1 then Lag = [Lag] ;Create a 1-element vector.

  if Double eq 0 then Auto = FLTARR(nLag) else Auto = DBLARR(nLag)

  if KEYWORD_SET(Covariance) eq 0 then begin ;Compute Autocorrelation.
      for k = 0, nLag-1 do $
        Auto[k] = Auto_Cov(X, ABS(Lag[k]), nX, Double = Double, XMEAN=xm)
      Auto = Auto / Auto_Cov(X, 0L, nX, Double = Double, XMEAN=xm)
  endif else begin              ;Compute Autocovariance.
      for k = 0, nLag-1 do $ 
        Auto[k] = Auto_Cov(X, ABS(Lag[k]), nX, Double = Double, XMEAN=xm) / nX
  endelse

  if Double eq 0 then RETURN, FLOAT(Auto) else $
  RETURN, Auto

END