C @(#)rfotdiaphr.for	19.1 (ES0-DMD) 02/25/03 13:30:14
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 Massachusetss Ave, Cambridge, 
C MA 02139, USA.
C
C Corresponding 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 DIAFRA    
C+++
C.IDENTIFICATION: RFOTDIAPHR
C.PURPOSE: Do aperture photometry at positions read from a registration table
C.AUTHOR:  R. Buonanno, G. Buscema, C. Corsi, I. Ferraro, G. Iannicola
C          Osservatorio Astronomico di Roma, and
C.         R.H. Warmels, ESO-IPG
C.VERSION: 881116 RB   First installation in MIDAS
C.VERSION: 890918 RHW  Rewritten for MIDAS tables
C.VERSION: 891222 RHW  Input obtained from registration table 
C----
      IMPLICIT       NONE
      INCLUDE        'MID_REL_INCL:RFOTDECL.INC'
C
      INTEGER        NCPAR
      PARAMETER      (NCPAR=12)
C ***
      INTEGER        NAXIS,NPIX(2)
      INTEGER        ICPAR(NCPAR)
      INTEGER        MADRID(1)
      INTEGER        IDNT
      INTEGER*8      IPNTR
      INTEGER        IMF
      INTEGER        TIDOUT
      INTEGER        KUN,KNUL
      INTEGER        NRREG,NCREG,NACREG,NARREG,NSREG
      INTEGER        OUTTYP
      INTEGER        TIDREG
      INTEGER        TINULL
      INTEGER        ICIDN
      INTEGER        IAV, ISTAT
      INTEGER        EC, ED, EL
      INTEGER        NPL, NL
      INTEGER        NCOUT, NROUT
      INTEGER        OUTCOL
      INTEGER        IC, IR, IROW
      INTEGER        IX, IY, IX0, IY0
      INTEGER        LF, LD, I1, I2, I1O
C
      DOUBLE PRECISION BEGIN(2),STEP(2)
      DOUBLE PRECISION TDNULL,TDTRUE,TDFALS
C
      REAL           TAB(NCPAR)
      REAL           RV(5000)
      REAL           TRNULL, TBLSEL
      REAL           FON, AP, S, DS 
C
      CHARACTER*60   FRAME
      CHARACTER*72   IDENT,CUNIT
      CHARACTER*60   REGFIL
      CHARACTER*60   OUTFIL
      CHARACTER*80   STRING
      CHARACTER*16   LABEL(NCPAR),OUTLAB
      CHARACTER*16   UNIT(NCPAR),OUTUNI
      CHARACTER*16   OUTFOR
      CHARACTER*16   FORMR4,FORMI4
C
      LOGICAL        LOG, NUL(NCPAR)
      LOGICAL        SFLAG
C
      INCLUDE        'MID_INCLUDE:TABLES.INC'
      COMMON         /VMR/MADRID
      INCLUDE        'MID_INCLUDE:TABLED.INC'
C ***
      DATA           ICIDN/1/
      DATA           FORMI4/'I6'/
      DATA           FORMR4/'E12.4'/
      DATA     LABEL /'X      ', 'Y      ', 'INT        ', 'LOC_BGD ', 
     2                'MAG1   ', 'MAG2   ', 'MAG3       ', 'MAG_CNV ', 
     3                'SIGMA  ', 'BETA   ', 'SIQ        ', 'CHI_SQ  '/
      DATA     UNIT / 'PIXEL  ', 'PIXEL  ', '           ', '        ',
     2                'MAG.   ', 'MAG.   ', 'MAG.       ', 'MAG.    ',
     3                '       ', '       ', '           ', '        '/
C
C *** Is MIDAS out there?
      CALL STSPRO('DIAFRAGM')
      CALL TBMNUL(TINULL,TRNULL,TDNULL)
      CALL TBMCON(TBLSEL,TDTRUE,TDFALS)
C
C *** get the frame name
      CALL STKRDC('IN_A',1,1,60,IAV,FRAME,KUN,KNUL,ISTAT)
      CALL STIGET(FRAME,D_R4_FORMAT,F_I_MODE,F_IMA_TYPE,2,NAXIS,
     2            NPIX,BEGIN,STEP,IDENT,CUNIT,IPNTR,IMF,ISTAT)
C 
      NPL = NPIX(1)
      NL  = NPIX(2)
C
C *** get the catalogue table
      CALL STKRDC('IN_B',1,1,60,IAV,REGFIL,KUN,KNUL,ISTAT)  !input register name
      IF (ISTAT.NE.0) THEN
         STRING = '*** FATAL: Error opening registration table'
         CALL STTPUT(STRING,ISTAT)
         CALL STSEPI
      ENDIF
C
C *** open the input register table
      CALL STECNT('GET',EC,ED,EL)
      CALL STECNT('PUT',1,0,0)
      CALL TBTOPN(REGFIL,F_IO_MODE,TIDREG,ISTAT)
      IF (ISTAT.NE.0) THEN
         STRING = '*** FATAL: Problems with opening registration'//
     2            ' table; Try again ...'
         CALL STTPUT(STRING,ISTAT)
         CALL STSEPI
      ENDIF
C
      CALL TBIGET(TIDREG,NCREG,NRREG,NSREG,NACREG,NARREG,ISTAT)
      IF (ISTAT.NE.0) THEN
         STRING = '*** FATAL: Problems with getting info for'//
     2            ' registration table; Try again ...'
         CALL STTPUT(STRING,ISTAT)
         CALL STSEPI
      ENDIF
      IF (NRREG.EQ.0) THEN
         STRING = '*** FATAL: No data in the registration table'
         CALL STTPUT(STRING,ISTAT)
         CALL STSEPI
      ENDIF
C
C *** get output table
      CALL STKRDC('OUT_A',1,1,60,IAV,OUTFIL,KUN,KNUL,ISTAT)    ! get output name
C
C *** create the output table
      CALL TBTINI(OUTFIL,0,F_O_MODE,NCOUT,NROUT,TIDOUT,ISTAT)
      IF (ISTAT.NE.0) THEN
         STRING = '*** FATAL:  Problems with opening a new '//
     2            'registration table'
         CALL STTPUT(STRING,ISTAT)
         CALL STSEPI
      ENDIF
      CALL STECNT('PUT',EC,ED,EL)
C
      OUTTYP = D_I4_FORMAT
      OUTFOR = FORMI4
      OUTUNI = '   '
      OUTLAB = 'IDENT'
      CALL TBCINI(TIDOUT,OUTTYP,1,OUTFOR,OUTUNI,
     2                   OUTLAB,OUTCOL,ISTAT)          ! create the ident column
C
      DO 101 IC = 1,NCPAR
         ICPAR(IC) = IC + 1
         OUTTYP    = D_R4_FORMAT
         OUTFOR    = FORMR4
         OUTUNI    = UNIT(IC)
         OUTLAB    = LABEL(IC)
         CALL TBCINI(TIDOUT,OUTTYP,1,OUTFOR,OUTUNI,
     2                      OUTLAB,OUTCOL,ISTAT)       ! create the data columns
  101 CONTINUE
C
C *** radius for integration
      CALL STKRDI('INPUTI',1,1,IAV,IR,KUN,KNUL,ISTAT)   ! get integration radius
C
C *** do the work
      DO 100 IROW=1,NRREG
         CALL TBSGET(TIDREG,IROW,SFLAG,ISTAT)
         IF (SFLAG) THEN
            CALL TBERDI(TIDREG,IROW,ICIDN,IDNT,NUL,ISTAT)
            CALL TBRRDR(TIDREG,IROW,NCPAR,ICPAR,TAB,NUL,ISTAT)
C
            IX  = INT((TAB(1)-BEGIN(1))/STEP(1)+1)
            IY  = INT((TAB(2)-BEGIN(2))/STEP(2)+1)
            FON = TAB(4)
            LOG = (IX+IR.LE.NPL) .AND. (IY+IR.LE.NL)
            IF (IX-IR.GT.0 .AND. IY-IR.GT.0 .AND. LOG) THEN
               LF  = 0
               IX0 = IX-IR
               IY0 = IY-IR
               LD  = IR*2+1
               S   = 0.
               IC  = LD/2+1
               DO 110 I1=IY0,IY0+LD-1
                  CALL REALIN(NPL,NL,I1,IX0,LD,MADRID(IPNTR),RV)
                  I1O = I1-IY0+1
                  DO 111 I2=1,LD
                     DS = (I1O-IC)**2+(I2-IC)**2
                     IF (SQRT(DS).LE.IR) THEN
                        AP = RV(I2)
                        S  = S+AP-FON
                     ENDIF
  111            CONTINUE
  110          CONTINUE
C
               IF (S.GT.0) THEN
                  S = -2.5*ALOG10(S)
               ELSE
                  S = 0.
               END IF
               TAB(8) =  S
C
               CALL TBEWRI(TIDOUT,IROW,ICIDN,IDNT,ISTAT)
               CALL TBRWRR(TIDOUT,IROW,NCPAR,ICPAR,TAB,ISTAT)
            END IF
         ENDIF
100   CONTINUE

      CALL TBSINI(TIDOUT,ISTAT)
      CALL TBTCLO(TIDOUT,ISTAT)
      CALL TBTCLO(TIDREG,ISTAT)
      CALL STSEPI
      END
