C @(#)raddet.for	17.1.1.1 (ESO-DMD) 01/25/02 17:15:44
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
C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C.IDENTIFICATION
C  Subroutine RADDET          version 1.0          870626
C  A. Kruszewski              ESO Garching
C.PURPOSE
C  Find active radius of detected object.
C.ALGORYTHM
C  Goes out from the object center in four directions looking when
C  intensity drop down below detection threshold. Active radius is
C  adopted as the second least distant drop-down point.
C.INPUT/OUTPUT
C  Input parameters
C  B_A       real*4 array          image data array
C  B_N       integer*4             number of elements in data array
C  B_PY      integer*4 array       line pointers
C  B_NYA     integer*4             number of active lines
C  A_NX      integer*4             number of line elements
C  A_YP      integer*4             offset of first active line
C  I         integer*4             x-coordinate
C  J         integer*4             y-coordinate
C  MINCR     integer*4             controls faint components
C  TRLM      real*4                detection threshold
C  Output parameters
C  AR        real*4                active radius
C 
C.Versions
C 
C 010118		last modif
C 
C-----------------------------------------------------------------------
      SUBROUTINE RADDET(A, JAPY, IBUF, I, J,
     &                  MINCR, TRLM, AVER, TRSH, AR) 
C
      IMPLICIT  NONE
C
      REAL      A(1)                    ! IN:  Image buffer
      INTEGER   JAPY(1)                 ! IN:  Pointers to image lines
      INTEGER   IBUF(4)                 ! IN:  Limits of image buffer
      INTEGER   I                       ! IN:  Object x-coordinate
      INTEGER   J                       ! IN:  Object y-coordinate
      INTEGER   MINCR                   ! IN:  Faint component control
      REAL      TRLM                    ! IN:  Detection level
      REAL      AVER                    ! IN:  9 pixels average
      REAL      TRSH                    ! IN:  Detection threshold
      REAL      AR                      ! OUT: Active radius
C
      INTEGER   ISTEP , IX
      INTEGER   JOF , JOFF , JY
      INTEGER   L
      INTEGER   MINCR2
      INTEGER   NX1 , NX2 , NY1 , NY2
C
      REAL   AC1 , AC2 , AC3
      REAL   CORR , DDR(8) , FCNTR
      REAL   SR2 , TEMP , TMPMX
C
      LOGICAL   UP , SKIP
C
      SR2 = SQRT(2.0)
      MINCR2 = NINT( FLOAT(MINCR) / SQRT(2.0) )
      MINCR2 = (MINCR+MINCR2) / 2
      NX2 = IBUF(3)
C 
      AC1 = 0.0				! was NOT initialized before...
      AC2 = 0.0
      Ac3 = 0.0
C
C *** JOF is the adrress of the start of central pixel line.
C
      JOFF = IBUF(2) - 1
      JOF = JAPY(J-JOFF)
C
C *** TEMP is a current pixel value.
C
      TEMP = A(JOF+I)
C
C *** FCNTR is the central pixel value.
C
      FCNTR = TEMP + 0.01 * TRSH
      IX = I + 1
      ISTEP = 0
      TMPMX = FCNTR
      UP = .TRUE.
      SKIP = .FALSE.
   10 CONTINUE
          IF ( TEMP .LT. TRLM .OR. IX .GE. NX2 ) GOTO 20
          IX = IX + 1
          TEMP = A(JOF+IX)
          IF ( UP .AND. IX-I .LE. MINCR .AND. TEMP .GE. TMPMX ) THEN
              TMPMX = TEMP
              ISTEP = ISTEP + 1
              SKIP = .FALSE.
          ELSE IF ( ( .NOT. SKIP ) .AND. TEMP .GT. TRLM .AND.
     &                                  TEMP .GE. TMPMX-TRSH ) THEN
              SKIP = .TRUE.
          ELSE
              UP = .FALSE.
          ENDIF
          GOTO 10
   20 CONTINUE
      CORR = A(JOF+IX-1) - TEMP
      IF ( CORR .GT. 0.0 ) THEN
          CORR = (TRLM-TEMP) / CORR
          IF ( CORR .LT. 0.0 ) THEN
              CORR = 0.0
          ELSE IF ( CORR .GT. 1.0 ) THEN
              CORR = 1.0
          ENDIF
      ELSE
          CORR = -1.0
      ENDIF
      DDR(1) = FLOAT(IX-I) - CORR
C
      TEMP = FCNTR
      IX = I - 1
      TMPMX = FCNTR
      NX1 = IBUF(1)
      UP = .TRUE.
      SKIP = .FALSE.
   30 CONTINUE
          IF ( TEMP .LT. TRLM .OR. IX .LE. NX1 ) GOTO 40
          IX = IX - 1
          TEMP = A(JOF+IX)
          IF ( UP .AND. I-IX .LE. MINCR .AND. TEMP .GE. TMPMX ) THEN
              TMPMX = TEMP
              ISTEP = ISTEP + 1
              SKIP = .FALSE.
          ELSE IF ( ( .NOT. SKIP ) .AND. TEMP .GE. TMPMX-TRSH
     &                           .AND. TEMP .GT. TRLM ) THEN
              SKIP = .TRUE.
          ELSE
              UP = .FALSE.
          ENDIF
          GOTO 30
   40 CONTINUE
      CORR = A(JOF+IX+1) - TEMP
      IF ( CORR .GT. 0.0 ) THEN
          CORR = (TRLM-TEMP) / CORR
          IF ( CORR .LT. 0.0 ) THEN
              CORR = 0.0
          ELSE IF ( CORR .GT. 1.0 ) THEN
              CORR = 1.0
          ENDIF
      ELSE
          CORR = -1.0
      ENDIF
      DDR(2) = FLOAT(I-IX) - CORR
C
      TEMP = FCNTR
      TMPMX = FCNTR
      NY2 = IBUF(4)
      JY = J + 1
      UP = .TRUE.
      SKIP = .FALSE.
   50 CONTINUE
          IF ( TEMP .LT. TRLM .OR. JY .GE. NY2 ) GOTO 60
          JY = JY + 1
          TEMP = A(JAPY(JY-JOFF)+I)
          IF ( UP .AND. JY-J .LE. MINCR .AND. TEMP .GE. TMPMX ) THEN
              TMPMX = TEMP
              ISTEP = ISTEP + 1
              SKIP = .FALSE.
          ELSE IF ( ( .NOT. SKIP ) .AND. TEMP .GE. TMPMX-TRSH
     &                          .AND. TEMP .GT. TRLM ) THEN
              SKIP = .TRUE.
          ELSE
              UP = .FALSE.
          ENDIF
          GOTO 50
   60 CONTINUE
      CORR = A(JAPY(JY-1-JOFF)+I) - TEMP
      IF ( CORR .GT. 0.0 ) THEN
          CORR = (TRLM-TEMP) / CORR
          IF ( CORR .LT. 0.0 ) THEN
              CORR = 0.0
          ELSE IF ( CORR .GT. 1.0 ) THEN
              CORR = 1.0
          ENDIF
      ELSE
          CORR = -1.0
      ENDIF
      DDR(3) = FLOAT(JY-J) - CORR
C
      JY = J - 1
      TEMP = FCNTR 
      TMPMX = FCNTR
      NY1 = IBUF(2)
      UP = .TRUE.
      SKIP = .FALSE.
   70 CONTINUE
          IF ( TEMP .LT. TRLM .OR. JY .LE. NY1 ) GOTO 80
          JY = JY - 1
          TEMP = A(JAPY(JY-JOFF)+I)
          IF ( UP .AND. J-JY .LE. MINCR .AND. TEMP .GE. TMPMX ) THEN
              TMPMX = TEMP
              ISTEP = ISTEP + 1
              SKIP = .FALSE.
          ELSE IF ( ( .NOT. SKIP ) .AND. TEMP .GE. TMPMX-TRSH
     &                      .AND. TEMP .GT. TRLM ) THEN
              SKIP = .TRUE.
          ELSE
              UP = .FALSE.
          ENDIF
          GOTO 70 
   80 CONTINUE
      CORR = A(JAPY(JY+1-JOFF)+I) - TEMP
      IF ( CORR .GT. 0.0 ) THEN
          CORR = (TRLM-TEMP) / CORR
          IF (CORR.LT.0.0) THEN
              CORR = 0.0
          ELSE IF ( CORR .GT. 1.0 ) THEN
              CORR = 1.0
          ENDIF
      ELSE
          CORR = -1.0
      ENDIF
      DDR(4) = FLOAT(J-JY) - CORR 
C
      IX = I + 1
      JY = J + 1
      TEMP = FCNTR
      TMPMX = FCNTR
      UP = .TRUE.
      SKIP = .FALSE.
   90 CONTINUE
          IF ( TEMP .LT. TRLM .OR. IX .GE. NX2 .OR. JY .GE. NY2 ) THEN
              GOTO 100
          ENDIF
          IX = IX + 1
          JY = JY + 1
          TEMP = A(JAPY(JY-JOFF)+IX)
          IF ( UP .AND. JY-J .LE. MINCR2 .AND. TEMP .GE. TMPMX ) THEN
              TMPMX = TEMP
              ISTEP = ISTEP + 1
              SKIP = .FALSE.
          ELSE IF ( ( .NOT. SKIP ) .AND. TEMP .GE. TMPMX-TRSH
     &                           .AND. TEMP .GT. TRLM ) THEN
              SKIP = .TRUE.
          ELSE
              UP = .FALSE.
          ENDIF
          GOTO 90
  100 CONTINUE
      CORR = A(JAPY(JY-1-JOFF)+IX-1) - TEMP
      IF ( CORR .GT. 0.0 ) THEN
          CORR = (TRLM-TEMP) / CORR
          IF ( CORR .LT. 0.0 ) THEN
              CORR = 0.0
          ELSE IF ( CORR .GT. 1.0 ) THEN
              CORR = 1.0
          ENDIF
      ELSE
          CORR = -1.0
      ENDIF
      DDR(5) = SR2 * ( FLOAT(JY-J) - CORR )
C
      IX = I - 1
      JY = J + 1
      TEMP = FCNTR
      TMPMX = FCNTR
      UP = .TRUE.
      SKIP = .FALSE.
  110 CONTINUE
          IF ( TEMP .LT. TRLM .OR. IX .LE. NX1 .OR. JY .GE. NY2 ) THEN
              GOTO 120
          ENDIF
          IX = IX - 1
          JY = JY + 1
          TEMP = A(JAPY(JY-JOFF)+IX)
          IF ( UP .AND. JY-J .LE. MINCR2 .AND. TEMP .GE. TMPMX ) THEN
              TMPMX = TEMP
              ISTEP = ISTEP + 1
              SKIP = .FALSE.
          ELSE IF ( ( .NOT. SKIP ) .AND. TEMP .GE. TMPMX-TRSH
     &                             .AND. TEMP .GT. TRLM ) THEN
              SKIP = .TRUE.
          ELSE
              UP = .FALSE.
          ENDIF
          GOTO 110
  120 CONTINUE
      CORR = A(JAPY(JY-1-JOFF)+IX+1) - TEMP
      IF ( CORR .GT. 0.0 ) THEN
          CORR = (TRLM-TEMP) / CORR
          IF ( CORR .LT. 0.0 ) THEN
              CORR = 0.0
          ELSE IF ( CORR .GT. 1.0 ) THEN
              CORR = 1.0
          ENDIF
      ELSE
          CORR = -1.0
      ENDIF
      DDR(6) = SR2 * ( FLOAT(JY-J) - CORR )
C
      IX = I - 1
      JY = J - 1
      TEMP = FCNTR
      TMPMX = FCNTR
      UP = .TRUE.
      SKIP = .FALSE.
  130 CONTINUE
          IF ( TEMP .LT. TRLM .OR. IX .LE. NX1 .OR. JY .LE. NY1 ) THEN
              GOTO 140
          ENDIF
          IX = IX - 1
          JY = JY - 1
          TEMP = A(JAPY(JY-JOFF)+IX)
          IF ( UP .AND. J-JY .LE. MINCR2 .AND. TEMP .GE. TMPMX ) THEN
              TMPMX = TEMP
              ISTEP = ISTEP + 1
              SKIP = .FALSE.
          ELSE IF ( ( .NOT. SKIP ) .AND. TEMP .GE. TMPMX-TRSH
     &                            .AND. TEMP .GT. TRLM ) THEN
              SKIP = .TRUE.
          ELSE
              UP = .FALSE.
          ENDIF
          GOTO 130
  140 CONTINUE
      CORR = A(JAPY(JY+1-JOFF)+IX+1) - TEMP
      IF ( CORR .GT. 0.0 ) THEN
          CORR = (TRLM-TEMP) / CORR
          IF ( CORR .LT. 0.0 ) THEN
              CORR = 0.0
          ELSE IF ( CORR .GT. 1.0 ) THEN
              CORR = 1.0
          ENDIF
      ELSE
          CORR = -1.0
      ENDIF
      DDR(7) = SR2 * ( FLOAT(J-JY) - CORR )
C
      IX = I + 1
      JY = J - 1
      TEMP = FCNTR
      TMPMX = FCNTR
      UP = .TRUE.
      SKIP = .FALSE.
  150 CONTINUE
          IF ( TEMP .LT. TRLM .OR. IX .GE. NX2 .OR. JY .LE. NY1 ) THEN
              GOTO 160
          ENDIF
          IX = IX + 1
          JY = JY - 1
          TEMP = A(JAPY(JY-JOFF)+IX)
          IF ( UP .AND. J-JY .LE. MINCR2 .AND. TEMP .GE. TMPMX ) THEN
              TMPMX = TEMP
              ISTEP = ISTEP + 1
              SKIP = .FALSE.
          ELSE IF ( ( .NOT. SKIP ) .AND. TEMP .GE. TMPMX-TRSH
     &                             .AND. TEMP .GT. TRLM ) THEN
              SKIP = .TRUE.
          ELSE
              UP = .FALSE.
          ENDIF
          GOTO 150
  160 CONTINUE
      CORR = A(JAPY(JY+1-JOFF)+IX-1) - TEMP
      IF ( CORR .GT. 0.0 ) THEN
          CORR = (TRLM-TEMP) / CORR
          IF ( CORR .LT. 0.0 ) THEN
              CORR = 0.0
          ELSE IF ( CORR .GT. 1.0 ) THEN
              CORR = 1.0
          ENDIF
      ELSE
          CORR = -1.0
      ENDIF
      DDR(8) = SR2 * ( FLOAT(J-JY) - CORR )
C
C *** Check if the detection is not spurious.
C
      IF ( ISTEP .GT. 2 ) THEN
          AR = -1.0
          RETURN
      ENDIF
C
C *** Calculate active radius.
C
      AC2 = FLOAT(NX2)
      AC3 = FLOAT(NX2)
      DO 170 L = 2 , 8
          IF ( DDR(L) .LE. AC1 ) THEN
              AC3 = AC2
              AC2 = AC1
              AC1 = DDR(L)
          ELSE IF ( DDR(L) .LE. AC2 ) THEN
              AC3 = AC2
              AC2 = DDR(L)
          ELSE IF ( DDR(L) .LT. AC3 ) THEN
              AC3 = DDR(L)
          ENDIF
  170 CONTINUE
      AR = (AC2+AC3) / 2.0
C
      RETURN
C
      END