C @(#)ccdillum.for 17.1.1.1 (ES0-DMD) 01/25/02 17:49:36 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 MKILLU C++++++++++++++++++++++++++++++++++++++++++++++++++ C.IDENTIFICATION: ccdillum.FOR C.KEYWORDS: box smoothing C.PURPOSE: box smoothing with growing boxcar size C.VERSION: 930810 RHW Creation C-------------------------------------------------- IMPLICIT NONE INTEGER IAV INTEGER STATUS,MADRID INTEGER NAXIS,NPIX(2) INTEGER IMNOC,IMNOA INTEGER*8 PNTRA,PNTRC INTEGER KUN,KNUL INTEGER XBMIN, XBMAX INTEGER YBMIN, YBMAX C CHARACTER*60 FRAMEA,FRAMEC CHARACTER*72 IDENT CHARACTER*48 CUNIT CHARACTER*3 ACTION C DOUBLE PRECISION START(2),STEP(2) C REAL AVER REAL XBOX(2), YBOX(2) REAL XBMINR, XBMAXR REAL YBMINR, YBMAXR REAL CLIP(2) REAL LOWSIG,HIGSIG C INCLUDE 'MID_INCLUDE:ST_DEF.INC' COMMON/VMR/MADRID(1) INCLUDE 'MID_INCLUDE:ST_DAT.INC' C C *** do the work CALL STSPRO('ILLUM') CALL STKRDC('IN_A',1,1,60,IAV,FRAMEA,KUN,KNUL,STATUS) CALL STIGET(FRAMEA,D_R4_FORMAT,F_I_MODE,F_IMA_TYPE, 1 2,NAXIS,NPIX,START,STEP,IDENT,CUNIT, 2 PNTRA,IMNOA,STATUS) CALL STKRDC('OUT_A',1,1,60,IAV,FRAMEC,KUN,KNUL,STATUS) CALL STIPUT(FRAMEC,D_R4_FORMAT,F_O_MODE,F_IMA_TYPE, 1 NAXIS,NPIX,START,STEP,IDENT,CUNIT,PNTRC, 2 IMNOC,STATUS) C C *** smoothing boxes CALL STKRDR('XBOX',1,2,IAV,XBOX,KUN,KNUL,STATUS) CALL STKRDR('YBOX',1,2,IAV,YBOX,KUN,KNUL,STATUS) XBMINR = XBOX(1) XBMAXR = XBOX(2) YBMINR = YBOX(1) YBMAXR = YBOX(2) C C *** check the smoothing boxes IF (XBMINR.LT.1.) THEN XBMINR = XBMINR*NPIX(1) ENDIF IF (XBMAXR.LT.1.) THEN XBMAXR = XBMAXR*NPIX(1) ENDIF IF (YBMINR.LT.1.) THEN YBMINR = YBMINR*NPIX(2) ENDIF IF (YBMAXR.LT.1.) THEN YBMAXR = YBMAXR*NPIX(2) ENDIF C XBMIN = INT(MAX(2,MIN(NPIX(1),NINT(MIN(XBMINR,XBMAXR))))) XBMAX = INT(MAX(2,MIN(NPIX(1),NINT(MAX(XBMINR,XBMAXR))))) YBMIN = INT(MAX(2,MIN(NPIX(2),NINT(MIN(YBMINR,YBMAXR))))) YBMAX = INT(MAX(2,MIN(NPIX(2),NINT(MAX(YBMINR,YBMAXR))))) C C *** to clip or not to clip CALL STKRDC('CLIP',1,1,3,IAV,ACTION,KUN,KNUL,STATUS) CALL UPCAS(ACTION,ACTION) IF (ACTION(1:1).EQ.'Y') THEN CALL STKRDR('SIGMA',1,2,IAV,CLIP,KUN,KNUL,STATUS) LOWSIG = CLIP(1) HIGSIG = CLIP(2) ENDIF C C *** test on clipping IF (ACTION(1:1).EQ.'Y') THEN CALL ILLUM(PNTRA,PNTRC,NPIX,XBMIN,XBMAX,YBMIN,YBMAX, 2 LOWSIG,HIGSIG,AVER) ELSE CALL QILLUM(PNTRA,PNTRC,NPIX,XBMIN,XBMAX,YBMIN,YBMAX,AVER) ENDIF C CALL STDCOP(IMNOA,IMNOC,1,' ',STATUS) CALL STDWRR(IMNOC,'CCDMEAN',AVER,1,1,KUN,STATUS) CALL STSEPI END SUBROUTINE ILLUM(IN,OUT,NPIX,XBMIN,XBMAX,YBMIN,YBMAX, 2 LOW,HIGH,MEAN) C IMPLICIT NONE INTEGER IN, OUT INTEGER NPIX(2) INTEGER XBMIN,XBMAX,YBMIN,YBMAX REAL LOW, HIGH REAL MEAN C INTEGER LINEIN, LINEOUT INTEGER IMNOS, IMNOG, IMNOP INTEGER*8 IPOINT, PTRS, SUM, AVG, OUTPUT INTEGER*8 PTR INTEGER I, NITER, NREJ INTEGER MADRID(1) INTEGER NCOLS,NLINES INTEGER STAT, IAV INTEGER YBOX2 REAL SCALE REAL ILLSUM EXTERNAL ILLSUM INTEGER BXCLN EXTERNAL BXCLN INCLUDE 'MID_INCLUDE:ST_DEF.INC' COMMON /VMR/MADRID INCLUDE 'MID_INCLUDE:ST_DAT.INC' DATA NITER/10/ C C *** here create vitual memory NCOLS = NPIX(1) NLINES = NPIX(2) C CALL STFCRE('dummsumm',D_R4_FORMAT,F_X_MODE,F_IMA_TYPE, 2 NCOLS,IMNOS,STAT) ! create sum array CALL STFMAP(IMNOS,F_X_MODE,1,NCOLS,IAV,SUM,STAT) C CALL STFCRE('dummavg',D_R4_FORMAT,F_X_MODE,F_IMA_TYPE, 2 NCOLS,IMNOG,STAT) !create output array CALL STFMAP(IMNOG,F_X_MODE,1,NCOLS,IAV,AVG,STAT) OUTPUT = AVG C CALL STFCRE('dummpntr',D_I4_FORMAT,F_X_MODE,F_IMA_TYPE, 2 YBMAX,IMNOP,STAT) ! create pointer array CALL STFMAP(IMNOP,F_X_MODE,1,YBMAX,IAV,PTRS,STAT) C C *** now start the work C IF (YBMAX.LT.NLINES) THEN C YBMAX = NLINES C ENDIF C C *** accumulate the minimum y box CALL ILLZER(MADRID(SUM),NCOLS) ! initialize sum array C C *** acumulate the minimum y box LINEIN = 0 10 CONTINUE IF (LINEIN.LT.YBMIN) THEN LINEIN = LINEIN + 1 IPOINT = IN + (LINEIN-1)*NCOLS CALL ILLADD(MADRID(IPOINT),MADRID(SUM),MADRID(SUM),NCOLS) PTR = PTRS + MOD(LINEIN,YBMAX) MADRID(PTR) = IPOINT GOTO 10 ENDIF C C *** output the minimum y box YBOX2 = YBMIN SCALE = FLOAT(YBMIN) CALL AGBCAR(MADRID(SUM),MADRID(AVG),NCOLS,XBMIN,XBMAX,SCALE) C C *** Iteratively clean the initial lines PTR = PTRS IF (YBOX2 .NE. YBMAX) THEN PTR = PTR+1 ENDIF C DO I = 1, NITER NREJ = 0 DO LINEOUT = 1,LINEIN IPOINT = MADRID(PTR+LINEOUT-1) NREJ = NREJ + BXCLN(MADRID(IPOINT),MADRID(AVG), 2 MADRID(SUM),NCOLS,LOW,HIGH) ENDDO IF (NREJ.GT.0) THEN CALL AGBCAR(MADRID(SUM),MADRID(AVG),NCOLS,XBMIN,XBMAX,SCALE) ELSE GOTO 20 ENDIF ENDDO C 20 CONTINUE YBOX2 = (YBMIN+1)/2 LINEOUT = 0 30 CONTINUE IF (LINEOUT.LT.YBOX2) THEN LINEOUT = LINEOUT + 1 IPOINT = OUT + (LINEOUT-1)*NCOLS CALL COPYF(MADRID(OUTPUT),MADRID(IPOINT),NCOLS) GO TO 30 ENDIF MEAN = YBOX2*ILLSUM(MADRID(OUTPUT),NCOLS) ! calculate the MEAN C C *** increase the y box size by steps of 2 until the maximum szie. 40 CONTINUE IF (LINEIN.LT.YBMAX) THEN LINEIN = LINEIN + 1 IPOINT = IN + (LINEIN-1)*NCOLS CALL ILLADD(MADRID(SUM),MADRID(IPOINT),MADRID(SUM),NCOLS) PTR = PTRS + MOD(LINEIN,YBMAX) MADRID(PTR) = IPOINT SCALE = SCALE + 1 C NREJ = BXCLN(MADRID(IPOINT),MADRID(AVG), 2 MADRID(SUM),NCOLS,LOW,HIGH) CALL AGBCAR(MADRID(SUM),MADRID(AVG),NCOLS,XBMIN,XBMAX,SCALE) C LINEIN = LINEIN + 1 IPOINT = IN + (LINEIN-1)*NCOLS CALL ILLADD(MADRID(SUM),MADRID(IPOINT),MADRID(SUM),NCOLS) PTR = PTRS + MOD(LINEIN,YBMAX) MADRID(PTR) = IPOINT C NREJ = BXCLN(MADRID(IPOINT),MADRID(AVG),MADRID(SUM), 2 NCOLS,LOW,HIGH) SCALE = SCALE + 1 CALL AGBCAR(MADRID(SUM),MADRID(AVG),NCOLS,XBMIN,XBMAX,SCALE) C LINEOUT = LINEOUT + 1 IPOINT = OUT + (LINEOUT-1)*NCOLS CALL COPYF(MADRID(AVG),MADRID(IPOINT),NCOLS) MEAN = MEAN + ILLSUM(MADRID(IPOINT),NCOLS) ! calculate MEAN GOTO 40 ENDIF C C *** for each line subtract the last line from the sum, add the C next line to the sum, and output a line 50 CONTINUE IF (LINEIN.LT.NLINES) THEN LINEIN = LINEIN + 1 PTR = PTRS + MOD(LINEIN,YBMAX) IPOINT = MADRID(PTR) CALL ILLSUB(MADRID(SUM),MADRID(IPOINT),MADRID(SUM),NCOLS) IPOINT = IN + (LINEIN-1)*NCOLS CALL ILLADD(MADRID(SUM),MADRID(IPOINT),MADRID(SUM),NCOLS) MADRID(PTR) = IPOINT C NREJ = BXCLN(MADRID(IPOINT),MADRID(AVG),MADRID(SUM), 2 NCOLS,LOW,HIGH) C LINEOUT = LINEOUT + 1 IPOINT = OUT + (LINEOUT-1)*NCOLS CALL AGBCAR(MADRID(SUM),MADRID(AVG),NCOLS,XBMIN,XBMAX,SCALE) CALL COPYF(MADRID(AVG),MADRID(IPOINT),NCOLS) MEAN = MEAN + ILLSUM(MADRID(IPOINT),NCOLS) ! calculate the MEAN GOTO 50 ENDIF C C *** decrease the y box in steps of 2 until mimimum y box 60 CONTINUE IF (LINEOUT.LT.(NLINES-YBOX2)) THEN LINEIN = LINEIN + 1 PTR = PTRS + MOD(LINEIN,YBMAX) IPOINT = MADRID(PTR) CALL ILLSUB(MADRID(SUM),MADRID(IPOINT),MADRID(SUM),NCOLS) LINEIN = LINEIN + 1 PTR = PTRS + MOD(LINEIN,YBMAX) IPOINT = MADRID(PTR) CALL ILLSUB(MADRID(SUM),MADRID(IPOINT),MADRID(SUM),NCOLS) SCALE = SCALE-2 C LINEOUT = LINEOUT + 1 IPOINT = OUT + (LINEOUT-1)*NCOLS CALL AGBCAR(MADRID(SUM),MADRID(IPOINT),NCOLS,XBMIN,XBMAX,SCALE) MEAN = MEAN + ILLSUM(MADRID(IPOINT),NCOLS) ! calculate the MEAN GOTO 60 ENDIF C C *** output the last liens of the minimum y box size CALL AGBCAR(MADRID(SUM),MADRID(AVG),NCOLS,XBMIN,XBMAX,SCALE) 70 CONTINUE IF (LINEOUT.LT.NLINES) THEN LINEOUT = LINEOUT + 1 IPOINT = OUT + (LINEOUT-1)*NCOLS CALL COPYF(MADRID(OUTPUT),MADRID(IPOINT),NCOLS) GOTO 70 ENDIF C C *** write the scaling vector MEAN = MEAN + (NLINES-LINEOUT)*ILLSUM(MADRID(AVG),NCOLS) MEAN = MEAN/(NCOLS*NLINES) C RETURN END SUBROUTINE QILLUM(IN,OUT,NPIX,XBMIN,XBMAX,YBMIN,YBMAX,MEAN) C IMPLICIT NONE INTEGER IN, OUT INTEGER NPIX(2) INTEGER XBMIN,XBMAX,YBMIN,YBMAX REAL MEAN C INTEGER LINEIN, LINEOUT INTEGER IMNOS, IMNOO, IMNOP INTEGER IPOINT, PTRS, SUM, OUTPUT INTEGER PTR INTEGER MADRID(1) INTEGER NCOLS,NLINES INTEGER STAT, IAV INTEGER YBOX1 REAL SCALE REAL ILLSUM EXTERNAL ILLSUM INCLUDE 'MID_INCLUDE:ST_DEF.INC' COMMON /VMR/MADRID INCLUDE 'MID_INCLUDE:ST_DAT.INC' C C *** here create vitual memory NCOLS = NPIX(1) NLINES = NPIX(2) C CALL STFCRE('dummsumm',D_R4_FORMAT,F_X_MODE,F_IMA_TYPE, 2 NCOLS,IMNOS,STAT) ! create sum array CALL STFMAP(IMNOS,F_X_MODE,1,NCOLS,IAV,SUM,STAT) C CALL STFCRE('dummoutp',D_R4_FORMAT,F_X_MODE,F_IMA_TYPE, 2 NCOLS,IMNOO,STAT) !create output array CALL STFMAP(IMNOO,F_X_MODE,1,NCOLS,IAV,OUTPUT,STAT) C CALL STFCRE('dummpntr',D_I4_FORMAT,F_X_MODE,F_IMA_TYPE, 2 YBMAX,IMNOP,STAT) ! create pointer array CALL STFMAP(IMNOP,F_X_MODE,1,YBMAX,IAV,PTRS,STAT) C C *** now start the work C IF (YBMAX.LT.NLINES) THEN C YBMAX = NLINES C ENDIF C C *** accumulate the minimum y box CALL ILLZER(MADRID(SUM),NCOLS) ! initialize sum array C C *** acumulate the minimum y box LINEIN = 0 10 CONTINUE IF (LINEIN.LT.YBMIN) THEN LINEIN = LINEIN + 1 IPOINT = IN + (LINEIN-1)*NCOLS CALL ILLADD(MADRID(IPOINT),MADRID(SUM),MADRID(SUM),NCOLS) PTR = PTRS + MOD(LINEIN,YBMAX) MADRID(PTR) = IPOINT GOTO 10 ENDIF C C *** output the minimum y box YBOX1 = (YBMIN+1)/2 SCALE = FLOAT(YBMIN) CALL AGBCAR(MADRID(SUM),MADRID(OUTPUT),NCOLS,XBMIN,XBMAX,SCALE) LINEOUT = 0 20 CONTINUE IF (LINEOUT.LT.YBOX1) THEN LINEOUT = LINEOUT + 1 IPOINT = OUT + (LINEOUT-1)*NCOLS CALL COPYF(MADRID(OUTPUT),MADRID(IPOINT),NCOLS) GO TO 20 ENDIF MEAN = YBOX1*ILLSUM(MADRID(OUTPUT),NCOLS) ! calculate the MEAN C C *** increase the y box size by steps of 2 until the maximum szie. 30 CONTINUE IF (LINEIN.LT.YBMAX) THEN LINEIN = LINEIN + 1 IPOINT = IN + (LINEIN-1)*NCOLS CALL ILLADD(MADRID(SUM),MADRID(IPOINT),MADRID(SUM),NCOLS) PTR = PTRS + MOD(LINEIN,YBMAX) MADRID(PTR) = IPOINT LINEIN = LINEIN + 1 IPOINT = IN + (LINEIN-1)*NCOLS CALL ILLADD(MADRID(SUM),MADRID(IPOINT),MADRID(SUM),NCOLS) PTR = PTRS + MOD(LINEIN,YBMAX) MADRID(PTR) = IPOINT C SCALE = SCALE + 2 LINEOUT = LINEOUT + 1 IPOINT = OUT + (LINEOUT-1)*NCOLS CALL AGBCAR(MADRID(SUM),MADRID(IPOINT),NCOLS,XBMIN,XBMAX,SCALE) MEAN = MEAN + ILLSUM(MADRID(IPOINT),NCOLS) ! calculate MEAN GOTO 30 ENDIF C C *** for each line subtract the last line from the sum, add the C next line to the sum, and output a line 40 CONTINUE IF (LINEIN.LT.NLINES) THEN LINEIN = LINEIN + 1 PTR = PTRS + MOD(LINEIN,YBMAX) IPOINT = MADRID(PTR) CALL ILLSUB(MADRID(SUM),MADRID(IPOINT),MADRID(SUM),NCOLS) IPOINT = IN + (LINEIN-1)*NCOLS CALL ILLADD(MADRID(SUM),MADRID(IPOINT),MADRID(SUM),NCOLS) MADRID(PTR) = IPOINT C LINEOUT = LINEOUT + 1 IPOINT = OUT + (LINEOUT-1)*NCOLS CALL AGBCAR(MADRID(SUM),MADRID(IPOINT),NCOLS,XBMIN,XBMAX,SCALE) MEAN = MEAN + ILLSUM(MADRID(IPOINT),NCOLS) ! calculate the MEAN GOTO 40 ENDIF C C *** decrease the y box in steps of 2 until mimimum y box 50 CONTINUE IF (LINEOUT.LT.(NLINES-YBOX1)) THEN LINEIN = LINEIN + 1 PTR = PTRS + MOD(LINEIN,YBMAX) IPOINT = MADRID(PTR) CALL ILLSUB(MADRID(SUM),MADRID(IPOINT),MADRID(SUM),NCOLS) LINEIN = LINEIN + 1 PTR = PTRS + MOD(LINEIN,YBMAX) IPOINT = MADRID(PTR) CALL ILLSUB(MADRID(SUM),MADRID(IPOINT),MADRID(SUM),NCOLS) C LINEOUT = LINEOUT + 1 SCALE = SCALE -2 IPOINT = OUT + (LINEOUT-1)*NCOLS CALL AGBCAR(MADRID(SUM),MADRID(IPOINT),NCOLS,XBMIN,XBMAX,SCALE) MEAN = MEAN + ILLSUM(MADRID(IPOINT),NCOLS) ! calculate the MEAN GOTO 50 ENDIF C C *** output the last liens of the minimum y box size CALL AGBCAR(MADRID(SUM),MADRID(OUTPUT),NCOLS,XBMIN,XBMAX,SCALE) 60 CONTINUE IF (LINEOUT.LT.NLINES) THEN LINEOUT = LINEOUT + 1 IPOINT = OUT + (LINEOUT-1)*NCOLS CALL COPYF(MADRID(OUTPUT),MADRID(IPOINT),NCOLS) GOTO 60 ENDIF C C *** write the scaling vector MEAN = MEAN + (NLINES-LINEOUT)*ILLSUM(MADRID(OUTPUT),NCOLS) MEAN = MEAN/(NCOLS*NLINES) C RETURN END C C C INTEGER FUNCTION BXCLN(DATA,BXAVG,SUM,NCOLS,LOW,HIGH) C C ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ C Reject data values from thr sum for the next boxcar average which C exceed the minimum and maximum residual values from the current C boxcar average. This excludes data from the moving avarge before it C enters the average. C----------------------------------------------------------------------- C IMPLICIT NONE REAL DATA(1) REAL BXAVG(1) REAL SUM(1) INTEGER NCOLS REAL LOW REAL HIGH C INTEGER I,NREJ REAL RMS,RESID,MINRES,MAXRES C RMS = 0.0 DO I = 1,NCOLS RMS = RMS + (DATA(I)-BXAVG(1))**2 ENDDO RMS = SQRT(RMS/NCOLS) MINRES = -LOW*RMS MAXRES = HIGH*RMS C NREJ = 0 DO I = 1,NCOLS RESID = DATA(I)-BXAVG(I) IF ((RESID.LT.MINRES) .OR. (RESID.GT.MAXRES)) THEN DATA(I) = BXAVG(I) SUM(I) = SUM(I)-RESID NREJ = NREJ+1 ENDIF ENDDO C BXCLN = NREJ RETURN END C C C SUBROUTINE ILLZER(A,NC) C IMPLICIT NONE REAL A(1) INTEGER NC INTEGER IC C DO IC = 1,NC A(IC) = 0.0 ENDDO RETURN END C C C SUBROUTINE ILLADD(A,B,C,NC) C IMPLICIT NONE REAL A(1),B(1),C(1) INTEGER NC INTEGER IC C DO IC = 1,NC C(IC) = A(IC)+B(IC) ENDDO RETURN END C C C SUBROUTINE ILLSUB(A,B,C,NC) C IMPLICIT NONE REAL A(1),B(1),C(1) INTEGER NC INTEGER IC C DO IC = 1,NC C(IC) = A(IC)-B(IC) ENDDO RETURN END C C C REAL FUNCTION ILLSUM(A,NC) C IMPLICIT NONE REAL A(1) INTEGER NC INTEGER IC C ILLSUM = 0.0 DO IC = 1,NC ILLSUM = ILLSUM + A(IC) ENDDO RETURN END C C C SUBROUTINE AGBCAR(INA,OUTA,NCOLS,XBMIN,XBMAX,YBOX) C C +++ Vector growing boxcar smooth. Taken from the IRAF CCD package C REAL INA(1) REAL OUTA(1) INTEGER NCOLS INTEGER XBMIN,XBMAX REAL YBOX C INTEGER COLIN, COLOUT, LASTCOL INTEGER NPIX, XBMIN2 REAL SUM,OUTPUT C XBMIN2 = (XBMIN+1)/2 COLIN = 0 SUM = 0.0 10 CONTINUE IF (COLIN.LT.XBMIN) THEN COLIN = COLIN + 1 SUM = SUM + INA(COLIN) GOTO 10 ENDIF C NPIX = XBMIN * YBOX OUTPUT = SUM/NPIX COLOUT = 0 20 CONTINUE IF (COLOUT.LT.XBMIN2) THEN COLOUT = COLOUT + 1 OUTA(COLOUT) = OUTPUT GO TO 20 ENDIF C 30 CONTINUE IF (COLIN.LT.XBMAX) THEN COLIN = COLIN + 1 SUM = SUM + INA(COLIN) COLIN = COLIN + 1 SUM = SUM + INA(COLIN) NPIX = NPIX + 2*YBOX COLOUT = COLOUT + 1 OUTA(COLOUT) = SUM/NPIX GO TO 30 ENDIF C LASTCOL = 0 40 CONTINUE IF (COLIN.LT.NCOLS) THEN COLIN = COLIN + 1 LASTCOL = LASTCOL +1 SUM = SUM + INA(COLIN) - INA(LASTCOL) COLOUT = COLOUT + 1 OUTA(COLOUT) = SUM/NPIX GO TO 40 ENDIF C 50 CONTINUE IF (COLOUT.LT. NCOLS - XBMIN2) THEN LASTCOL = LASTCOL + 1 SUM = SUM - INA(LASTCOL) LASTCOL = LASTCOL + 1 SUM = SUM - INA(LASTCOL) NPIX = NPIX - 2*YBOX COLOUT = COLOUT + 1 OUTA(COLOUT)= SUM/NPIX GO TO 50 ENDIF C OUTPUT = SUM/NPIX 70 CONTINUE IF (COLOUT.LT.NCOLS) THEN COLOUT = COLOUT + 1 OUTA(COLOUT) = OUTPUT GO TO 70 ENDIF RETURN END