C @(#)updat.for 17.1.1.1 (ES0-DMD) 01/25/02 17:19:02 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.IDENT: UPDAT C.PURPOSE: Update the mask frame and the buffer setting to a flag value points C with residual greater than kappa*sigma C.INPUT: MASK : INPUT_B is the frame masking bright objects in the image C ATEMP : a buffer containing the rebinned input image - used for speed C FLMSK : flag = .TRUE. if MASK exists C REJ : value assigned to 'bad' pixels in the buffer ATEMP C NMX,... usual descriptors of the frame MASK C IDEG : polinomial degree C RKAPPA: parameter for kappa-sigma clipping C SIGM1 : standard deviation of residuals image-fit.pol. for the C 'good' pixels C C : up to 21 coefficients of the polinomial surface C NSX,... usual descriptors for the temporary frame ATEMP C.OUTPUT: none C ------------------------------------------------------------------ SUBROUTINE UPDAT(MASK,ATEMP,FLMSK,REJ,NMX,NMY,RXSTR,RYSTR,RXSTP, + RYSTP,NSX,NSY,SXSTR,SYSTR,SXSTP,SYSTP,IDEG, + RKAPPA,SIGM1,C) IMPLICIT NONE C INTEGER NMX INTEGER NMY REAL MASK(NMX,NMY) REAL ATEMP(128,128) LOGICAL FLMSK REAL REJ DOUBLE PRECISION RXSTR DOUBLE PRECISION RYSTR DOUBLE PRECISION RXSTP DOUBLE PRECISION RYSTP INTEGER NSX INTEGER NSY DOUBLE PRECISION SXSTR DOUBLE PRECISION SYSTR DOUBLE PRECISION SXSTP DOUBLE PRECISION SYSTP INTEGER IDEG REAL RKAPPA REAL SIGM1 DOUBLE PRECISION C(21) C INTEGER I INTEGER I1MSK INTEGER IX, IY INTEGER J INTEGER J1MSK INTEGER NCV(5) REAL DCX, DCY REAL POLIN REAL RES REAL X1 REAL Y1 REAL V DOUBLE PRECISION CC(6) DOUBLE PRECISION XD,YD REAL BTSIG C EXTERNAL POLIN C DATA NCV/3,6,10,15,21/ C C *** calculate residuals and update mask and temporary array C DCX = (NSX-1.)/2. DCY = (NSY-1.)/2. C C *** loop : update mask and set -rej- flag in temporary array C BTSIG = RKAPPA*SIGM1 DO 70 J = 1,NSY IY = J - 1 Y1 = SYSTR + SYSTP*FLOAT(IY) IF (FLMSK) THEN J1MSK = (Y1-RYSTR)/RYSTP + 1.5 END IF YD = DCY - IY CALL COEFY(YD,C,CC,IDEG) C C *** check mask C IF ( .NOT. FLMSK) GO TO 30 ! y outofbounds IF (J1MSK.LT.1 .OR. J1MSK.GT.NMY) GO TO 30 C DO 20 I = 1,NSX IX = I - 1 XD = IX - DCX RES = ATEMP(I,J) - POLIN(XD,CC,IDEG) X1 = SXSTR + SXSTP*FLOAT(IX) I1MSK = (X1-RXSTR)/RXSTP + 1.5 IF (I1MSK.LT.1 .OR. I1MSK.GT.NMX) THEN IF (ATEMP(I,J).EQ.REJ) THEN GO TO 10 END IF IF (ABS(RES).GT.BTSIG) THEN ATEMP(I,J) = REJ END IF ELSE V = MASK(I1MSK,J1MSK) IF (V.LE.0.) THEN GO TO 10 END IF IF (ABS(RES).LE.BTSIG) THEN GO TO 10 END IF MASK(I1MSK,J1MSK) = -1 ATEMP(I,J) = REJ END IF 10 CONTINUE 20 CONTINUE GO TO 60 C 30 DO 50 I = 1,NSX IF (ATEMP(I,J).EQ.REJ) THEN GO TO 40 END IF IX = I - 1 XD = IX - DCX RES = ATEMP(I,J) - POLIN(XD,CC,IDEG) IF (ABS(RES).GT.BTSIG) THEN ATEMP(I,J) = REJ END IF 40 CONTINUE 50 CONTINUE C 60 CONTINUE 70 CONTINUE C RETURN END