;+
;
; @file_comments
; 
;
; @categories
; Statistics
;
; @param X {in}{required}
; An 2 dimension Array [nx,ny]
;
; @param LAG {in}{required}
; 2-element vector, in the intervals [-(nx-2), (nx-2)],[-(ny-2), (ny-2)],
; of type integer that specifies the absolute distance(s) between 
; indexed elements of X.
;
; @keyword ZERO2NAN
;
; @keyword DOUBLE
; If set to a non-zero value, computations are done in double precision arithmetic.
;
; @history
; 28/2/2000 Sebastien Masson (smasson\@lodyc.jussieu.fr)
; Based on the A_CORRELATE procedure of IDL
;
; @version
; $Id$
;
;-
;
FUNCTION auto_cov2d, x, lag, DOUBLE = double, ZERO2NAN = zero2nan
;
  compile_opt idl2, strictarrsubs
;
   XDim = SIZE(X, /dimensions)
   nx = XDim[0]
   ny = XDim[1]
;Sample autocovariance function
   Xmean = TOTAL(X, Double = Double) / (1.*nx*ny)
;
   res = TOTAL( (X[0:nx-1-lag[0], 0:ny-1-lag[1]] - Xmean) * $
                (X[lag[0]:nx-1, lag[1]:ny-1] - Xmean) $
                , Double = Double )
   if keyword_set(zero2nan) AND res EQ 0 then res = !values.f_nan
   RETURN, res

END
;+
;
; @file_comments
; This function computes the autocorrelation Px(K,L) or
; autocovariance Rx(K,L) of a sample population X[nx,ny] as a
; function of the lag (K,L).
;
; @categories
; Statistics
;
; @param X {in}{required}
; An 2 dimension Array [nx,ny]
;
; @param LAG {in}{required}
; 2-element vector, in the intervals [-(nx-2), (nx-2)],[-(ny-2), (ny-2)],
; of type integer that specifies the absolute distance(s) between 
; indexed elements of X.
;
; @keyword COVARIANCE
; If set to a non-zero value, the sample autocovariance is computed.
;
; @keyword DOUBLE
; If set to a non-zero value, computations are done in double precision arithmetic.
;
; @history
; 28/2/2000 Sebastien Masson (smasson\@lodyc.jussieu.fr)
; Based on the A_CORRELATE procedure of IDL
;
; @version
; $Id$
;
;-
;
FUNCTION a_correlate2d, x, lag, COVARIANCE = covariance, DOUBLE = double
;
  compile_opt idl2, strictarrsubs
;

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

   ON_ERROR, 2

   XDim = SIZE(X, /dimensions)
   XNDim = SIZE(X, /n_dimensions)
   nx = XDim[0]
   ny = XDim[1]
   if XNDim NE 2 then $
    ras = report("X array must contain 2 dimensions.")
;Check length.
   if nx lt 2 then $
    ras = report("first dimension of X array must contain 2 or more elements.")
   if ny lt 2 then $
    ras = report("second dimension of X array must contain 2 or more elements.")
   if n_elements(Lag) NE 2 THEN $
    ras = report("Lag array must contain 2 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 = (SIZE(X, /type) eq 5)

   if KEYWORD_SET(Covariance) eq 0 then begin ;Compute Autocorrelation.
      Auto = Auto_Cov2d(X, ABS(Lag), Double = Double) / $
          Auto_Cov2d(X, [0L, 0L], Double = Double, /zero2nan)
   endif else begin             ;Compute Autocovariance.
      Auto = Auto_Cov2d(X, ABS(Lag), Double = Double) / n_elements(X) 
   endelse

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

END