C @(#)rlgrnt.for	19.1 (ES0-DMD) 02/25/03 13:25:38
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
C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C.IDENTIFICATION
C  subroutine RLGRNT         version 2.3       830916
C  A. Kruszewski             ESO Garching
C  modified                  version 2.4       870304
C  A. Kruszewski             Obs. de Geneve
C.KEYWORDS
C  profile gradient
C.PURPOSE
C  calculates gradient in an observed profile relative to the
C  standard stellar profile
C  this relative gradient should be close to zero for stars
C  negative for galaxies, and positive for majority of defects
C.ALGORYTHM
C  quantity  -LOG10(AVPR)-SPRF , or in other words, deviations of
C  an observed profile from the standard profile, is subjected to
C  the regression with respect to a ring number K expressed
C  in pixels
C  the coefficient of regression is accepted as relative gradient
C.INPUT/OUTPUT
C  input arguments
C  AVPR        real*4 array      observed one-dimensional profile
C  SPRF        real*4 array      cumulative standard profile
C  HHCUT       real*4            upper limit of usable data
C  TRSH        real*4            lower limiting treshold
C  output arguments
C  AVGR        real*4            relative gradient
C  SIGMA       real*4            sigma of single data point
C-----------------------------------------------------------------------
      SUBROUTINE RLGRNT(AVPR, SPRF, HHCUT, TRSH, AVGR, SIGMA)
C
      IMPLICIT  NONE
      INCLUDE  'MID_REL_INCL:INVENT.INC/NOLIST'
C
      REAL      AVPR(0:MAXSUB)
      REAL      SPRF(0:MAXSUB)
      REAL      HHCUT
      REAL      TRSH
      REAL      AVGR
      REAL      SIGMA
C
      INTEGER   K
      INTEGER   M1
      INTEGER   N, NN
      INTEGER   L
      REAL      ALIM, AMEAN
      REAL      BLIM
      REAL      GRSL
      REAL      ZRLV
      REAL      SGGR
      REAL      IC, ID
      INTEGER   IW(0:MAXSUB )
      REAL      X(0:MAXSUB), Y(0:MAXSUB)
      REAL      DVGR(0:MAXSUB)
      REAL      YY(MAXSUB+1)
C
      DO 10 K = 0 , MAXSUB
          IW(K) = 0
          X(K) = 0.0
          Y(K) = 0.0
   10 CONTINUE
      M1 = MAXSUB + 1
      DO 20 K = 1 , M1
          YY(K) = 0.0
   20 CONTINUE
      NN = 0
      BLIM = MIN( AVPR(0) , HHCUT )
      ALIM = MAX( 0.01*BLIM , 0.5*TRSH )
      L = 0
C
C ******      Select usable part of an observed profile.
C
   31 CONTINUE
          IF ( AVPR(L) .LT. ALIM .OR. L .GT. MAXSUB ) THEN
              N = L - 1
              GOTO 30
          ENDIF
          IF ( AVPR(L) .LT. HHCUT ) THEN
              X(L) = FLOAT(L)
              Y(L) = -ALOG10( AVPR(L) ) - SPRF(L)
              IW(L) = 1
              NN = NN + 1
              YY(NN) = Y(L)
          ELSE
              IW(L) = 0
          ENDIF
          L = L + 1
          GOTO 31
   30 CONTINUE
C
C ******      Calculate root mean square deviation SIGMA of used data.
C
      IF ( NN .GT. 1 ) THEN
          CALL MEAN( YY , NN , AMEAN , SIGMA )
      ELSE
          SIGMA = 1.0
      ENDIF
C
C ******      Calculate first approximation for relative gradient.
C
      IF ( NN .GE. 2 ) THEN
          CALL GRADET( L , X , Y , IW , AVGR , ZRLV , DVGR , SGGR )
          IF ( NN .EQ. 2 ) THEN
              IC = 0
              SGGR = 1.0
          ELSE
              IC = 3
          ENDIF
      ELSE
          AVGR = 1.0
          SGGR = 1.0
          IC = 0
      ENDIF
C
C ******      Repeat calculations after rejecting
C ******      points deviating more than 1.4*SGGR.
C
   41 CONTINUE
      IF ( IC .EQ. 0 ) GOTO 40
          GRSL = 1.4 * SGGR
          NN = 0
          ID = 0
          DO 50 L = 0 , N
              IF ( IW(L) .EQ. 1 ) THEN
                  IF ( ABS( DVGR(L) ) .GT. GRSL ) THEN
                      IW(L) = 0
                      ID = ID + 1
                  ELSE
                      NN = NN + 1
                  ENDIF
              ENDIF
   50     CONTINUE
          IF ( ID .EQ. 0 .OR. NN .LT. 2 ) THEN
              IC = 0
          ELSE
              CALL GRADET( MAXSUB , X , Y , IW , AVGR ,
     &                                    ZRLV , DVGR , SGGR )
              IF ( NN .GT. 2 ) THEN
                  IC = IC - 1
              ELSE 
                  IC = 0
                  SGGR = 1.0
              ENDIF
          ENDIF
          GOTO 41
   40 CONTINUE
C
      RETURN
      END
