C++*******************************************************************
C
C  AVERG.F
C
C **********************************************************************
C=* FROM: SPIDER - MODULAR IMAGE PROCESSING SYSTEM.   AUTHOR: J.FRANK  *
C=* Copyright (C) 1985-2005  Health Research Inc.                      *
C=*                                                                    *
C=* HEALTH RESEARCH INCORPORATED (HRI),                                *   
C=* ONE UNIVERSITY PLACE, RENSSELAER, NY 12144-3455.                   *
C=*                                                                    *
C=* Email:  spider@wadsworth.org                                       *
C=*                                                                    *
C=* This program is free software; you can redistribute it and/or      *
C=* modify it under the terms of the GNU General Public License as     *
C=* published by the Free Software Foundation; either version 2 of the *
C=* License, or (at your option) any later version.                    *
C=*                                                                    *
C=* This program is distributed in the hope that it will be useful,    *
C=* but WITHOUT ANY WARRANTY; without even the implied warranty of     *
C=* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU  *
C=* General Public License for more details.                           *
C=*                                                                    *
C=* You should have received a copy of the GNU General Public License  *
C=* along with this program; if not, write to the                      *
C=* Free Software Foundation, Inc.,                                    *
C=* 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.      *
C=*                                                                    *
C **********************************************************************
C
C      AVERG(LUN1,XAVG,J,K,NSAM,NROW,DS,DR,MOVWAY)
C      
C  XAVG: RECURSIVE LOCAL AVERAGE OF PICTURE AREA
C  X( ): DR PICTURE LINES STORED  X[J,K]<>X((MOD(J-1,DR))*NSAM+K)
C  J,K PICTURE COORD
C  NSAM,NROW: PICTURE SIZE
C  DS,DR: LOCAL BOX  SIZE
C  MOVWAY= 1: READ IN DR LINE AND COMPUTE 1ST AVERAGE VALUE FROM SCRATCH
C             IF J,K IS IN A BORDER AREA COMPUTE XAVG FOR NEAREST
C             ACTIVE POINT
C        = 2: MOVE RT ONE PIXEL
C        = 3: MOVE LEFT ONE PIXEL
C        = 4: MOVE DOWN ONE PIXEL, READ IN NEXT LINE
C
C IMAGE_PROCESSING_ROUTINE
C
C        1         2         3         4         5         6         7
C23456789012345678901234567890123456789012345678901234567890123456789012
C--*******************************************************************

        SUBROUTINE AVERG(LUN1,XAVG,J,K,NSAM,NROW,DS,DR,MOVWAY)



        INTEGER DS,DR,HDS,HDR,EDGFST,EDGLST,EDGLFT,EDGRT
        COMMON /UNITS/LUN,NIN,NOUT,NECHO,IFOUND,NPROC,NDAT
        DIMENSION X(1)
        COMMON X

        SAVE

        IF (MOVWAY .EQ. 1) THEN       
C          INITIALIZE PARM AND XAVG
           HDR=DR/2
           HDS=DS/2
           EDGFST=HDR+1
           EDGLST=NROW-HDR
           EDGLFT=HDS+1
           EDGRT=NSAM-HDS

C          CHECK IF IN BORDER
           XAVG=0.
           JT=J
           KT=K
           IF(J.LT.EDGFST)JT=EDGFST
           IF(J.GT.EDGLST)JT=EDGLST
           IF(K.LT.EDGLFT)KT=EDGLFT
           IF(K.GT.EDGRT)KT=EDGRT

C          READ IN DR LINES AND COMPUTE AVERAGE FROM SCRATCH
           DO  JJ=JT-HDR,JT+HDR
              LSTART=(MOD(JJ-1,DR))*NSAM
              CALL REDLIN(LUN1,X(LSTART+1),NSAM,JJ)
              DO  KK=KT-HDS,KT+HDS
                 XAVG=XAVG+X(LSTART+KK)
              ENDDO
           ENDDO

           XAVG=XAVG/FLOAT(DR*DS)
           RETURN

        ELSEIF (MOVWAY .EQ. 2) THEN
C          MOVE RT ONE PIXEL, UPDATE XAVG
           JT=J
           IF (K.EQ.EDGLFT.OR.K.LT.EDGLFT.OR.K.GT.EDGRT)RETURN
           IF (J.LT.EDGFST)JT=EDGFST
           IF (J.GT.EDGLST)JT=EDGLST
           CORECT=0.
           DO  JJ=JT-HDR,JT+HDR
              LSTART=(MOD(JJ-1,DR))*NSAM
              CORECT=CORECT+X(LSTART+K+HDS)-X(LSTART+K-HDS-1)
           ENDDO

           XAVG=CORECT/FLOAT(DR*DS)+XAVG
           RETURN

        ELSEIF (MOVWAY .EQ. 3) THEN
C          MOVE LEFT ONE PIXEL UPDATE XAVG
           JT=J
           IF(K.LT.EDGLFT.OR.K.GT.EDGRT.OR.K.EQ.EDGRT)RETURN
           IF(J.LT.EDGFST)JT=EDGFST
           IF(J.GT.EDGLST)JT=EDGLST
           CORECT=0.
           DO  JJ=JT-HDR,JT+HDR
              LSTART=(MOD(JJ-1,DR))*NSAM
              CORECT=CORECT+X(LSTART+K-HDS)-X(LSTART+K+HDS+1)
           ENDDO
           XAVG=CORECT/FLOAT(DR*DS)+XAVG
           RETURN

        ELSEIF (MOVWAY .EQ. 4) THEN
C          J INCREMENTED ONE
           IF (J.LE.EDGFST.OR.J.GT.EDGLST)RETURN
           KT=K
           IF (K.LT.EDGLFT)KT=EDGLFT
           IF (K.GT.EDGRT)KT=EDGRT
           CORECT=0.
           JOLD=J-HDR-1
           LSTART=(MOD(JOLD-1,DR))*NSAM
           DO  KK=KT-HDS,KT+HDS
              CORECT=CORECT-X(LSTART+KK)
           ENDDO
           CALL REDLIN(LUN1,X(LSTART+1),NSAM,J+HDR)
           DO  KK=KT-HDS,KT+HDS
              CORECT=CORECT+X(LSTART+KK)
           ENDDO
           XAVG=XAVG+CORECT/FLOAT(DS*DR)

           RETURN
        ENDIF
        END
