C @(#)calfon.for	19.1 (ES0-DMD) 02/25/03 13:19:39
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
      SUBROUTINE CALFON(ETAB,IDIM,FONDM,FOND,R,IERR)                            
C +++                                                                           
C.PURPOSE: PROGRAMME DE CALCUL DU FOND AU VOISINAGE D'UN OBJET                  
C          PAR METHODE SEPARATION DE ZONE EN NIVEAUX(2)                         
C.INPUT/OUTPUT:  ETAB:  [E] tableau des valeurs entrees                         
C                IDIM:  [E] dimension du tableau                                
C                FONDM: [E] val. max. estimee du fond                           
C                FOND:  [S] valeur du fond                                      
C                R:     [S] dimension de l'objet(circulaire..)en pixels         
C                IERR:  [S] error parameter: =0:  ok                            
C                                            =-1: excessive background          
C                                            =-2: indefinite diameter           
C                                            =-3: the max. is out object        
C.VERSION: 80015  ???                          Created                          
C ---                                                                           
      IMPLICIT NONE	
      REAL     ETAB(1)
      INTEGER  IDIM
      REAL     FONDM
      REAL     FOND
      REAL     R
      INTEGER  IERR
C
      INTEGER  I
      INTEGER  IFIN
      INTEGER  ITEMP
      INTEGER  J
      INTEGER  K
      INTEGER  MADRID(1)
      INTEGER  MM
      INTEGER  N
      INTEGER  TWO
C
      REAL     DR
      REAL     ITAB(256)
      REAL     SIGMA
      REAL     SOM
      REAL     TAB(256)                                        
      REAL     TFOND
      REAL     TEMP
      REAL     VAL
      REAL     VTEST
      REAL     XMOY
      REAL     XMOYC
C
      COMMON  /VMR/MADRID                                                    
                                                                                
      FOND   = 0.                                                               
      R      = 0.                                                               
      IERR   = 0                                                                
      IFIN   = IDIM - 1                                                         
      DO 10 K = 1,IFIN                                                          
         ITAB(K) = K                                                            
   10 CONTINUE                                                                  
      DO 20 I = 2,IDIM                                                          
          TAB(I-1) = ETAB(I)                                                    
   20 CONTINUE                                                                  
C                                                                               
C *** sort by increassing order                                                 
C                                                                               
      DO 50 I = 1,IFIN - 1                                                      
         DO 40 J = I + 1,IFIN                                                   
            IF (TAB(I).LT.TAB(J)) GO TO 40                                      
            ITEMP   = ITAB(I)                                                   
            ITAB(I) = ITAB(J)                                                  
            ITAB(J) = ITEMP                                                    
            TEMP    = TAB(I)                                                    
            TAB(I)  = TAB(J)                                                    
            TAB(J)  = TEMP                                                      
   40    CONTINUE                                                               
   50 CONTINUE                                                                  
                                                                                
      IF (ITAB(IFIN).GT.3) IERR   = -3                                          
      IF (ITAB(IFIN).GT.3) RETURN                        ! max.hors de l'objet  
C                                                                               
C **** 1st background & sigma estimation with the 4 lesser values               
      VAL    = 0.                                                               
      SOM    = 0.                                                               
      DO 60 N = 1,4                                                             
         VAL    = VAL + TAB(N)*TAB(N)                                           
         SOM    = SOM + TAB(N)                                                  
   60 CONTINUE                                                                  
C                                                                               
      XMOY   = SOM/4.                                                           
      XMOYC  = SOM*SOM/4.                                                       
      SIGMA  = (VAL-XMOYC)/3.                                                   
      SIGMA  = SQRT(SIGMA)                                                      
      FOND   = XMOY                                                             
C                                                                               
C *** compute background & test                                                 
      MM     = 3                                                                
      VTEST  = XMOY + 5.*SIGMA                                                  
   70 MM     = MM + 1                                                           
      IF (TAB(MM).GT.VTEST) GO TO 80                                            
      GO TO 70                                                                  
C                                                                               
   80 CONTINUE                                                                  
      IF (MM.EQ.4) GO TO 100                                                    
      DO 90 K = 5,MM                                                            
          VAL    = VAL + TAB(K)*TAB(K)                                          
          SOM    = SOM + TAB(K)                                                 
   90 CONTINUE                                                                  
C                                                                               
      FOND   = SOM/FLOAT(MM)                                                    
      SIGMA  = (VAL-SOM*FOND)/FLOAT(MM-1)                                       
      SIGMA  = SQRT(SIGMA)                                                      
  100 IF (FOND.GT.FONDM) IERR   = -1                                            
      IF (FOND.GT.FONDM) RETURN                                                 
C                                                                               
C *** compute radius (to background+6*sigma level) & test                       
      TFOND  = FOND + 6.*SIGMA                                                  
      TWO    = 2
      DR     = ETAB(TWO)

      DO 110 J = 2,IDIM                                                         
         IF (ETAB(J).LE.TFOND) GO TO 120     ! test remontee(diam>480 microns)  
         IF (ETAB(J).GT.DR .AND. J.GT.8) GO TO 130                              
         DR     = ETAB(J)                                                       
  110 CONTINUE                                                                  
C                                                                               
      IERR = -2                                                                 
      RETURN                                                                    
C                                                                               
  120 R    = FLOAT(J-1)                                                         
      RETURN                                                                    
C                                                                               
  130 R    = FLOAT(J-2)                                                         
      RETURN                                                                    
      END                                                                       
