C @(#)averow.for	19.1 (ESO-DMD) 02/25/03 14:01:40
C===========================================================================
C Copyright (C) 1995 European Southern Observatory (ESO)
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 
C the 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
C GNU General Public License for more details.
C
C You should have received a copy of the GNU General Public 
C License along with this program; if not, write to the Free 
C Software Foundation, Inc., 675 Massachusetts Ave, Cambridge, 
C MA 02139, USA.
C
C Correspondence concerning ESO-MIDAS should be addressed as follows:
C	Internet e-mail: midas@eso.org
C	Postal address: European Southern Observatory
C			Data Management Division 
C			Karl-Schwarzschild-Strasse 2
C			D 85748 Garching bei Muenchen 
C			GERMANY
C===========================================================================
C
      PROGRAM AVEROW
C --------------------------------------------------------------------
C
C.COPYRIGHT: Copyright (c) 1987 European Southern Observatory,
C                                         all rights reserved
C
C.VERSION: 1.0  ESO-FORTRAN Conversion, AA  13:37 - 5 JAN 1988
C
C.LANGUAGE: F77+ESOext
C
C.AUTHOR: J.D.Ponz
C         900202   KB  take care of options in lower case
C
C.IDENTIFICATION
C
C AVEROW.FOR
C
C.KEYWORDS
C
C  average
C
C.PURPOSE
C
C  produce a 1d image from a 2d by averaging over rows or columns
C
C  COMMAND
C  AVERAGE/ROW    output = input start,end [SUM]
C  AVERAGE/COLUMN output = input start,end [SUM]
C
C.ALGORITHM
C
C  add rows/columns
C
C.INPUT/OUTPUT
C
C KEYWORDS out_a, in_a, p4, p5, action
C
C 001207			last modif
C 
C ------------------------------------------------------------------
C
      IMPLICIT NONE
C
      INTEGER MADRID(1)
      INTEGER NDIM,KUN,KNUL
      INTEGER I1,I2,II1,II2
      INTEGER I,IFIRST,INUM,SUBLO(3)
      INTEGER NAXISA,NAXISB,NN,IMNOA,IMNOB
      INTEGER NPIXA(2),NPIXB(2),STATUS
      INTEGER*8 PNTRA,PNTRB
C
      REAL IAV,LHCUTS(4),FACT
C
      DOUBLE PRECISION STEPA(2),STEPB(2)
      DOUBLE PRECISION STARTA(2),STARTB(2)
C
      CHARACTER*80    FRAMEA, FRAMEB
      CHARACTER*72    IDENT,CUNITA,RANGE
      CHARACTER*72    CUNITB 
      CHARACTER*1     IOP,COMLIN
      CHARACTER       NEWSTR*80
C
      INCLUDE 'MID_INCLUDE:ST_DEF.INC'
C
      COMMON  /VMR/ MADRID
C
      INCLUDE 'MID_INCLUDE:ST_DAT.INC'
C
      DATA   NDIM    /2/
      DATA   IDENT   /' '/,   CUNITA  /' '/
      DATA   CUNITB  /' '/
      DATA   LHCUTS  /0.,0.,0.,0./
C
C  connect to Midas environment
C
      CALL STSPRO('AVEROW')
C
C  read params
C
      CALL STKRDC('OUT_A',1,1,80,IAV,FRAMEB,KUN,KNUL,STATUS)
      CALL STKRDC('IN_A',1,1,80,IAV,FRAMEA,KUN,KNUL,STATUS)
      CALL STKRDC('P4',1,1,72,IAV,RANGE,KUN,KNUL,STATUS)
      CALL STKRDC('P5',1,1,1,IAV,IOP,KUN,KNUL,STATUS)
      CALL STKRDC('ACTION',1,2,1,IAV,COMLIN,KUN,KNUL,STATUS)
C
C  open input frame and read descriptors
C
      CALL STFOPN(FRAMEA,D_R4_FORMAT,F_I_MODE,F_IMA_TYPE,IMNOA,STATUS)
      CALL STDRDI(IMNOA,'NAXIS',1,1,IAV,NAXISA,KUN,KNUL,STATUS)
      CALL STDRDI(IMNOA,'NPIX',1,2,IAV,NPIXA,KUN,KNUL,STATUS)
      CALL STDRDD(IMNOA,'START',1,2,IAV,STARTA,KUN,KNUL,STATUS)
      CALL STDRDD(IMNOA,'STEP',1,2,IAV,STEPA,KUN,KNUL,STATUS)
      CALL STDRDC(IMNOA,'IDENT',1,1,72,IAV,IDENT,KUN,KNUL,STATUS)
      CALL STDRDC(IMNOA,'CUNIT',1,1,72,IAV,CUNITA,KUN,KNUL,STATUS)    
C
C  setup descriptors of result frame
C
      NAXISB = 1
      NPIXB(2) = 1
      STARTB(2) = 1.D0
      STEPB(2) = 0.D0
      IF (RANGE(1:1).EQ.'[') THEN
         II1 = 2
      ELSE
         II1 = 1
      ENDIF
      II2 = INDEX(RANGE,',')
      IF (II2.LT.2) CALL STETER(22,'missing comma in range spec.')
C
C  get range of pixels I1,I2
C
      IF (COMLIN.EQ.'R') THEN
         NPIXB(1) = NPIXA(1)
         STARTB(1) = STARTA(1)
         STEPB(1) = STEPA(1)
         CUNITB(1:16) = CUNITA(1:16)
         CUNITB(17:32) = CUNITA(17:32)
         NEWSTR(1:) = '[@1, '
         NEWSTR(5:) = RANGE(II1:II2-1)
         NN = INDEX(NEWSTR,' ')
         NEWSTR(NN:) = '] '
         CALL EXTCO1(IMNOA,NEWSTR,2,NN,SUBLO,STATUS)
         IF (STATUS.NE.0) CALL STSEPI
         I1 = SUBLO(2)
         NEWSTR(1:) = '[@1, '
         NEWSTR(5:) = RANGE(II2+1:)
         IF (II1.EQ.2) THEN
            NN = INDEX(NEWSTR,']')
         ELSE
            NN = INDEX(NEWSTR,' ')
         ENDIF
         NEWSTR(NN:) = '] '
         CALL EXTCO1(IMNOA,NEWSTR,2,NN,SUBLO,STATUS)
         IF (STATUS.NE.0) CALL STSEPI
         I2 = SUBLO(2)
      ELSE
         NPIXB(1) = NPIXA(2)
         STARTB(1) = STARTA(2)
         STEPB(1)= STEPA(2)
         CUNITB(1:16) = CUNITA(1:16)
         CUNITB(17:32) = CUNITA(33:49)	 
         NEWSTR(1:) = '[ '
         NEWSTR(2:) = RANGE(II1:II2-1)
         NN = INDEX(NEWSTR,' ')
         NEWSTR(NN:) = ',@1] '
         CALL EXTCO1(IMNOA,NEWSTR,2,NN,SUBLO,STATUS)
         IF (STATUS.NE.0) CALL STSEPI
         I1 = SUBLO(1)
         NEWSTR(1:) = '[ '
         NEWSTR(2:) = RANGE(II2+1:)
         IF (II1.EQ.2) THEN
            NN = INDEX(NEWSTR,']')
         ELSE
            NN = INDEX(NEWSTR,' ')
         ENDIF
         NEWSTR(NN:) = ',@1] '
         CALL EXTCO1(IMNOA,NEWSTR,2,NN,SUBLO,STATUS)
         IF (STATUS.NE.0) CALL STSEPI
         I2 = SUBLO(1)
      ENDIF
      IF (I1.GT.I2) CALL STETER(22,'Invalid range')
C
C  map output frame + initialize data
C
      CALL STIPUT(FRAMEB,D_R4_FORMAT,F_O_MODE,F_IMA_TYPE,NAXISB,
     +            NPIXB,STARTB,STEPB,IDENT,CUNITA,PNTRB,IMNOB,STATUS)
      CALL WORK(0,MADRID(PNTRB),NPIXB,0.0,LHCUTS(3))
C
C  do it
C
      IF (COMLIN.EQ.'R') THEN
         CALL STFXMP(NPIXB(1),D_R4_FORMAT,PNTRA,STATUS)
         IFIRST = (I1-1)*NPIXA(1)+1
         DO 100, I=I1,I2
C
C read row after row and sum up
            CALL STFGET(IMNOA,IFIRST,NPIXA(1),IAV,MADRID(PNTRA),STATUS)
            CALL AVER1(MADRID(PNTRB),MADRID(PNTRA),NPIXB(1))
            IFIRST = IFIRST + NPIXA(1)
100      CONTINUE
C 
      ELSE
         INUM = I2-I1+1
         CALL STFXMP(INUM,D_R4_FORMAT,PNTRA,STATUS)
         IFIRST = I1
         DO 200, I=1,NPIXA(2)
C
C read the columns which are to be  summed up
            CALL STFGET(IMNOA,IFIRST,INUM,IAV,MADRID(PNTRA),STATUS)
            CALL AVER2(MADRID(PNTRB),MADRID(PNTRA),
     +           I,INUM,NPIXA(1),NPIXA(2))
            IFIRST = IFIRST + NPIXA(1)
200      CONTINUE
      ENDIF
      FACT = 1.0/(I2-I1+1)
      IF ((IOP.NE.'S').AND.(IOP.NE.'s')) THEN
         CALL WORK(1,MADRID(PNTRB),NPIXB(1),FACT,
     +             LHCUTS(3))                      !average + find minmax
      ELSE
         CALL WORK(2,MADRID(PNTRB),NPIXB(1),FACT,
     +             LHCUTS(3))                      !find minmax only
      ENDIF
      CALL STDWRR(IMNOB,'LHCUTS',LHCUTS,1,4,KUN,STATUS)
      CALL STDWRC(IMNOB,'CUNIT',1,CUNITB,1,72,KUN,STATUS)
      CALL DSCUPT(IMNOA,IMNOB,' ',STATUS)
C
      CALL STSEPI
      END

      SUBROUTINE WORK(FLAG,Y,NPIX,FACT,AUX)
C
      IMPLICIT NONE
C 
      INTEGER  FLAG,NPIX,J
C 
      REAL     Y(NPIX),FACT,AUX(2)
C
      IF (FLAG.EQ.0) THEN
         DO 20, J=1,NPIX
            Y(J) = 0.
20       CONTINUE
      ELSE IF (FLAG.EQ.1) THEN
         AUX(1) = Y(1)*FACT
         AUX(2) = AUX(1)
         DO 30, J=1,NPIX
            Y(J) = Y(J)*FACT
            IF (Y(J).LT.AUX(1)) AUX(1) = Y(J)
            IF (Y(J).GT.AUX(2)) AUX(2) = Y(J)
30       CONTINUE
      ELSE
         AUX(1) = Y(1)
         AUX(2) = AUX(1)
         DO 40, J=1,NPIX
            IF (Y(J).LT.AUX(1)) AUX(1) = Y(J)
            IF (Y(J).GT.AUX(2)) AUX(2) = Y(J)
40       CONTINUE
      ENDIF
C
      RETURN
      END

      SUBROUTINE AVER1(Y,Y1,NPIX1)
C
C sum over rows
C
      IMPLICIT NONE
      INTEGER  NPIX1,J
      REAL     Y(NPIX1), Y1(NPIX1)
C
      DO 20, J=1,NPIX1
         Y(J) = Y(J) + Y1(J)
20    CONTINUE
C
      RETURN
      END

      SUBROUTINE AVER2(Y,Y1,ICOL,INUM,NPXA,NPXB)
C
C sum over columns
C
      IMPLICIT NONE
      INTEGER  ICOL,J,INUM,NPXA,NPXB
      REAL     Y(NPXA), Y1(NPXB)
C
      DO 20, J=1,INUM
         Y(ICOL) = Y(ICOL)+Y1(J)
20    CONTINUE
C    
      RETURN
      END

