C @(#)smosub2.for 17.1.1.1 (ES0-DMD) 01/25/02 17:40: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 SUBROUTINE AVER2A(METH,A,B,NX,NY,RADIUS,TR,ALL,V,SUM,CUTS) C C ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ C C.IDENTIFICATION C subroutine AVER2A version 1.00 861216 C K. Banse ESO - Garching C 1.10 870227 1.20 870302 1.50 880224 C 1.60 910521 C C.KEYWORDS C smoothing, neighbourhood C C.PURPOSE C remove spot noise C C.ALGORITHM C a) smooth image by taking the average A(x,y)of all points in neighbourhood C of x,y with given radius, and replace original value I(x,y) by A(x,y), C if average differs more than THRESH*I(x,y) from original value C b) as a) but take median instead of average C C.INPUT/OUTPUT C call as AVER2A(METH,A,B,NX,NY,RADIUS,TR,ALL,V,SUM) C C input par: C METH: char.exp. (1:1) = method to use: always A for average C (2:2) = A for absolute, C R for relative threshold, C F for fixed threshold interval C A: R*4 array input array C B: R*4 array output array C NX: I*4 no. of pixels per line C NY: I*4 no. of lines in input array C RADIUS: I*4 array "radius" of neighbourhood in x and y, C i.e. no. of pixels in horizontal and vertical C direction to include C TR: R*4 (array) threshold(s) for center point to determine C replacement C ALL: I*4 flag for inclusion of center point in averaging C = 1: include center point C = 0: do not include center point C V: R*8 working buff. of initial rowsums (RADIUS(2)*2+1) C SUM: R*8 working buff. of sums (NPIX(1)) C C.VERSIONS C 1.00 has been built on 2.80 of old module AVERAGE C 1.10 fix bug with index DIFX in relative threshold option C (it was 1 too low...) C 1.20 add fixed threshold interval option C 1.50 move to FORTRAN 77 C 1.60 also filter the borders (input frame larger now) C C -------------------------------------------------------------------- C IMPLICIT NONE C INTEGER NX,NY,ALL INTEGER DIFX,ICOUNT INTEGER INDX,BINDX INTEGER N,NN,NL,NO,NOZ INTEGER NP,NRL,NXMRX INTEGER OFF,RX,RX1,RXM1,RY,RYNX,TOP INTEGER RADIUS(2),NW1,NW2,NOSUM C REAL A(*),B(*),TR(*),CUTS(*) REAL W,TEST,COUNT,FACTOR,THRESH C DOUBLE PRECISION V(*),SUM(NX) DOUBLE PRECISION VNEW C CHARACTER*(*) METH C C initialize C BINDX = 1 RX = RADIUS(1) RY = RADIUS(2) NW1 = RX*2 + 1 !no. of values on row of filter NW2 = RY*2 + 1 !no. of values on column of filter ICOUNT = NW1 * NW2 COUNT = FLOAT(ICOUNT) !no. of pixels in neighbourhood with center C RYNX = RY * NX NXMRX = NX - RX RX1 = RX + 1 RXM1 = RX - 1 OFF = RYNX !offset for center row C IF (ALL.EQ.0) THEN FACTOR = 1. / (COUNT-1.) ELSE FACTOR = 1. / COUNT ENDIF C C fork according to test criterion (for speed reasons) IF (METH(2:2).EQ.'F') THEN GOTO 4000 ELSE THRESH = TR(1) ENDIF IF (THRESH.GT.0.) GOTO 2000 C C C -------------------------- C C here we handle boxcar filtering without thresholding C C -------------------------- C C C set up first row sums NOZ = 0 !offset for lower left C DO 300 N=1,NW2 !we have NW2 row sums V(N) = 0.D0 DO 200 NP=NOZ+1,NOZ+NW1 V(N) = V(N) + A(NP) 200 CONTINUE NOZ = NOZ + NX !follow with offset... 300 CONTINUE C C *** C do average calculations for first "center" line C *** C C average in the inner part TOP = 1 !init to lower left corner of matrix NOSUM = 1 !point for array of average values C IF (ALL.EQ.0) THEN DO 700 INDX=OFF+RX1,OFF+NXMRX DIFX = TOP !moving offset along line SUM(NOSUM)= 0.D0 C C sum over neighbourhood + update the row sums on the fly C [ if RX = 0, there should be no out-of-range problem, C just spill over into the next lines (RX = 0 and RY = 0 is not allowed!) ] C DO 500 NRL=1,NW2 SUM(NOSUM) = SUM(NOSUM) + V(NRL) V(NRL) = V(NRL) - A(DIFX) + A(DIFX+NW1) DIFX = DIFX + NX 500 CONTINUE C C replace pixel in question by average value TEST = (SUM(NOSUM) - A(INDX)) * FACTOR B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 700 CONTINUE C ELSE DO 900 INDX=OFF+RX1,OFF+NXMRX DIFX = TOP !moving offset along line SUM(NOSUM) = 0.D0 C DO 800 NRL=1,NW2 SUM(NOSUM) = SUM(NOSUM) + V(NRL) V(NRL) = V(NRL) - A(DIFX) + A(DIFX+NW1) DIFX = DIFX + NX 800 CONTINUE C TEST = SUM(NOSUM) * FACTOR B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 900 CONTINUE ENDIF C OFF = OFF + NX !update offset of center line C C *** C now handle average calculations for all following lines C *** C NO = 0 !pointer to elements which are out... DO 1600 NL=RY+2,NY-RY C C copy first row sums + set up the diff-row sum VNEW = 0.D0 NN = NO - NOZ DO 1200 NP=NOZ+1,NOZ+NW1 VNEW = VNEW + A(NP) - A(NP+NN) 1200 CONTINUE TOP = NOZ + 1 !init to lower left corner of matrix NOZ = NOZ + NX !follow with offset... C C average in the inner part NOSUM = 1 DIFX = NO + 1 C IF (ALL.EQ.0) THEN DO 1300 INDX=OFF+RX1,OFF+NXMRX-1 C C sum over neighbourhood + update the diff-row sum on the fly SUM(NOSUM) = SUM(NOSUM) + VNEW VNEW = VNEW - A(TOP) + A(TOP+NW1) + + A(DIFX) - A(DIFX+NW1) TEST = (SUM(NOSUM) - A(INDX)) * FACTOR C C replace pixel in question by average value B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 DIFX = DIFX + 1 NOSUM = NOSUM + 1 1300 CONTINUE C C handle last pixel on line specially to avoid out_of_bounds for (TOP+NW1) SUM(NOSUM) = SUM(NOSUM) + VNEW TEST = (SUM(NOSUM) - A(OFF+NXMRX)) * FACTOR ELSE DO 1400 INDX=OFF+RX1,OFF+NXMRX-1 C SUM(NOSUM) = SUM(NOSUM) + VNEW VNEW = VNEW - A(TOP) + A(TOP+NW1) + + A(DIFX) - A(DIFX+NW1) TEST = SUM(NOSUM) * FACTOR B(BINDX) = TEST C C replace pixel in question by average value BINDX = BINDX + 1 TOP = TOP + 1 DIFX = DIFX + 1 NOSUM = NOSUM + 1 1400 CONTINUE C C handle last pixel on line specially to avoid out_of_bounds for (TOP+NW1) SUM(NOSUM) = SUM(NOSUM) + VNEW TEST = SUM(NOSUM) * FACTOR ENDIF C B(BINDX) = TEST !for last pixel... BINDX = BINDX + 1 OFF = OFF + NX !update offset of center line NO = NO + NX C 1600 CONTINUE C C finally update the cuts CALL CHCUTS(B,BINDX,CUTS) RETURN C C check, if absolute or relative thresholding C 2000 IF (METH(2:2).NE.'A') GOTO 3000 C C C C -------------------------- C C here we handle boxcar filtering with absolute thresholding C C -------------------------- C C C C set up first row sums NOZ = 0 !offset for lower left C DO 2100 N=1,NW2 !we have NW2 row sums V(N) = 0.D0 DO 2060 NP=NOZ+1,NOZ+NW1 V(N) = V(N) + A(NP) 2060 CONTINUE NOZ = NOZ + NX !follow with offset... 2100 CONTINUE C C *** C do average calculations for first "center" line C *** C C average in the inner part TOP = 1 !init to lower left corner of matrix NOSUM = 1 !point for array of average values C IF (ALL.EQ.0) THEN DO 2200 INDX=OFF+RX1,OFF+NXMRX DIFX = TOP !moving offset along line W = A(INDX) SUM(NOSUM)= 0.D0 C DO 2180 NRL=1,NW2 SUM(NOSUM) = SUM(NOSUM) + V(NRL) V(NRL) = V(NRL) - A(DIFX) + A(DIFX+NW1) DIFX = DIFX + NX 2180 CONTINUE C C replace pixel in question by average value TEST = (SUM(NOSUM) - W) * FACTOR IF (ABS(TEST-W).LE.THRESH) TEST = W B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 2200 CONTINUE C ELSE DO 2300 INDX=OFF+RX1,OFF+NXMRX DIFX = TOP !moving offset along line W = A(INDX) SUM(NOSUM) = 0.0 C DO 2250 NRL=1,NW2 SUM(NOSUM) = SUM(NOSUM) + V(NRL) V(NRL) = V(NRL) - A(DIFX) + A(DIFX+NW1) DIFX = DIFX + NX 2250 CONTINUE C TEST = SUM(NOSUM) * FACTOR IF (ABS(TEST-W).LE.THRESH) TEST = W B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 2300 CONTINUE ENDIF C OFF = OFF + NX !update offset of center line C C *** C now handle average calculations for all following lines C *** C NO = 0 !offset for elements which are out... DO 2900 NL=RY+2,NY-RY C C copy first row sums + set up the diff-row sum VNEW = 0.D0 NN = NO - NOZ DO 2450 NP=NOZ+1,NOZ+NW1 VNEW = VNEW + A(NP) - A(NP+NN) 2450 CONTINUE TOP = NOZ + 1 !init to lower left corner of matrix NOZ = NOZ + NX !follow with offset... C C average in the inner part NOSUM = 1 DIFX = NO + 1 C IF (ALL.EQ.0) THEN DO 2600 INDX=OFF+RX1,OFF+NXMRX-1 W = A(INDX) SUM(NOSUM) = SUM(NOSUM) + VNEW VNEW = VNEW - A(TOP) + A(TOP+NW1) + + A(DIFX) - A(DIFX+NW1) TEST = (SUM(NOSUM) - W) * FACTOR C C replace pixel in question by average value IF (ABS(TEST-W).LE.THRESH) TEST = W B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 DIFX = DIFX + 1 2600 CONTINUE C C handle last pixel on line specially to avoid out_of_bounds for (TOP+NW1) W = A(OFF+NXMRX) SUM(NOSUM) = SUM(NOSUM) + VNEW TEST = (SUM(NOSUM) - W) * FACTOR IF (ABS(TEST-W).LE.THRESH) TEST = W ELSE DO 2700 INDX=OFF+RX1,OFF+NXMRX-1 W = A(INDX) SUM(NOSUM) = SUM(NOSUM) + VNEW VNEW = VNEW - A(TOP) + A(TOP+NW1) + + A(DIFX) - A(DIFX+NW1) TEST = SUM(NOSUM) * FACTOR C C replace pixel in question by average value IF (ABS(TEST-W).LE.THRESH) TEST = W B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 DIFX = DIFX + 1 2700 CONTINUE C C handle last pixel on line specially to avoid out_of_bounds for (TOP+NW1) W = A(OFF+NXMRX) SUM(NOSUM) = SUM(NOSUM) + VNEW TEST = SUM(NOSUM) * FACTOR IF (ABS(TEST-W).LE.THRESH) TEST = W ENDIF C B(BINDX) = TEST !for last pixel... BINDX = BINDX + 1 OFF = OFF + NX !update offset of center line NO = NO + NX C 2900 CONTINUE C C finally update the cuts CALL CHCUTS(B,BINDX,CUTS) RETURN C C C -------------------------- C C here we handle boxcar filtering with relative thresholding C C -------------------------- C C C set up first row sums 3000 NOZ = 0 !offset for lower left C DO 3100 N=1,NW2 !we have NW2 row sums V(N) = 0.D0 DO 3060 NP=NOZ+1,NOZ+NW1 V(N) = V(N) + A(NP) 3060 CONTINUE NOZ = NOZ + NX !follow with offset... 3100 CONTINUE C C *** C do average calculations for first "center" line C *** C C average in the inner part TOP = 1 !init to lower left corner of matrix NOSUM = 1 !point for array of average values C IF (ALL.EQ.0) THEN DO 3300 INDX=OFF+RX1,OFF+NXMRX W = A(INDX) DIFX = TOP !moving offset along line SUM(NOSUM)= 0.D0 C DO 3200 NRL=1,NW2 SUM(NOSUM) = SUM(NOSUM) + V(NRL) V(NRL) = V(NRL) - A(DIFX) + A(DIFX+NW1) DIFX = DIFX + NX 3200 CONTINUE C C replace pixel in question by average value TEST = (SUM(NOSUM) - W) * FACTOR IF (ABS(TEST-W).LE.ABS(THRESH*W)) TEST = W B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 3300 CONTINUE C ELSE DO 3500 INDX=OFF+RX1,OFF+NXMRX W = A(INDX) DIFX = TOP !moving offset along line SUM(NOSUM) = 0.D0 C DO 3400 NRL=1,NW2 SUM(NOSUM) = SUM(NOSUM) + V(NRL) V(NRL) = V(NRL) - A(DIFX) + A(DIFX+NW1) DIFX = DIFX + NX 3400 CONTINUE C TEST = SUM(NOSUM) * FACTOR IF (ABS(TEST-W).LE.ABS(THRESH*W)) TEST = W B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 3500 CONTINUE ENDIF C OFF = OFF + NX !update offset of center line C C *** C now handle average calculations for all following lines C *** C NO = 0 !offset for elements which are out... DO 3900 NL=RY+2,NY-RY C C copy first row sums + set up the diff-row sum VNEW = 0.D0 NN = NO - NOZ DO 3650 NP=NOZ+1,NOZ+NW1 VNEW = VNEW + A(NP) - A(NP+NN) 3650 CONTINUE TOP = NOZ + 1 !init to lower left corner of matrix NOZ = NOZ + NX !follow with offset... C C average in the inner part NOSUM = 1 DIFX = NO + 1 C IF (ALL.EQ.0) THEN DO 3700 INDX=OFF+RX1,OFF+NXMRX-1 W = A(INDX) SUM(NOSUM) = SUM(NOSUM) + VNEW VNEW = VNEW - A(TOP) + A(TOP+NW1) + + A(DIFX) - A(DIFX+NW1) TEST = (SUM(NOSUM) - W) * FACTOR C C replace pixel in question by average value IF (ABS(TEST-W).LE.ABS(THRESH*W)) TEST = W B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 DIFX = DIFX + 1 3700 CONTINUE C C handle last pixel on line specially to avoid out_of_bounds for (TOP+NW1) W = A(OFF+NXMRX) SUM(NOSUM) = SUM(NOSUM) + VNEW TEST = (SUM(NOSUM) - W) * FACTOR IF (ABS(TEST-W).LE.ABS(THRESH*W)) TEST = W ELSE DO 3800 INDX=OFF+RX1,OFF+NXMRX-1 W = A(INDX) SUM(NOSUM) = SUM(NOSUM) + VNEW VNEW = VNEW - A(TOP) + A(TOP+NW1) + + A(DIFX) - A(DIFX+NW1) TEST = SUM(NOSUM) * FACTOR C C replace pixel in question by average value IF (ABS(TEST-W).LE.ABS(THRESH*W)) TEST = W B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 DIFX = DIFX + 1 3800 CONTINUE C C handle last pixel on line specially to avoid out_of_bounds for (TOP+NW1) W = A(OFF+NXMRX) SUM(NOSUM) = SUM(NOSUM) + VNEW TEST = SUM(NOSUM) * FACTOR IF (ABS(TEST-W).LE.ABS(THRESH*W)) TEST = W ENDIF C B(BINDX) = TEST !for last pixel... BINDX = BINDX + 1 OFF = OFF + NX !update offset of center line NO = NO + NX C 3900 CONTINUE C C finally update the cuts CALL CHCUTS(B,BINDX,CUTS) RETURN C C C -------------------------- C C here we handle boxcar filtering with fixed thresholding interval C C -------------------------- C C C set up first row sums 4000 NOZ = 0 !offset for lower left C DO 4100 N=1,NW2 !we have NW2 row sums V(N) = 0.D0 DO 4060 NP=NOZ+1,NOZ+NW1 V(N) = V(N) + A(NP) 4060 CONTINUE NOZ = NOZ + NX !follow with offset... 4100 CONTINUE C C *** C do average calculations for first "center" line C *** C C average in the inner part TOP = 1 !init to lower left corner of matrix NOSUM = 1 !point for array of average values C IF (ALL.EQ.0) THEN DO 4300 INDX=OFF+RX1,OFF+NXMRX TEST = A(INDX) DIFX = TOP !moving offset along line SUM(NOSUM)= 0.D0 C DO 4200 NRL=1,NW2 SUM(NOSUM) = SUM(NOSUM) + V(NRL) V(NRL) = V(NRL) - A(DIFX) + A(DIFX+NW1) DIFX = DIFX + NX 4200 CONTINUE C C replace pixel in question by average value IF ( (TEST.GE.TR(1)) .AND. (TEST.LE.TR(2)) ) + TEST = (SUM(NOSUM) - TEST) * FACTOR B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 4300 CONTINUE C ELSE DO 4500 INDX=OFF+RX1,OFF+NXMRX TEST = A(INDX) DIFX = TOP !moving offset along line SUM(NOSUM) = 0.D0 C DO 4400 NRL=1,NW2 SUM(NOSUM) = SUM(NOSUM) + V(NRL) V(NRL) = V(NRL) - A(DIFX) + A(DIFX+NW1) DIFX = DIFX + NX 4400 CONTINUE C IF ( (TEST.GE.TR(1)) .AND. (TEST.LE.TR(2)) ) + TEST = SUM(NOSUM) * FACTOR B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 4500 CONTINUE ENDIF C OFF = OFF + NX !update offset of center line C C *** C now handle average calculations for all following lines C *** C NO = 0 !offset for elements which are out... DO 4900 NL=RY+2,NY-RY C C copy first row sums + set up the diff-row sum VNEW = 0.D0 NN = NO - NOZ DO 4650 NP=NOZ+1,NOZ+NW1 VNEW = VNEW + A(NP) - A(NP+NN) 4650 CONTINUE TOP = NOZ + 1 !init to lower left corner of matrix NOZ = NOZ + NX !follow with offset... C C average in the inner part NOSUM = 1 DIFX = NO + 1 C IF (ALL.EQ.0) THEN DO 4700 INDX=OFF+RX1,OFF+NXMRX-1 TEST = A(INDX) SUM(NOSUM) = SUM(NOSUM) + VNEW VNEW = VNEW - A(TOP) + A(TOP+NW1) + + A(DIFX) - A(DIFX+NW1) C C replace pixel in question by average value IF ( (TEST.GE.TR(1)) .AND. (TEST.LE.TR(2)) ) + TEST = (SUM(NOSUM) - TEST) * FACTOR B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 DIFX = DIFX + 1 4700 CONTINUE C C handle last pixel on line specially to avoid out_of_bounds for (TOP+NW1) TEST = A(OFF+NXMRX) SUM(NOSUM) = SUM(NOSUM) + VNEW IF ( (TEST.GE.TR(1)) .AND. (TEST.LE.TR(2)) ) + TEST = (SUM(NOSUM) - TEST) * FACTOR ELSE DO 4800 INDX=OFF+RX1,OFF+NXMRX-1 TEST = A(INDX) SUM(NOSUM) = SUM(NOSUM) + VNEW VNEW = VNEW - A(TOP) + A(TOP+NW1) + + A(DIFX) - A(DIFX+NW1) C C replace pixel in question by average value IF ( (TEST.GE.TR(1)) .AND. (TEST.LE.TR(2)) ) + TEST = SUM(NOSUM) * FACTOR B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 DIFX = DIFX + 1 4800 CONTINUE C C handle last pixel on line specially to avoid out_of_bounds for (TOP+NW1) TEST = A(OFF+NXMRX) SUM(NOSUM) = SUM(NOSUM) + VNEW IF ( (TEST.GE.TR(1)) .AND. (TEST.LE.TR(2)) ) + TEST = SUM(NOSUM) * FACTOR ENDIF C B(BINDX) = TEST !for last pixel... BINDX = BINDX + 1 OFF = OFF + NX !update offset of center line NO = NO + NX C 4900 CONTINUE C C finally update the cuts CALL CHCUTS(B,BINDX,CUTS) RETURN C END SUBROUTINE BVER2A(METH,A,B,NX,NY,RADIUS,TR,ALL,V,SUM,CUTS) C IMPLICIT NONE C INTEGER NX,NY,ALL INTEGER DIFX,ICOUNT INTEGER INDX,BINDX INTEGER N,NN,NL,NO,NOZ INTEGER NP,NRL,NXMRX INTEGER OFF,RX,RX1,RXM1,RY,RYNX,TOP INTEGER RADIUS(2),NW1,NW2,NOSUM C REAL A(*),B(*),TR(*),CUTS(*) REAL W,TEST,COUNT,FACTOR,THRESH C DOUBLE PRECISION V(*),SUM(*) DOUBLE PRECISION VNEW C CHARACTER*(*) METH C C initialize C BINDX = 1 RX = RADIUS(1) RY = RADIUS(2) NW1 = RX*2 + 1 !no. of values on row of filter NW2 = RY*2 + 1 !no. of values on column of filter ICOUNT = NW1 * NW2 COUNT = FLOAT(ICOUNT) !no. of pixels in neighbourhood with center C RYNX = RY * NX NXMRX = NX - RX RX1 = RX + 1 RXM1 = RX - 1 OFF = RYNX !offset for center row C IF (ALL.EQ.0) THEN FACTOR = 1. / (COUNT-1.) ELSE FACTOR = 1. / COUNT ENDIF C C fork according to test criterion (for speed reasons) IF (METH(2:2).EQ.'F') THEN GOTO 4000 ELSE THRESH = TR(1) ENDIF IF (THRESH.GT.0.) GOTO 2000 C C C -------------------------- C C here we handle boxcar filtering without thresholding C C -------------------------- C C C set up first row sums NOZ = 0 !offset for lower left C DO 300 N=1,NW2 !we have NW2 row sums V(N) = 0.D0 DO 200 NP=NOZ+1,NOZ+NW1 V(N) = V(N) + A(NP) 200 CONTINUE NOZ = NOZ + NX !follow with offset... 300 CONTINUE C C *** C do average calculations for first "center" line C *** C C average in the inner part TOP = 1 !init to lower left corner of matrix NOSUM = 1 !point for array of average values C IF (ALL.EQ.0) THEN DO 700 INDX=OFF+RX1,OFF+NXMRX DIFX = TOP !moving offset along line SUM(NOSUM)= 0.D0 C C sum over neighbourhood + update the row sums on the fly C [ if RX = 0, there should be no out-of-range problem, C just spill over into the next lines (RX = 0 and RY = 0 is not allowed!) ] C DO 500 NRL=1,NW2 SUM(NOSUM) = SUM(NOSUM) + V(NRL) V(NRL) = V(NRL) - A(DIFX) + A(DIFX+NW1) DIFX = DIFX + NX 500 CONTINUE C C replace pixel in question by average value TEST = (SUM(NOSUM) - A(INDX)) * FACTOR B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 700 CONTINUE C ELSE DO 900 INDX=OFF+RX1,OFF+NXMRX DIFX = TOP !moving offset along line SUM(NOSUM) = 0.D0 C DO 800 NRL=1,NW2 SUM(NOSUM) = SUM(NOSUM) + V(NRL) V(NRL) = V(NRL) - A(DIFX) + A(DIFX+NW1) DIFX = DIFX + NX 800 CONTINUE C TEST = SUM(NOSUM) * FACTOR B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 900 CONTINUE ENDIF C OFF = OFF + NX !update offset of center line C C *** C now handle average calculations for all following lines C *** C NO = 0 !pointer to elements which are out... DO 1600 NL=RY+2,NY-RY C C copy first row sums + set up the diff-row sum VNEW = 0.D0 NN = NO - NOZ DO 1200 NP=NOZ+1,NOZ+NW1 VNEW = VNEW + A(NP) - A(NP+NN) 1200 CONTINUE TOP = NOZ + 1 !init to lower left corner of matrix NOZ = NOZ + NX !follow with offset... C C average in the inner part NOSUM = 1 DIFX = NO + 1 C IF (ALL.EQ.0) THEN DO 1300 INDX=OFF+RX1,OFF+NXMRX C C sum over neighbourhood + update the diff-row sum on the fly SUM(NOSUM) = SUM(NOSUM) + VNEW VNEW = VNEW - A(TOP) + A(TOP+NW1) + + A(DIFX) - A(DIFX+NW1) TEST = (SUM(NOSUM) - A(INDX)) * FACTOR C C replace pixel in question by average value B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 DIFX = DIFX + 1 NOSUM = NOSUM + 1 1300 CONTINUE C ELSE DO 1400 INDX=OFF+RX1,OFF+NXMRX C SUM(NOSUM) = SUM(NOSUM) + VNEW VNEW = VNEW - A(TOP) + A(TOP+NW1) + + A(DIFX) - A(DIFX+NW1) TEST = SUM(NOSUM) * FACTOR B(BINDX) = TEST BINDX = BINDX + 1 C C replace pixel in question by average value TOP = TOP + 1 DIFX = DIFX + 1 NOSUM = NOSUM + 1 1400 CONTINUE ENDIF C OFF = OFF + NX !update offset of center line NO = NO + NX C 1600 CONTINUE C C finally update the cuts CALL CHCUTS(B,BINDX,CUTS) RETURN C C C check, if absolute or relative thresholding C 2000 IF (METH(2:2).NE.'A') GOTO 3000 C C C C -------------------------- C C here we handle boxcar filtering with absolute thresholding C C -------------------------- C C C C set up first row sums NOZ = 0 !offset for lower left C DO 2100 N=1,NW2 !we have NW2 row sums V(N) = 0.D0 DO 2060 NP=NOZ+1,NOZ+NW1 V(N) = V(N) + A(NP) 2060 CONTINUE NOZ = NOZ + NX !follow with offset... 2100 CONTINUE C C *** C do average calculations for first "center" line C *** C C average in the inner part TOP = 1 !init to lower left corner of matrix NOSUM = 1 !point for array of average values C IF (ALL.EQ.0) THEN DO 2200 INDX=OFF+RX1,OFF+NXMRX DIFX = TOP !moving offset along line W = A(INDX) SUM(NOSUM)= 0.D0 C DO 2180 NRL=1,NW2 SUM(NOSUM) = SUM(NOSUM) + V(NRL) V(NRL) = V(NRL) - A(DIFX) + A(DIFX+NW1) DIFX = DIFX + NX 2180 CONTINUE C C replace pixel in question by average value TEST = (SUM(NOSUM) - W) * FACTOR IF (ABS(TEST-W).LE.THRESH) TEST = W B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 2200 CONTINUE C ELSE DO 2300 INDX=OFF+RX1,OFF+NXMRX DIFX = TOP !moving offset along line W = A(INDX) SUM(NOSUM) = 0.D0 C DO 2250 NRL=1,NW2 SUM(NOSUM) = SUM(NOSUM) + V(NRL) V(NRL) = V(NRL) - A(DIFX) + A(DIFX+NW1) DIFX = DIFX + NX 2250 CONTINUE C TEST = SUM(NOSUM) * FACTOR IF (ABS(TEST-W).LE.THRESH) TEST = W B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 2300 CONTINUE ENDIF C OFF = OFF + NX !update offset of center line C C *** C now handle average calculations for all following lines C *** C NO = 0 !offset for elements which are out... DO 2900 NL=RY+2,NY-RY C C copy first row sums + set up the diff-row sum VNEW = 0.D0 NN = NO - NOZ DO 2450 NP=NOZ+1,NOZ+NW1 VNEW = VNEW + A(NP) - A(NP+NN) 2450 CONTINUE TOP = NOZ + 1 !init to lower left corner of matrix NOZ = NOZ + NX !follow with offset... C C average in the inner part NOSUM = 1 DIFX = NO + 1 C IF (ALL.EQ.0) THEN DO 2600 INDX=OFF+RX1,OFF+NXMRX W = A(INDX) SUM(NOSUM) = SUM(NOSUM) + VNEW VNEW = VNEW - A(TOP) + A(TOP+NW1) + + A(DIFX) - A(DIFX+NW1) TEST = (SUM(NOSUM) - W) * FACTOR C C replace pixel in question by average value IF (ABS(TEST-W).LE.THRESH) TEST = W B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 DIFX = DIFX + 1 2600 CONTINUE C ELSE DO 2700 INDX=OFF+RX1,OFF+NXMRX W = A(INDX) SUM(NOSUM) = SUM(NOSUM) + VNEW VNEW = VNEW - A(TOP) + A(TOP+NW1) + + A(DIFX) - A(DIFX+NW1) TEST = SUM(NOSUM) * FACTOR C C replace pixel in question by average value IF (ABS(TEST-W).LE.THRESH) TEST = W B(BINDX) = TEST TOP = TOP + 1 NOSUM = NOSUM + 1 DIFX = DIFX + 1 2700 CONTINUE ENDIF C OFF = OFF + NX !update offset of center line NO = NO + NX C 2900 CONTINUE C C finally update the cuts CALL CHCUTS(B,BINDX,CUTS) RETURN C C C -------------------------- C C here we handle boxcar filtering with relative thresholding C C -------------------------- C C C set up first row sums 3000 NOZ = 0 !offset for lower left C DO 3100 N=1,NW2 !we have NW2 row sums V(N) = 0.D0 DO 3060 NP=NOZ+1,NOZ+NW1 V(N) = V(N) + A(NP) 3060 CONTINUE NOZ = NOZ + NX !follow with offset... 3100 CONTINUE C C *** C do average calculations for first "center" line C *** C C average in the inner part TOP = 1 !init to lower left corner of matrix NOSUM = 1 !point for array of average values C IF (ALL.EQ.0) THEN DO 3300 INDX=OFF+RX1,OFF+NXMRX W = A(INDX) DIFX = TOP !moving offset along line SUM(NOSUM)= 0.D0 C DO 3200 NRL=1,NW2 SUM(NOSUM) = SUM(NOSUM) + V(NRL) V(NRL) = V(NRL) - A(DIFX) + A(DIFX+NW1) DIFX = DIFX + NX 3200 CONTINUE C C replace pixel in question by average value TEST = (SUM(NOSUM) - W) * FACTOR IF (ABS(TEST-W).LE.ABS(THRESH*TEST)) TEST = W B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 3300 CONTINUE C ELSE DO 3500 INDX=OFF+RX1,OFF+NXMRX W = A(INDX) DIFX = TOP !moving offset along line SUM(NOSUM) = 0.D0 C DO 3400 NRL=1,NW2 SUM(NOSUM) = SUM(NOSUM) + V(NRL) V(NRL) = V(NRL) - A(DIFX) + A(DIFX+NW1) DIFX = DIFX + NX 3400 CONTINUE C TEST = SUM(NOSUM) * FACTOR IF (ABS(TEST-W).LE.ABS(THRESH*TEST)) TEST = W B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 3500 CONTINUE ENDIF C OFF = OFF + NX !update offset of center line C C *** C now handle average calculations for all following lines C *** C NO = 0 !offset for elements which are out... DO 3900 NL=RY+2,NY-RY C C copy first row sums + set up the diff-row sum VNEW = 0.D0 NN = NO - NOZ DO 3650 NP=NOZ+1,NOZ+NW1 VNEW = VNEW + A(NP) - A(NP+NN) 3650 CONTINUE TOP = NOZ + 1 !init to lower left corner of matrix NOZ = NOZ + NX !follow with offset... C C average in the inner part NOSUM = 1 DIFX = NO + 1 C IF (ALL.EQ.0) THEN DO 3700 INDX=OFF+RX1,OFF+NXMRX W = A(INDX) SUM(NOSUM) = SUM(NOSUM) + VNEW VNEW = VNEW - A(TOP) + A(TOP+NW1) + + A(DIFX) - A(DIFX+NW1) TEST = (SUM(NOSUM) - W) * FACTOR C C replace pixel in question by average value IF (ABS(TEST-W).LE.ABS(THRESH*TEST)) TEST = W B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 DIFX = DIFX + 1 3700 CONTINUE C ELSE DO 3800 INDX=OFF+RX1,OFF+NXMRX W = A(INDX) SUM(NOSUM) = SUM(NOSUM) + VNEW VNEW = VNEW - A(TOP) + A(TOP+NW1) + + A(DIFX) - A(DIFX+NW1) TEST = SUM(NOSUM) * FACTOR C C replace pixel in question by average value IF (ABS(TEST-W).LE.ABS(THRESH*TEST)) TEST = W B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 DIFX = DIFX + 1 3800 CONTINUE ENDIF C OFF = OFF + NX !update offset of center line NO = NO + NX C 3900 CONTINUE C C finally update the cuts CALL CHCUTS(B,BINDX,CUTS) RETURN C C C C -------------------------- C C here we handle boxcar filtering with fixed thresholding interval C C -------------------------- C C C set up first row sums 4000 NOZ = 0 !offset for lower left C DO 4100 N=1,NW2 !we have NW2 row sums V(N) = 0.D0 DO 4060 NP=NOZ+1,NOZ+NW1 V(N) = V(N) + A(NP) 4060 CONTINUE NOZ = NOZ + NX !follow with offset... 4100 CONTINUE C C *** C do average calculations for first "center" line C *** C C average in the inner part TOP = 1 !init to lower left corner of matrix NOSUM = 1 !point for array of average values C IF (ALL.EQ.0) THEN DO 4300 INDX=OFF+RX1,OFF+NXMRX TEST = A(INDX) DIFX = TOP !moving offset along line SUM(NOSUM)= 0.D0 C DO 4200 NRL=1,NW2 SUM(NOSUM) = SUM(NOSUM) + V(NRL) V(NRL) = V(NRL) - A(DIFX) + A(DIFX+NW1) DIFX = DIFX + NX 4200 CONTINUE C C replace pixel in question by average value IF ( (TEST.GE.TR(1)) .AND. (TEST.LE.TR(2)) ) + TEST = (SUM(NOSUM) - TEST) * FACTOR B(BINDX) = TEST TOP = TOP + 1 NOSUM = NOSUM + 1 4300 CONTINUE C ELSE DO 4500 INDX=OFF+RX1,OFF+NXMRX TEST = A(INDX) DIFX = TOP !moving offset along line SUM(NOSUM) = 0.D0 C DO 4400 NRL=1,NW2 SUM(NOSUM) = SUM(NOSUM) + V(NRL) V(NRL) = V(NRL) - A(DIFX) + A(DIFX+NW1) DIFX = DIFX + NX 4400 CONTINUE C IF ( (TEST.GE.TR(1)) .AND. (TEST.LE.TR(2)) ) + TEST = SUM(NOSUM) * FACTOR B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 4500 CONTINUE ENDIF C OFF = OFF + NX !update offset of center line C C *** C now handle average calculations for all following lines C *** C NO = 0 !offset for elements which are out... DO 4900 NL=RY+2,NY-RY C C copy first row sums + set up the diff-row sum VNEW = 0.D0 NN = NO - NOZ DO 4650 NP=NOZ+1,NOZ+NW1 VNEW = VNEW + A(NP) - A(NP+NN) 4650 CONTINUE TOP = NOZ + 1 !init to lower left corner of matrix NOZ = NOZ + NX !follow with offset... C C average in the inner part NOSUM = 1 DIFX = NO + 1 C IF (ALL.EQ.0) THEN DO 4700 INDX=OFF+RX1,OFF+NXMRX TEST = A(INDX) SUM(NOSUM) = SUM(NOSUM) + VNEW VNEW = VNEW - A(TOP) + A(TOP+NW1) + + A(DIFX) - A(DIFX+NW1) C C replace pixel in question by average value IF ( (TEST.GE.TR(1)) .AND. (TEST.LE.TR(2)) ) + TEST = (SUM(NOSUM) - TEST) * FACTOR B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 DIFX = DIFX + 1 4700 CONTINUE C ELSE DO 4800 INDX=OFF+RX1,OFF+NXMRX TEST = A(INDX) SUM(NOSUM) = SUM(NOSUM) + VNEW VNEW = VNEW - A(TOP) + A(TOP+NW1) + + A(DIFX) - A(DIFX+NW1) C C replace pixel in question by average value IF ( (TEST.GE.TR(1)) .AND. (TEST.LE.TR(2)) ) + TEST= SUM(NOSUM) * FACTOR B(BINDX) = TEST BINDX = BINDX + 1 TOP = TOP + 1 NOSUM = NOSUM + 1 DIFX = DIFX + 1 4800 CONTINUE ENDIF C OFF = OFF + NX !update offset of center line NO = NO + NX C 4900 CONTINUE C C finally update the cuts CALL CHCUTS(B,BINDX,CUTS) RETURN C END