;------------------------------------------------------------------
;+
; NAME:
;	ROTLATAVG
; PURPOSE:
;
;       This function takes a grid of rotation frequencies found
;       using Louis Strous' gridmatch function and produces an
;       average value for each latitude.
;
;       Points with a rotation frequency of exactly 0.0 are not
;       included since they most likely points from a part of the
;       map with no data or there was no magnetic feature present
;       to give a valid cross-correlation.
;
;       Points which have a rotation frequency a factor of three
;       different from the initial mean are not included since
;       they are most likely either points next to a no-data area
;       on the map (the LCT algorithm is influenced by such areas).
;
; CATEGORY:
;	POLAR
; CALLING SEQUENCE:
;	latavg = ROTLATAVG(gridimg)
; INPUTS:
;       gridimg = array output from gridmatch.  This is a 2 x nx x ny
;                 array.
; KEYWORD PARAMETERS:
; OUTPUTS:
;       A vector containing latitude averages.
; COMMON BLOCKS:
; NOTES:
; MODIFICATION HISTORY:
;       J. Varsik,      20 May  1999
;-
;-------------------------------------------------------------------
FUNCTION ROTLATAVG,gridimg,HELP=help

    ; Display IDL header if help is required.
    IF (KEYWORD_SET(help)) THEN BEGIN
      GET_IDLHDR,'rotlatavg.pro'
      GOTO,finishup
    ENDIF

    gsz = SIZE(gridimg)
    IF (gsz(0) NE 3) THEN BEGIN
        PRINT,"Wrong size grid"
        GOTO,finishup
    ENDIF

    IF (gsz(1) NE 2) THEN BEGIN
        PRINT,"Does not appear to be a LCT grid"
        GOTO,finishup
    ENDIF

    nx = gsz(2)
    ny = gsz(3)
    grid = reform(gridimg(0,*,*))
;   initial guess
    initg = TOTAL(grid,1)/FLOAT(nx)
    latavg = fltarr(ny)

;   Now check on each row
    FOR i = 0, ny - 1 DO BEGIN
        imean = initg(i)
        ipts = 0
        FOR j = 0, nx - 1 DO BEGIN
            IF (grid(j,i) NE 0.0) THEN BEGIN
                IF (grid(j,i) LT imean/3.0) THEN BEGIN
                    grid(j,i) = 0.0
                ENDIF ELSE BEGIN
                    ipts = ipts + 1
                ENDELSE
            ENDIF
        ENDFOR
        IF (ipts NE 0) THEN BEGIN
            latavg(i) = TOTAL(reform(grid(*,i)),1)/FLOAT(ipts)
        ENDIF ELSE BEGIN
            latavg(i) = 0.0
        ENDELSE
    ENDFOR

    finishup:

    RETURN,latavg
end