C @(#)genzz1.for 13.1.1.3 (ESO-DMD) 02/12/99 19:14:53 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 Massachusetts Ave, Cambridge, C MA 02139, USA. C C Correspondence 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 GENZZ1 C C++++++++++++++++++++++++++++++++++++++++++++++++++ C C.IDENTIFICATION C program GENZZ1 version 1.00 890705 C K. Banse ESO - Garching C 1.10 910204 C C.KEYWORDS C general image processing functions C C.PURPOSE C merge FNDPIX.FOR, POLYGO.FOR, ROTATE.FOR, ROTGEN.FOR C SCALE.FOR, WEIGHT.FOR C C.ALGORITHM C straight forward C C.INPUT/OUTPUT C the following keys are used: C C ACTION/C/1/2 option, C = FN, for FNDPIX C = PO, for POLYGO C C.VERSIONS C 010410 last modif C C-------------------------------------------------- C IMPLICIT NONE C INTEGER IAV,STAT INTEGER UNIT(1),NULLO C CHARACTER ACTION*2 C C get into MIDAS CALL STSPRO('GENZZ1') C C get action flag CALL STKRDC('ACTION',1,1,2,IAV,ACTION,UNIT,NULLO,STAT) C IF (ACTION.EQ.'FN') THEN CALL SUBFNP ELSE IF (ACTION.EQ.'PO') THEN CALL SUBPOL ELSE CALL STETER(11,'Invalid task for module "genzz1.for"...') ENDIF C C That's it folks... CALL STSEPI C END SUBROUTINE SUBFNP C C++++++++++++++++++++++++++++++++++++++++++++++++++ C C.IDENTIFICATION C subroutine SUBFNP version 1.40 850822 C K. Banse ESO - Garching C 1.50 890123 1.60 891124 C C.KEYWORDS C bulk data frame, pixel range C C.PURPOSE C 1) find 1. pixel which falls inside or outside given interval C 2) find pixel no's. of min + max C C.ALGORITHM C straight forward C C.INPUT/OUTPUT C the following keys are used: C C IN_A/C/1/60 name of input frame C INPUTR/R/1/2 XL, XH for interval [XL,XH] C ACTION/C/3/4 IN or OUT, to indicate search for inside or C outside of given interval C or MI, if we search for min + max C OUTPUTI/I/1/2 x-, y-coordinate of 1. pixel which C matches, C = 0,0 if no match C OUTPUTR/R/1/1 value, if matching pixel found C or C OUTPUTI/I/1/4 x-, y-coords. of minimum and maximum C OUTPUTR/R/1/2 min + max value, C also stored into descriptor LHCUTS/R/3/2 C OUTPUTR/R/3/4 world coords. of OUTPUTI C C.VERSIONS C 1.40 creation C 1.50 move to FORTRAN 77 + new ST interfaces C and add table option to record all pixels in/out of interval C 1.60 really put 0,0 into OUTPUTI if nothing found ... C-------------------------------------------------- C IMPLICIT NONE C INTEGER IAV,IDUM INTEGER N,NAXIS,NPIX(3),SPIX(3),EPIX(3) INTEGER*8 PNTR INTEGER IMNO,STAT INTEGER UNIT(1),NULLO,MADRID(1) INTEGER TID,TCOLS(4),NCOLS,ROWMAX,XRFLAG INTEGER RETBUF(3),CHUNK C REAL XRNG(2) C DOUBLE PRECISION DDUM C CHARACTER FRAME*60,ACTION*4,TABLE*72,AREA*80 CHARACTER TLABL(4)*16,TUNIT(4)*16,CXRNG*60 C INCLUDE 'MID_INCLUDE:ST_DEF.INC' C COMMON /VMR/ MADRID C DATA TLABL /'X_PIX ','Y_PIX ','Z_PIX ','VALUE '/ DATA TUNIT /'FRAME_PIXEL','FRAME_PIXEL','FRAME_PIXEL',' '/ DATA XRNG /0.0,0.0/ DATA NPIX /1,1,1/, SPIX /1,1,1/, EPIX /1,1,1/ C INCLUDE 'MID_INCLUDE:ST_DAT.INC' C CALL STPSET(F_FITS_PARM,1,STAT) C C get name of input image, interval + action flags CALL STKRDC('IN_A',1,1,60,IAV,FRAME,UNIT,NULLO,STAT) CALL STKRDC('P2',1,1,60,IAV,CXRNG,UNIT,NULLO,STAT) CALL STKRDC('ACTION',1,3,4,IAV,ACTION,UNIT,NULLO,STAT) CALL UPCAS(ACTION,ACTION) C C find max. no. of lines we can use CALL STFOPN(FRAME,D_R4_FORMAT,0,F_IMA_TYPE,IMNO,STAT) CALL STDRDI(IMNO,'NAXIS',1,1,IAV,NAXIS,UNIT,NULLO,STAT) CALL STDRDI(IMNO,'NPIX',1,3,IAV,NPIX,UNIT,NULLO,STAT) CALL HACKUP(NPIX(1),D_R4_FORMAT,RETBUF) C CHUNK = RETBUF(1)*NPIX(1) !size of memory to allocate CALL STFXMP(CHUNK,D_R4_FORMAT,PNTR,STAT) C C handle FIND/MINMAX C IF (ACTION(1:2).EQ.'MI') THEN CALL FNDMX(RETBUF,IMNO,MADRID(PNTR),NAXIS,NPIX) ELSE C C handle FIND/PIXEL C XRFLAG = 0 !check interval IF (CXRNG(1:2).EQ.'<,') THEN IF (CXRNG(3:3).EQ.'>') THEN XRFLAG = 3 ! <,> ELSE XRFLAG = 1 ! <,m CALL GENCNV(CXRNG(3:),2,1,IDUM,XRNG(2),DDUM,N) IF (N.NE.1) CALL STETER(3,'Bad syntax for interval...') ENDIF ELSE N = INDEX(CXRNG,',') IF (N.LT.2) CALL STETER(3,'Bad syntax for interval...') N = N + 1 IF (CXRNG(N:N).EQ.'>') THEN XRFLAG = 2 ! n,> CALL GENCNV(CXRNG,2,1,IDUM,XRNG(1),DDUM,N) IF (N.NE.1) CALL STETER(3,'Bad syntax for interval...') ELSE CALL GENCNV(CXRNG,2,2,IDUM,XRNG,DDUM,N) IF (N.NE.2) CALL STETER(3,'Bad syntax for interval...') ENDIF ENDIF C CALL STKRDC('P7',1,1,60,IAV,AREA,UNIT,NULLO,STAT) IF (AREA(1:1).NE.'?') !check sub-area option + CALL EXTCOO(IMNO,AREA,3,IAV,SPIX,EPIX,STAT) C IF (ACTION(4:4).NE.'A') THEN !only look for first one IF (AREA(1:1).NE.'?') THEN CALL FNDPXA(RETBUF,IMNO,MADRID(PNTR),NAXIS,NPIX, + ACTION(1:3),XRNG,XRFLAG,CXRNG,SPIX,EPIX) ELSE CALL FNDPX(RETBUF,IMNO,MADRID(PNTR),NAXIS,NPIX, + ACTION(1:3),XRNG,XRFLAG,CXRNG) ENDIF C ELSE !look for all pixels CALL STKRDI('INPUTI',1,1,IAV,ROWMAX,UNIT,NULLO,STAT) CALL STKRDC('IN_B',1,1,72,IAV,TABLE,UNIT,NULLO,STAT) N = INDEX(TABLE,',') IF (TABLE(1:1).EQ.'+') THEN TID = -99 !no table or image output IF (AREA(1:1).NE.'?') THEN CALL FNDPXS(RETBUF,IMNO,MADRID(PNTR),NAXIS,NPIX, + ACTION(1:3),XRNG,XRFLAG,CXRNG, + TID,TCOLS,ROWMAX,SPIX,EPIX) ELSE CALL FNDPXT(RETBUF,IMNO,MADRID(PNTR),NAXIS,NPIX, + ACTION(1:3),XRNG,XRFLAG,CXRNG, + TID,TCOLS,ROWMAX) ENDIF C ELSE IF (N.GT.0) THEN !we have name,image TABLE(N:) = ' ' N = NPIX(1)*NPIX(2)*NPIX(3) CALL STFCRE('middummprov',D_R4_FORMAT,F_O_MODE, + F_IMA_TYPE,N,TID,STAT) TID = -TID !so the subroutine knows IF (AREA(1:1).NE.'?') THEN CALL FNDPXS(RETBUF,IMNO,MADRID(PNTR),NAXIS,NPIX, + ACTION(1:3),XRNG,XRFLAG,CXRNG, + TID,TCOLS,ROWMAX,SPIX,EPIX) ELSE CALL FNDPXT(RETBUF,IMNO,MADRID(PNTR),NAXIS,NPIX, + ACTION(1:3),XRNG,XRFLAG,CXRNG, + TID,TCOLS,ROWMAX) ENDIF CALL STFCLO(IMNO,STAT) !close input frame IF (RETBUF(2).LT.1) THEN !no pixel there... CALL STFCLO(TID,STAT) CALL STFDEL('middummprov.bdf',STAT) CALL STTPUT('no output file created',STAT) RETURN ENDIF IF (ROWMAX.NE.0) THEN IF (ROWMAX.LT.0) ROWMAX = -ROWMAX IF (RETBUF(2).GT.ROWMAX) RETBUF(2) = ROWMAX ENDIF C C if more than 5 percent different in size, crop out image IDUM = (N / 100) * 5 IF ((N-RETBUF(2)).GT.IDUM) THEN NPIX(1) = RETBUF(2) !real no. of pixels CALL STFCRE(TABLE,D_R4_FORMAT,F_O_MODE,F_IMA_TYPE, + RETBUF(2),IMNO,STAT) IF (NPIX(1).LE.CHUNK) THEN CALL STFGET(TID,1,NPIX(1),IAV,MADRID(PNTR),STAT) CALL STFPUT(IMNO,1,NPIX(1),MADRID(PNTR),STAT) ELSE IDUM = NPIX(1)/CHUNK TCOLS(1) = 1 DO 770, N=1,IDUM CALL STFGET(TID,TCOLS(1),CHUNK,IAV,MADRID(PNTR), + STAT) CALL STFPUT(IMNO,TCOLS(1),CHUNK,MADRID(PNTR), + STAT) TCOLS(1) = TCOLS(1) + CHUNK 770 CONTINUE IDUM = NPIX(1) - (IDUM*CHUNK) !look for remainder IF (IDUM.GT.0) THEN CALL STFGET(TID,TCOLS(1),IDUM,IAV,MADRID(PNTR), + STAT) CALL STFPUT(IMNO,TCOLS(1),IDUM,MADRID(PNTR), + STAT) ENDIF ENDIF CALL STFCLO(TID,STAT) CALL STFDEL('middummprov.bdf',STAT) TID = -99 !switch for later ELSE NPIX(1) = N IMNO = TID ENDIF TCOLS(1) = 1 !use 1-dim result image CALL STDWRI(IMNO,'NAXIS',TCOLS,1,1,UNIT,STAT) CALL STDWRI(IMNO,'NPIX',NPIX,1,1,UNIT,STAT) DDUM = 0.0 CALL STDWRD(IMNO,'START',DDUM,1,1,UNIT,STAT) DDUM = 1.0 CALL STDWRD(IMNO,'STEP',DDUM,1,1,UNIT,STAT) AREA(1:) = 'FIND/PIX outfile ' CALL STDWRC(IMNO,'IDENT',1,AREA,1,72,UNIT,STAT) AREA(1:) = 'data ' CALL STDWRC(IMNO,'CUNIT',1,AREA,1,32,UNIT,STAT) IF (TID.NE.-99) THEN N = INDEX(TABLE,'.') !already file type? IF (N.LT.1) THEN N = INDEX(TABLE,' ') TABLE(N:) = '.bdf' ENDIF CALL STFRNM('middummprov.bdf',TABLE,STAT) ENDIF C ELSE NCOLS = NAXIS + 1 CALL TBTINI(TABLE,0,F_O_MODE,NCOLS+5,200,TID,STAT) DO 800, N=1,NCOLS-1 CALL TBCINI(TID,D_R4_FORMAT,1,'F10.1',TUNIT(N), + TLABL(N),TCOLS(N),STAT) 800 CONTINUE CALL TBCINI(TID,D_R4_FORMAT,1,'G12.6',TUNIT(4), + TLABL(4),TCOLS(NCOLS),STAT) C IF (AREA(1:1).NE.'?') THEN CALL FNDPXS(RETBUF,IMNO,MADRID(PNTR),NAXIS,NPIX, + ACTION(1:3),XRNG,XRFLAG,CXRNG, + TID,TCOLS,ROWMAX,SPIX,EPIX) ELSE CALL FNDPXT(RETBUF,IMNO,MADRID(PNTR),NAXIS,NPIX, + ACTION(1:3),XRNG,XRFLAG,CXRNG, + TID,TCOLS,ROWMAX) ENDIF CALL TBSINI(TID,STAT) CALL DSCUPT(TID,TID,' ',STAT) !create descr. HISTORY CALL TBTCLO(TID,STAT) ENDIF ENDIF ENDIF C RETURN END SUBROUTINE FNDMX(INFO,IMNO,A,NAXIS,NPIX) C IMPLICIT NONE C INTEGER INFO(*),IMNO,NAXIS,NPIX(*) INTEGER K,M,N,UNIT(1),NOLOOP INTEGER KOFF,KLAST,CHUNK,IAV,STAT INTEGER MIMPIX(2),MINPIX,MMIN,LMIN,KMIN INTEGER OUTI(6),MAXPIX,MMAX,LMAX,KMAX C DOUBLE PRECISION FINDXY(4),WINDXY(4) C REAL A(*) REAL CUTS(2),FCUTS(2),OUTR(8) C CHARACTER CBUF*80 C CHUNK = INFO(1) * NPIX(1) !no. of lines * x-dim KLAST = NPIX(2) * NPIX(3) !total no. of lines NOLOOP = KLAST / INFO(1) IF ((NOLOOP*INFO(1)).LT.KLAST) NOLOOP = NOLOOP + 1 KOFF = 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) CALL MNMX(A,IAV,CUTS(1),MIMPIX) !do it once MINPIX = MIMPIX(1) MAXPIX = MIMPIX(2) C IF (NOLOOP.GT.1) THEN DO 100, N=2,NOLOOP KOFF = KOFF + IAV CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) CALL MNMX(A,IAV,FCUTS(1),MIMPIX) IF (FCUTS(1).LT.CUTS(1)) THEN MINPIX = MIMPIX(1) + KOFF - 1 CUTS(1) = FCUTS(1) ENDIF IF (FCUTS(2).GT.CUTS(2)) THEN MAXPIX = MIMPIX(2) + KOFF - 1 CUTS(2) = FCUTS(2) ENDIF 100 CONTINUE ENDIF C C now we have the final data IF (NAXIS.EQ.3) THEN K = NPIX(1)*NPIX(2) M = MINPIX/K MMIN = MINPIX - M*K IF (MMIN.EQ.0) THEN MMIN = M LMIN = NPIX(2) KMIN = NPIX(1) ELSE MINPIX = MMIN MMIN = M + 1 N = MINPIX/NPIX(1) KMIN = MINPIX - N*NPIX(1) IF (KMIN.EQ.0) THEN LMIN = N KMIN = NPIX(1) ELSE LMIN = N + 1 ENDIF ENDIF ELSE IF (NAXIS.EQ.2) THEN N = MINPIX/NPIX(1) KMIN = MINPIX - N*NPIX(1) IF (KMIN.EQ.0) THEN LMIN = N KMIN = NPIX(1) ELSE LMIN = N + 1 ENDIF ELSE KMIN = MINPIX ENDIF C IF (NAXIS.EQ.3) THEN K = NPIX(1)*NPIX(2) M = MAXPIX/K MMAX = MAXPIX - M*K IF (MMAX.EQ.0) THEN MMAX = M LMAX = NPIX(2) KMAX = NPIX(1) ELSE MAXPIX = MMAX MMAX = M + 1 N = MAXPIX/NPIX(1) KMAX = MAXPIX - N*NPIX(1) IF (KMAX.EQ.0) THEN LMAX = N KMAX = NPIX(1) ELSE LMAX = N + 1 ENDIF ENDIF ELSE IF (NAXIS.EQ.2) THEN N = MAXPIX/NPIX(1) KMAX = MAXPIX - N*NPIX(1) IF (KMAX.EQ.0) THEN LMAX = N KMAX = NPIX(1) ELSE LMAX = N + 1 ENDIF ELSE KMAX = MAXPIX ENDIF C C use the correct WCS conversion routines C DO 440,N=1,4 FINDXY(N) = 0.0 WINDXY(N) = 0.0 440 CONTINUE OUTR(1) = CUTS(1) OUTR(2) = CUTS(2) DO 500,N=1,6 OUTI(N) = 1 OUTR(N+2) = 0.0 500 CONTINUE CALL FPXWCO(0,IMNO,FINDXY,WINDXY,STAT) !init conversion C OUTI(1) = KMIN OUTI(3) = KMAX FINDXY(1) = KMIN C IF (NAXIS.EQ.1) THEN CALL FPXWCO(1,IMNO,FINDXY,WINDXY,STAT) !frame pix -> wc OUTR(3) = WINDXY(1) FINDXY(1) = KMAX CALL FPXWCO(1,IMNO,FINDXY,WINDXY,STAT) OUTR(5) = WINDXY(1) WRITE(CBUF,20000) CUTS(1),CUTS(2),KMIN,KMAX C ELSE IF (NAXIS.EQ.2) THEN FINDXY(2) = LMIN OUTI(2) = LMIN CALL FPXWCO(1,IMNO,FINDXY,WINDXY,STAT) !frame pix -> wc OUTR(3) = WINDXY(1) OUTR(4) = WINDXY(2) FINDXY(1) = KMAX FINDXY(2) = LMAX OUTI(4) = LMAX CALL FPXWCO(1,IMNO,FINDXY,WINDXY,STAT) !frame pix -> wc OUTR(5) = WINDXY(1) OUTR(6) = WINDXY(2) WRITE(CBUF,20001) CUTS(1),CUTS(2) CALL STTPUT(CBUF,STAT) WRITE(CBUF,20002) KMIN,LMIN,KMAX,LMAX C ELSE IF (NAXIS.EQ.3) THEN FINDXY(2) = LMIN FINDXY(3) = MMIN OUTI(2) = LMIN OUTI(5) = MMIN CALL FPXWCO(1,IMNO,FINDXY,WINDXY,STAT) !frame pix -> wc OUTR(3) = WINDXY(1) OUTR(4) = WINDXY(2) OUTR(7) = WINDXY(3) FINDXY(1) = KMAX FINDXY(2) = LMAX FINDXY(3) = MMAX OUTI(4) = LMAX OUTI(6) = MMAX CALL FPXWCO(1,IMNO,FINDXY,WINDXY,STAT) !frame pix -> wc OUTR(5) = WINDXY(1) OUTR(6) = WINDXY(2) OUTR(8) = WINDXY(3) WRITE(CBUF,20001) CUTS(1),CUTS(2) CALL STTPUT(CBUF,STAT) WRITE(CBUF,20003) KMIN,LMIN,MMIN,KMAX,LMAX,MMAX ENDIF C CALL STTPUT(CBUF,STAT) CALL STKWRI('OUTPUTI',OUTI,1,6,UNIT,STAT) CALL STKWRR('OUTPUTR',OUTR,1,8,UNIT,STAT) C C finally update descr. LHCUTS CALL STDWRR(IMNO,'LHCUTS',CUTS,3,2,UNIT,STAT) C C that's it folks... RETURN C C Formats 20000 FORMAT('min =',G15.5,', max =',G15.5, + ' at frame pixel (',I6,'), (',I6,')') 20001 FORMAT('min =',G15.5,', max =',G15.5) 20002 FORMAT('at frame pixel (',I6,',',I6,'), (',I6,',',I6,')') 20003 FORMAT('at frame pixel (',I6,',',I6,',',I6,'), (', + I6,',',I6,',',I6,')') END SUBROUTINE FNDPX(INFO,IMNO,A,NAXIS,NPIX,ACTION, + XRNG,XRFLAG,CXB) C IMPLICIT NONE C INTEGER INFO(*),IMNO,NAXIS,NPIX(*),XRFLAG INTEGER N,IX,IY,IZ,KX,KY,KZ,INDXY(3),UNIT(1) INTEGER KOFF,KLAST,CHUNK,IAV,STAT C REAL A(*),XRNG(*) REAL XL,XH,VAL C CHARACTER*(*) ACTION,CXB CHARACTER SIDE*8,CBUF*80 C C init KX = 0 KY = 0 KZ = 0 XL = XRNG(1) XH = XRNG(2) CBUF(1:) = ' ' CHUNK = INFO(1) * NPIX(1) !no. of lines * x-dim KLAST = NPIX(1) * NPIX(2) * NPIX(3) !total no. of pixels KOFF = 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 C C branch according to action IF (ACTION(1:3).EQ.'OUT') GOTO 2000 C SIDE(1:) = 'inside ' C C find 1. pixel inside [XL,XH] IF (XRFLAG.EQ.0) THEN DO 500, IZ=1,NPIX(3) DO 450, IY=1,NPIX(2) DO 400, IX=1,NPIX(1) IF ((A(N).GE.XL).AND.(A(N).LE.XH)) THEN KX = IX KY = IY KZ = IZ VAL = A(N) GOTO 8000 ENDIF N = N + 1 400 CONTINUE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (KOFF.LE.KLAST)) THEN CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 ENDIF 450 CONTINUE 500 CONTINUE C ELSE IF (XRFLAG.EQ.1) THEN ! <,XH DO 700, IZ=1,NPIX(3) DO 650, IY=1,NPIX(2) DO 600, IX=1,NPIX(1) IF (A(N).LE.XH) THEN KX = IX KY = IY KZ = IZ VAL = A(N) GOTO 8000 ENDIF N = N + 1 600 CONTINUE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (KOFF.LE.KLAST)) THEN CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 ENDIF 650 CONTINUE 700 CONTINUE C ELSE IF (XRFLAG.EQ.2) THEN ! XL,> DO 900, IZ=1,NPIX(3) DO 850, IY=1,NPIX(2) DO 800, IX=1,NPIX(1) IF (A(N).GE.XL) THEN KX = IX KY = IY KZ = IZ VAL = A(N) GOTO 8000 ENDIF N = N + 1 800 CONTINUE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (KOFF.LE.KLAST)) THEN CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 ENDIF 850 CONTINUE 900 CONTINUE C ELSE KX = 1 KY = 1 KZ = 1 VAL = A(N) GOTO 8000 ENDIF C C no pixel found N = INDEX(CXB,' ') - 1 IF (N.LT.1) N = LEN(CXB) WRITE(CBUF,10000) SIDE,CXB(1:N) CALL STTPUT(CBUF,STAT) GOTO 8100 C C find 1. pixel outside [XL,XH] C 2000 SIDE(1:) = 'outside ' C IF (XRFLAG.EQ.0) THEN DO 2200, IZ=1,NPIX(3) DO 2150, IY=1,NPIX(2) DO 2100, IX=1,NPIX(1) IF ((A(N).LT.XL).OR.(A(N).GT.XH)) THEN KX = IX KY = IY KZ = IZ VAL = A(N) GOTO 8000 ENDIF N = N + 1 2100 CONTINUE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (KOFF.LE.KLAST)) THEN CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 ENDIF 2150 CONTINUE 2200 CONTINUE C ELSE IF (XRFLAG.EQ.1) THEN ! <,XH DO 2400, IZ=1,NPIX(3) DO 2350, IY=1,NPIX(2) DO 2300, IX=1,NPIX(1) IF (A(N).GT.XH) THEN KX = IX KY = IY KZ = IZ VAL = A(N) GOTO 8000 ENDIF N = N + 1 2300 CONTINUE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (KOFF.LE.KLAST)) THEN CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 ENDIF 2350 CONTINUE 2400 CONTINUE C ELSE IF (XRFLAG.EQ.2) THEN ! XL,> DO 2600, IZ=1,NPIX(3) DO 2550, IY=1,NPIX(2) DO 2500, IX=1,NPIX(1) IF (A(N).LT.XL) THEN KX = IX KY = IY KZ = IZ VAL = A(N) GOTO 8000 ENDIF N = N + 1 2500 CONTINUE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (KOFF.LE.KLAST)) THEN CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 ENDIF 2550 CONTINUE 2600 CONTINUE ENDIF C C no pixel found N = INDEX(CXB,' ') - 1 IF (N.LT.1) N = LEN(CXB) WRITE(CBUF,10000) SIDE,CXB(1:N) CALL STTPUT(CBUF,STAT) GOTO 8100 C C we found a matching pixel 8000 N = INDEX(CXB,' ') - 1 IF (N.LT.1) N = LEN(CXB) C IF (NAXIS.EQ.1) THEN WRITE(CBUF,20000) KX,VAL ELSE IF (NAXIS.EQ.2) THEN WRITE(CBUF,20001) KX,KY,VAL ELSE WRITE(CBUF,20002) KX,KY,KZ,VAL ENDIF CALL STTPUT(CBUF,STAT) WRITE(CBUF,30001) SIDE,CXB(1:N) CALL STTPUT(CBUF,STAT) CALL STKWRR('OUTPUTR',VAL,1,1,UNIT,STAT) C 8100 INDXY(1) = KX INDXY(2) = KY INDXY(3) = KZ CALL STKWRI('OUTPUTI',INDXY,1,3,UNIT,STAT) RETURN C 10000 FORMAT('No pixel ',A,'[',A,'] ...') 20000 FORMAT('frame pixel (',I5,'), value = ',G12.5) 20001 FORMAT('frame pixel (',I5,',',I5,'), value = ',G12.5) 20002 FORMAT('frame pixel (',I5,',',I5,',',I5,'), value = ',G12.5) 30001 FORMAT('is 1. pixel ',A, '[',A,']') C END SUBROUTINE FNDPXT(INFO,IMNO,A,NAXIS,NPIX,ACTION,XRNG,XRFLAG, + CXB,TID,TCOLS,MMAX) C C work on full frame C IMPLICIT NONE C INTEGER INFO(*),NAXIS,NPIX(*),XRFLAG,TID,TCOLS(*),MMAX,RMAX INTEGER N,IX,IY,IZ,NCOLS,KCOUNT,UNIT(1),STAT INTEGER KOFF,KLAST,CHUNK,IAV,IMNO INTEGER TFLAG,NDUM,NDOFF,DUMLIM C REAL A(*),XRNG(*) REAL XL,XH,RVALS(4),DUMBUF(1000) C CHARACTER*(*) ACTION,CXB CHARACTER SIDE*8,CBUF*80 C C init KCOUNT = 0 NCOLS = NAXIS + 1 XL = XRNG(1) XH = XRNG(2) CBUF(1:) = ' ' IF (MMAX.LT.0) THEN RMAX = -MMAX ELSE IF (MMAX.GT.0) THEN RMAX = MMAX ELSE RMAX = 99999999 !99 999 999 ENDIF CHUNK = INFO(1) * NPIX(1) !no. of lines * x-dim N = 1 KLAST = NPIX(1) * NPIX(2) * NPIX(3) !total no. of pixels KOFF = 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV C IF (TID.LT.0) THEN TID = -TID IF (TID.EQ.99) THEN TFLAG = -1 !only display output w. coords ELSE TFLAG = 0 !image output pixels only, no coords NDOFF = 1 NDUM = 1 DUMLIM = 1000 !synchronize with array DUMBUF !!! ENDIF ELSE TFLAG = 1 !table output w. coords ENDIF C C branch according to action IF (ACTION(1:3).EQ.'OUT') GOTO 2000 C SIDE(1:) = 'inside ' C C find all pixels inside [XL,XH] IF (XRFLAG.EQ.0) THEN DO 200, IZ=1,NPIX(3) DO 150, IY=1,NPIX(2) DO 100, IX=1,NPIX(1) IF ((A(N).GE.XL).AND.(A(N).LE.XH)) THEN KCOUNT = KCOUNT + 1 RVALS(NCOLS) = A(N) IF (KCOUNT.GT.RMAX) THEN IF (MMAX.LT.0) THEN KCOUNT = RMAX GOTO 3000 !stop already ENDIF ELSE IF (TFLAG.EQ.1) THEN RVALS(1) = IX IF (NAXIS.GE.2) RVALS(2) = IY IF (NAXIS.GE.3) RVALS(3) = IZ CALL TBRWRR(TID,KCOUNT,NCOLS,TCOLS, + RVALS,STAT) ELSE IF (TFLAG.EQ.0) THEN DUMBUF(NDUM) = A(N) NDUM = NDUM + 1 IF (NDUM.GT.DUMLIM) THEN CALL STFPUT(TID,NDOFF,DUMLIM,DUMBUF,STAT) NDOFF = NDOFF + DUMLIM NDUM = 1 ENDIF ELSE IF (NAXIS.EQ.1) THEN WRITE(CBUF,10050) IX,RVALS(NCOLS) ELSE IF (NAXIS.EQ.2) THEN WRITE(CBUF,10051) IX,IY,RVALS(NCOLS) ELSE WRITE(CBUF,10052) IX,IY,IZ,RVALS(NCOLS) ENDIF CALL STTPUT(CBUF,STAT) ENDIF ENDIF ENDIF N = N + 1 100 CONTINUE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (KOFF.LE.KLAST)) THEN CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 ENDIF 150 CONTINUE 200 CONTINUE C ELSE IF (XRFLAG.EQ.1) THEN ! <,XH DO 400, IZ=1,NPIX(3) DO 350, IY=1,NPIX(2) DO 300, IX=1,NPIX(1) IF (A(N).LE.XH) THEN KCOUNT = KCOUNT + 1 RVALS(NCOLS) = A(N) IF (KCOUNT.GT.RMAX) THEN IF (MMAX.LT.0) THEN KCOUNT = RMAX GOTO 3000 !stop already ENDIF ELSE IF (TFLAG.EQ.1) THEN RVALS(1) = IX IF (NAXIS.GE.2) RVALS(2) = IY IF (NAXIS.GE.3) RVALS(3) = IZ CALL TBRWRR(TID,KCOUNT,NCOLS,TCOLS, + RVALS,STAT) ELSE IF (TFLAG.EQ.0) THEN DUMBUF(NDUM) = A(N) NDUM = NDUM + 1 IF (NDUM.GT.DUMLIM) THEN CALL STFPUT(TID,NDOFF,DUMLIM,DUMBUF,STAT) NDOFF = NDOFF + DUMLIM NDUM = 1 ENDIF ELSE IF (NAXIS.EQ.1) THEN WRITE(CBUF,10050) IX,RVALS(NCOLS) ELSE IF (NAXIS.EQ.2) THEN WRITE(CBUF,10051) IX,IY,RVALS(NCOLS) ELSE WRITE(CBUF,10052) IX,IY,IZ,RVALS(NCOLS) ENDIF CALL STTPUT(CBUF,STAT) ENDIF ENDIF ENDIF N = N + 1 300 CONTINUE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (KOFF.LE.KLAST)) THEN CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 ENDIF 350 CONTINUE 400 CONTINUE C ELSE IF (XRFLAG.EQ.2) THEN ! XL,> DO 600, IZ=1,NPIX(3) DO 550, IY=1,NPIX(2) DO 500, IX=1,NPIX(1) IF (A(N).GE.XL) THEN KCOUNT = KCOUNT + 1 RVALS(NCOLS) = A(N) IF (KCOUNT.GT.RMAX) THEN IF (MMAX.LT.0) THEN KCOUNT = RMAX GOTO 3000 !stop already ENDIF ELSE IF (TFLAG.EQ.1) THEN RVALS(1) = IX IF (NAXIS.GE.2) RVALS(2) = IY IF (NAXIS.GE.3) RVALS(3) = IZ CALL TBRWRR(TID,KCOUNT,NCOLS,TCOLS, + RVALS,STAT) ELSE IF (TFLAG.EQ.0) THEN DUMBUF(NDUM) = A(N) NDUM = NDUM + 1 IF (NDUM.GT.DUMLIM) THEN CALL STFPUT(TID,NDOFF,DUMLIM,DUMBUF,STAT) NDOFF = NDOFF + DUMLIM NDUM = 1 ENDIF ELSE IF (NAXIS.EQ.1) THEN WRITE(CBUF,10050) IX,RVALS(NCOLS) ELSE IF (NAXIS.EQ.2) THEN WRITE(CBUF,10051) IX,IY,RVALS(NCOLS) ELSE WRITE(CBUF,10052) IX,IY,IZ,RVALS(NCOLS) ENDIF CALL STTPUT(CBUF,STAT) ENDIF ENDIF ENDIF N = N + 1 500 CONTINUE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (KOFF.LE.KLAST)) THEN CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 ENDIF 550 CONTINUE 600 CONTINUE ENDIF C GOTO 3000 C C find 1. pixel outside [XL,XH] C 2000 SIDE(1:) = 'outside ' C IF (XRFLAG.EQ.0) THEN DO 2200, IZ=1,NPIX(3) DO 2150, IY=1,NPIX(2) DO 2100, IX=1,NPIX(1) IF ((A(N).LT.XL).OR.(A(N).GT.XH)) THEN KCOUNT = KCOUNT + 1 RVALS(NCOLS) = A(N) IF (KCOUNT.GT.RMAX) THEN IF (MMAX.LT.0) THEN KCOUNT = RMAX GOTO 3000 !stop already ENDIF ELSE IF (TFLAG.EQ.1) THEN RVALS(1) = IX IF (NAXIS.GE.2) RVALS(2) = IY IF (NAXIS.GE.3) RVALS(3) = IZ CALL TBRWRR(TID,KCOUNT,NCOLS,TCOLS, + RVALS,STAT) ELSE IF (TFLAG.EQ.0) THEN DUMBUF(NDUM) = A(N) NDUM = NDUM + 1 IF (NDUM.GT.DUMLIM) THEN CALL STFPUT(TID,NDOFF,DUMLIM,DUMBUF,STAT) NDOFF = NDOFF + DUMLIM NDUM = 1 ENDIF ELSE IF (NAXIS.EQ.1) THEN WRITE(CBUF,10050) IX,RVALS(NCOLS) ELSE IF (NAXIS.EQ.2) THEN WRITE(CBUF,10051) IX,IY,RVALS(NCOLS) ELSE WRITE(CBUF,10052) IX,IY,IZ,RVALS(NCOLS) ENDIF CALL STTPUT(CBUF,STAT) ENDIF ENDIF ENDIF N = N + 1 2100 CONTINUE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (KOFF.LE.KLAST)) THEN CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 ENDIF 2150 CONTINUE 2200 CONTINUE C ELSE IF (XRFLAG.EQ.1) THEN ! <,XH DO 2400, IZ=1,NPIX(3) DO 2350, IY=1,NPIX(2) DO 2300, IX=1,NPIX(1) IF (A(N).GT.XH) THEN KCOUNT = KCOUNT + 1 RVALS(NCOLS) = A(N) IF (KCOUNT.GT.RMAX) THEN IF (MMAX.LT.0) THEN KCOUNT = RMAX GOTO 3000 !stop already ENDIF ELSE IF (TFLAG.EQ.1) THEN RVALS(1) = IX IF (NAXIS.GE.2) RVALS(2) = IY IF (NAXIS.GE.3) RVALS(3) = IZ CALL TBRWRR(TID,KCOUNT,NCOLS,TCOLS, + RVALS,STAT) ELSE IF (TFLAG.EQ.0) THEN DUMBUF(NDUM) = A(N) NDUM = NDUM + 1 IF (NDUM.GT.DUMLIM) THEN CALL STFPUT(TID,NDOFF,DUMLIM,DUMBUF,STAT) NDOFF = NDOFF + DUMLIM NDUM = 1 ENDIF ELSE IF (NAXIS.EQ.1) THEN WRITE(CBUF,10050) IX,RVALS(NCOLS) ELSE IF (NAXIS.EQ.2) THEN WRITE(CBUF,10051) IX,IY,RVALS(NCOLS) ELSE WRITE(CBUF,10052) IX,IY,IZ,RVALS(NCOLS) ENDIF CALL STTPUT(CBUF,STAT) ENDIF ENDIF ENDIF N = N + 1 2300 CONTINUE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (KOFF.LE.KLAST)) THEN CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 ENDIF 2350 CONTINUE 2400 CONTINUE C ELSE IF (XRFLAG.EQ.2) THEN ! XL,> DO 2600, IZ=1,NPIX(3) DO 2550, IY=1,NPIX(2) DO 2500, IX=1,NPIX(1) IF (A(N).LT.XL) THEN KCOUNT = KCOUNT + 1 RVALS(NCOLS) = A(N) IF (KCOUNT.GT.RMAX) THEN IF (MMAX.LT.0) THEN KCOUNT = RMAX GOTO 3000 !stop already ENDIF ELSE IF (TFLAG.EQ.1) THEN RVALS(1) = IX IF (NAXIS.GE.2) RVALS(2) = IY IF (NAXIS.GE.3) RVALS(3) = IZ CALL TBRWRR(TID,KCOUNT,NCOLS,TCOLS, + RVALS,STAT) ELSE IF (TFLAG.EQ.0) THEN DUMBUF(NDUM) = A(N) NDUM = NDUM + 1 IF (NDUM.GT.DUMLIM) THEN CALL STFPUT(TID,NDOFF,DUMLIM,DUMBUF,STAT) NDOFF = NDOFF + DUMLIM NDUM = 1 ENDIF ELSE IF (NAXIS.EQ.1) THEN WRITE(CBUF,10050) IX,RVALS(NCOLS) ELSE IF (NAXIS.EQ.2) THEN WRITE(CBUF,10051) IX,IY,RVALS(NCOLS) ELSE WRITE(CBUF,10052) IX,IY,IZ,RVALS(NCOLS) ENDIF CALL STTPUT(CBUF,STAT) ENDIF ENDIF ENDIF N = N + 1 2500 CONTINUE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (KOFF.LE.KLAST)) THEN CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 ENDIF 2550 CONTINUE 2600 CONTINUE ENDIF C C test, if no pixel found 3000 N = INDEX(CXB,' ') - 1 IF (N.LT.1) N = LEN(CXB) IF (KCOUNT.LE.0) THEN WRITE(CBUF,10000) SIDE,CXB(1:N) ELSE IF (KCOUNT.EQ.1) THEN WRITE(CBUF,20001) SIDE,CXB(1:N) ELSE WRITE(CBUF,20000) KCOUNT,SIDE,CXB(1:N) ENDIF ENDIF CALL STTPUT(CBUF,STAT) CALL STKWRI('OUTPUTI',KCOUNT,1,1,UNIT,STAT) INFO(2) = KCOUNT C C for image output we have to check more IF ((TFLAG.EQ.0) .AND. (NDUM.GT.1)) THEN DUMLIM = NDUM - 1 CALL STFPUT(TID,NDOFF,DUMLIM,DUMBUF,STAT) ENDIF RETURN C 10000 FORMAT('No pixel ',A,'[',A,'] ...') 10050 FORMAT('frame pixel (',I6,'), value = ',G12.5) 10051 FORMAT('frame pixel (',I6,',',I6,'), value = ',G12.5) 10052 FORMAT('frame pixel (',I6,',',I6,',',I6,'), value = ',G12.5) 20000 FORMAT(I6,' pixels ',A,'[',A,']') 20001 FORMAT('1 pixel ',A,'[',A,']') C END SUBROUTINE SUBPOL C C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ C C.IDENTIFICATION: C subroutine SUBPOL version 1.00 881028 C K. Banse ESO - Garching C C.KEYWORDS: C polygon fill C C.PURPOSE: C calculate the average value (and rms) according to different methods C of an area defined by the cursor rectangle or a table or by coordinates C C.ALGORITHM: C This info may also be stored either in a descriptor or a table... C C.INPUT/OUTPUT: C C The following keywords are used: C IN_A/C/1/60 input frame/mask C IN_B/C/1/60 test frame/mask C INPUTC/C/1/60 name of table with subimages of IN_A C OUT_A/C/1/60 output frame/mask C P4/C/1/80 thresholds for polygon detection + C right hand side of replacement C.VERSIONS C 1.00 created from version 1.00 as of 860701 C use FORTRAN 77 + new ST interfaces C C-------------------------------------------------------------------------- C IMPLICIT NONE C INTEGER IAV,N,NROW,TABFLG,STAT INTEGER*8 PNTRA,PNTRB,PNTRC,PNTRT INTEGER NAXIS,NPIX(3) INTEGER IMNOA,IMNOB,IMNOC,IMNOT,TID INTEGER TNAXIS,TNPIX(2) INTEGER SUBPIX(4) INTEGER TABCLN(4) INTEGER UNI(1),NULO,MADRID(1) INTEGER GENNUM,INDEX C CHARACTER*60 INFRAM,OUTFRA,TABLE,TSTFRA CHARACTER CUNIT*64,IDENT*72 CHARACTER OUTLAB(4)*16,CBUF*80 CHARACTER CASE*1,INFRA2*60,P4*80 C REAL THRESH(2),FILL,RPIX(4) REAL CUTS(4) C DOUBLE PRECISION START(3),STEP(3),TSTART(3),TSTEP(3),DD C LOGICAL TABNUL(4),SELFLG C COMMON /VMR/ MADRID C INCLUDE 'MID_INCLUDE:ST_DEF.INC' C DATA OUTLAB /'XSTART ','YSTART ','XEND ','YEND '/ DATA CUNIT /' '/, IDENT /' '/ C INCLUDE 'MID_INCLUDE:ST_DAT.INC' C C get input, test, output frame as well as intable (optional) C CALL STKRDC('IN_A',1,1,60,IAV,INFRAM,UNI,NULO,STAT) CALL STKRDC('IN_B',1,1,60,IAV,TSTFRA,UNI,NULO,STAT) CALL STKRDC('OUT_A',1,1,60,IAV,OUTFRA,UNI,NULO,STAT) C C process threshold values + right hand side CALL STKRDC('P4',1,1,80,IAV,P4,UNI,NULO,STAT) IAV = INDEX(P4,'=') IF (IAV.LE.1) + CALL STETER(7,'invalid replacement string...') CALL GENCNV(P4(1:IAV-1),2,2,UNI,THRESH,DD,STAT) IF (STAT.NE.2) + CALL STETER(7,'invalid replacement string...') INFRA2(1:) = P4(IAV+1:)//' ' C IF (GENNUM(INFRA2).EQ.1) THEN !constant CASE = 'C' CALL GENCNV(INFRA2,2,1,UNI,FILL,DD,STAT) ELSE !second input frame CASE = 'I' CALL CLNFRA(INFRA2,INFRA2,0) ENDIF C C look if we use a table CALL STKRDC('INPUTC',1,1,60,IAV,TABLE,UNI,NULO,STAT) IF (TABLE(1:1).EQ.'+') THEN TABFLG = 0 ELSE CALL CLNTAB(TABLE,TABLE,0) TABFLG = 1 C C test, if labels there CALL TBTOPN(TABLE,F_I_MODE,TID,STAT) CALL TBIGET(TID,N,NROW,N,N,N,STAT) !get total no. of rows DO 600, N=1,4 CALL TBLSER(TID,OUTLAB(N),TABCLN(N),STAT) !without the `:' IF (TABCLN(N).EQ.-1) THEN CBUF(1:) = 'column labelled '//CBUF(1:10)// + 'not found... ' CALL STETER(1,CBUF) ENDIF 600 CONTINUE ENDIF C C map input + check dimensions CALL STIGET(INFRAM,D_R4_FORMAT,F_I_MODE,F_IMA_TYPE, + 3,NAXIS,NPIX,START,STEP, + IDENT,CUNIT,PNTRA,IMNOA,STAT) CALL STDRDR(IMNOA,'LHCUTS',1,4,IAV,CUTS,UNI,NULO,STAT) C IF (NAXIS.NE.2) THEN IF (NAXIS.GE.3) THEN IF (NPIX(3).GT.1) THEN CALL STETER(2,'only 2-dim frames supported...') ELSE NAXIS = 2 ENDIF ELSE CALL STETER(2,'only 2-dim frames supported...') ENDIF ENDIF C C map test frame IF (TSTFRA.NE.INFRAM) THEN CALL STIGET(TSTFRA,D_R4_FORMAT,F_I_MODE,F_IMA_TYPE, + 2,TNAXIS,TNPIX,TSTART,TSTEP, + IDENT,CUNIT,PNTRT,IMNOT,STAT) IF ( TNPIX(1)*TNPIX(2) .NE. NPIX(1)*NPIX(2) ) + CALL STETER(6,'frames must have same dimensions...') ELSE IMNOT = IMNOA PNTRT = PNTRA ENDIF C C map second input frame IF ( (CASE.NE.'C') .AND. (TSTFRA.NE.INFRA2) ) THEN CALL STIGET(INFRA2,D_R4_FORMAT,F_I_MODE,F_IMA_TYPE, + 2,TNAXIS,TNPIX,TSTART,TSTEP, + IDENT,CUNIT,PNTRB,IMNOB,STAT) IF ( TNPIX(1)*TNPIX(2) .NE. NPIX(1)*NPIX(2) ) + CALL STETER(6,'frames must have same dimensions...') ELSE IMNOB = IMNOT PNTRB = PNTRT ENDIF C C o.k. - map output frame + copy input frame IF (OUTFRA.NE.INFRAM) THEN CALL STIPUT(OUTFRA,D_R4_FORMAT,F_O_MODE,F_IMA_TYPE, + NAXIS,NPIX,START,STEP, + IDENT,CUNIT,PNTRC,IMNOC,STAT) CALL STDWRR(IMNOC,'LHCUTS',CUTS,1,4,UNI,STAT) CALL COPYF(MADRID(PNTRA),MADRID(PNTRC),NPIX(1)*NPIX(2)) ELSE IMNOC = IMNOA PNTRC = PNTRA ENDIF C C if complete frame we just have to do it once IF (TABFLG.EQ.0) THEN SUBPIX(1) = 1 SUBPIX(2) = 1 SUBPIX(3) = NPIX(1) SUBPIX(4) = NPIX(2) C C and do it via POLY_FILL CALL POLFIL(CASE,MADRID(PNTRA),MADRID(PNTRT), + MADRID(PNTRB),MADRID(PNTRC), + NPIX,SUBPIX,THRESH,FILL) C ELSE C C loop through table DO 2000, N=1,NROW CALL TBSGET(TID,N,SELFLG,STAT) !only use selected rows IF (.NOT.SELFLG) GOTO 2000 C C get next row of values CALL TBRRDR(TID,N,4,TABCLN,RPIX,TABNUL,STAT) SUBPIX(1) = NINT( (RPIX(1)-START(1))/STEP(1) ) + 1 !xbeg SUBPIX(2) = NINT( (RPIX(2)-START(2))/STEP(2) ) + 1 !ybeg SUBPIX(3) = NINT( (RPIX(3)-START(1))/STEP(1) ) + 1 !xend SUBPIX(4) = NINT( (RPIX(4)-START(2))/STEP(2) ) + 1 !yend C C and do it via POLFIL CALL POLFIL(CASE,MADRID(PNTRA),MADRID(PNTRT), + MADRID(PNTRB),MADRID(PNTRC), + NPIX,SUBPIX,THRESH,FILL) 2000 CONTINUE ENDIF C C That's it folks... IF (TABFLG.EQ.1) CALL TBTCLO(TID,STAT) C RETURN END SUBROUTINE FNDPXA(INFO,IMNO,A,NAXIS,NPIX,ACTION, + XRNG,XRFLAG,CXB,SPIX,EPIX) C IMPLICIT NONE C INTEGER INFO(*),IMNO,NAXIS,NPIX(*),XRFLAG INTEGER SPIX(3),EPIX(3),XOFF,YOFF,ZOFF,XSAV,YSAV INTEGER N,IX,IY,IZ,KX,KY,KZ,INDXY(3),UNIT(1) INTEGER KOFF,KLAST,KPLANE,CHUNK,IAV,STAT C REAL A(*),XRNG(*) REAL XL,XH,VAL C CHARACTER*(*) ACTION,CXB CHARACTER SIDE*8,CBUF*80 C C init KX = 0 KY = 0 KZ = 0 XL = XRNG(1) XH = XRNG(2) CBUF(1:) = ' ' CHUNK = INFO(1) * NPIX(1) !no. of lines * x-dim KPLANE = NPIX(1) * NPIX(2) KLAST = KPLANE * NPIX(3) !total no. of pixels XOFF = SPIX(1) - 1 YOFF = (SPIX(2)-1) * NPIX(1) ZOFF = (SPIX(3)-1) * KPLANE KOFF = ZOFF + YOFF + 1 !start at beginning of line CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) N = 1 + XOFF !first valid pixel in A C C branch according to action IF (ACTION(1:3).EQ.'OUT') GOTO 2000 C SIDE(1:) = 'inside ' C C find 1. pixel inside [XL,XH] IF (XRFLAG.EQ.0) THEN DO 500, IZ=SPIX(3),EPIX(3) YSAV = N DO 450, IY=SPIX(2),EPIX(2) XSAV = N DO 400, IX=SPIX(1),EPIX(1) IF ((A(N).GE.XL).AND.(A(N).LE.XH)) THEN KX = IX KY = IY KZ = IZ VAL = A(N) GOTO 8000 ENDIF N = N + 1 400 CONTINUE N = XSAV + NPIX(1) IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IY.NE.EPIX(2))) THEN !if not last line of plane KOFF = ((IZ-1)*KPLANE) + (IY*NPIX(1)) + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 450 CONTINUE N = YSAV + KPLANE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IZ.NE.EPIX(3))) THEN !if not last plane KOFF = (IZ*KPLANE) + YOFF + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 500 CONTINUE C ELSE IF (XRFLAG.EQ.1) THEN ! <,XH DO 700, IZ=SPIX(3),EPIX(3) YSAV = N DO 650, IY=1,NPIX(2) XSAV = N DO 600, IX=1,NPIX(1) IF (A(N).LE.XH) THEN KX = IX KY = IY KZ = IZ VAL = A(N) GOTO 8000 ENDIF N = N + 1 600 CONTINUE N = XSAV + NPIX(1) IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IY.NE.EPIX(2))) THEN !if not last line of plane KOFF = ((IZ-1)*KPLANE) + (IY*NPIX(1)) + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 650 CONTINUE N = YSAV + KPLANE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IZ.NE.EPIX(3))) THEN !if not last plane KOFF = (IZ*KPLANE) + YOFF + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 700 CONTINUE C ELSE IF (XRFLAG.EQ.2) THEN ! XL,> DO 900, IZ=SPIX(3),EPIX(3) YSAV = N DO 850, IY=1,NPIX(2) XSAV = N DO 800, IX=1,NPIX(1) IF (A(N).GE.XL) THEN KX = IX KY = IY KZ = IZ VAL = A(N) GOTO 8000 ENDIF N = N + 1 800 CONTINUE N = XSAV + NPIX(1) IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IY.NE.EPIX(2))) THEN !if not last line of plane KOFF = ((IZ-1)*KPLANE) + (IY*NPIX(1)) + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 850 CONTINUE N = YSAV + KPLANE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IZ.NE.EPIX(3))) THEN !if not last plane KOFF = (IZ*KPLANE) + YOFF + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 900 CONTINUE C ELSE KX = 1 KY = 1 KZ = 1 VAL = A(N) GOTO 8000 ENDIF C C no pixel found N = INDEX(CXB,' ') - 1 IF (N.LT.1) N = LEN(CXB) WRITE(CBUF,10000) SIDE,CXB(1:N) CALL STTPUT(CBUF,STAT) GOTO 8100 C C find 1. pixel outside [XL,XH] C 2000 SIDE(1:) = 'outside ' C IF (XRFLAG.EQ.0) THEN DO 2200, IZ=SPIX(3),EPIX(3) YSAV = N DO 2150, IY=1,NPIX(2) XSAV = N DO 2100, IX=1,NPIX(1) IF ((A(N).LT.XL).OR.(A(N).GT.XH)) THEN KX = IX KY = IY KZ = IZ VAL = A(N) GOTO 8000 ENDIF N = N + 1 2100 CONTINUE N = XSAV + NPIX(1) IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IY.NE.EPIX(2))) THEN !if not last line of plane KOFF = ((IZ-1)*KPLANE) + (IY*NPIX(1)) + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 2150 CONTINUE N = YSAV + KPLANE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IZ.NE.EPIX(3))) THEN !if not last plane KOFF = (IZ*KPLANE) + YOFF + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 2200 CONTINUE C ELSE IF (XRFLAG.EQ.1) THEN ! <,XH DO 2400, IZ=SPIX(3),EPIX(3) YSAV = N DO 2350, IY=1,NPIX(2) XSAV = N DO 2300, IX=1,NPIX(1) IF (A(N).GT.XH) THEN KX = IX KY = IY KZ = IZ VAL = A(N) GOTO 8000 ENDIF N = N + 1 2300 CONTINUE N = XSAV + NPIX(1) IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IY.NE.EPIX(2))) THEN !if not last line of plane KOFF = ((IZ-1)*KPLANE) + (IY*NPIX(1)) + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 2350 CONTINUE N = YSAV + KPLANE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IZ.NE.EPIX(3))) THEN !if not last plane KOFF = (IZ*KPLANE) + YOFF + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 2400 CONTINUE C ELSE IF (XRFLAG.EQ.2) THEN ! XL,> DO 2600, IZ=SPIX(3),EPIX(3) YSAV = N DO 2550, IY=1,NPIX(2) XSAV = N DO 2500, IX=1,NPIX(1) IF (A(N).LT.XL) THEN KX = IX KY = IY KZ = IZ VAL = A(N) GOTO 8000 ENDIF N = N + 1 2500 CONTINUE N = XSAV + NPIX(1) IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IY.NE.EPIX(2))) THEN !if not last line of plane KOFF = ((IZ-1)*KPLANE) + (IY*NPIX(1)) + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 2550 CONTINUE N = YSAV + KPLANE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IZ.NE.EPIX(3))) THEN !if not last plane KOFF = (IZ*KPLANE) + YOFF + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 2600 CONTINUE ENDIF C C no pixel found N = INDEX(CXB,' ') - 1 IF (N.LT.1) N = LEN(CXB) WRITE(CBUF,10000) SIDE,CXB(1:N) CALL STTPUT(CBUF,STAT) GOTO 8100 C C we found a matching pixel 8000 N = INDEX(CXB,' ') - 1 IF (N.LT.1) N = LEN(CXB) C IF (NAXIS.EQ.1) THEN WRITE(CBUF,20000) KX,VAL ELSE IF (NAXIS.EQ.2) THEN WRITE(CBUF,20001) KX,KY,VAL ELSE WRITE(CBUF,20002) KX,KY,KZ,VAL ENDIF CALL STTPUT(CBUF,STAT) WRITE(CBUF,30001) SIDE,CXB(1:N) CALL STTPUT(CBUF,STAT) CALL STKWRR('OUTPUTR',VAL,1,1,UNIT,STAT) C 8100 INDXY(1) = KX INDXY(2) = KY INDXY(3) = KZ CALL STKWRI('OUTPUTI',INDXY,1,3,UNIT,STAT) RETURN C 10000 FORMAT('No pixel ',A,'[',A,'] ...') 20000 FORMAT('frame pixel (',I5,'), value = ',G12.5) 20001 FORMAT('frame pixel (',I5,',',I5,'), value = ',G12.5) 20002 FORMAT('frame pixel (',I5,',',I5,',',I5,'), value = ',G12.5) 30001 FORMAT('is 1. pixel ',A, '[',A,']') C END SUBROUTINE FNDPXS(INFO,IMNO,A,NAXIS,NPIX,ACTION,XRNG,XRFLAG, + CXB,TID,TCOLS,MMAX,SPIX,EPIX) C C as FNDPXT but work on window of frame C IMPLICIT NONE C INTEGER INFO(*),NAXIS,NPIX(*),XRFLAG,TID,TCOLS(*),MMAX,RMAX INTEGER SPIX(3),EPIX(3),XOFF,YOFF,ZOFF,XSAV,YSAV INTEGER N,IX,IY,IZ,NCOLS,KCOUNT,UNIT(1),STAT INTEGER KOFF,KLAST,KPLANE,CHUNK,IAV,IMNO INTEGER TFLAG,NDUM,NDOFF,DUMLIM C REAL A(*),XRNG(*) REAL XL,XH,RVALS(4),DUMBUF(1000) C CHARACTER*(*) ACTION,CXB CHARACTER SIDE*8,CBUF*80 C C init KCOUNT = 0 NCOLS = NAXIS + 1 XL = XRNG(1) XH = XRNG(2) CBUF(1:) = ' ' IF (MMAX.LT.0) THEN RMAX = -MMAX ELSE IF (MMAX.GT.0) THEN RMAX = MMAX ELSE RMAX = 99999999 !99 999 999 ENDIF CHUNK = INFO(1) * NPIX(1) !no. of lines * x-dim KPLANE = NPIX(1) * NPIX(2) KLAST = KPLANE * NPIX(3) !total no. of pixels XOFF = SPIX(1) - 1 YOFF = (SPIX(2)-1) * NPIX(1) ZOFF = (SPIX(3)-1) * KPLANE KOFF = ZOFF + YOFF + 1 !start at beginning of line CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) N = 1 + XOFF !first valid pixel in A C IF (TID.LT.0) THEN TID = -TID IF (TID.EQ.99) THEN TFLAG = -1 !only display output w. coords ELSE TFLAG = 0 !image output pixels only, no coords NDOFF = 1 NDUM = 1 DUMLIM = 1000 !synchronize with array DUMBUF !!! ENDIF ELSE TFLAG = 1 !table output w. coords ENDIF C C branch according to action IF (ACTION(1:3).EQ.'OUT') GOTO 2000 C SIDE(1:) = 'inside ' C C find all pixels inside [XL,XH] IF (XRFLAG.EQ.0) THEN DO 200, IZ=SPIX(3),EPIX(3) YSAV = N DO 150, IY=SPIX(2),EPIX(2) XSAV = N DO 100, IX=SPIX(1),EPIX(1) IF ((A(N).GE.XL).AND.(A(N).LE.XH)) THEN KCOUNT = KCOUNT + 1 RVALS(NCOLS) = A(N) IF (KCOUNT.GT.RMAX) THEN IF (MMAX.LT.0) THEN KCOUNT = RMAX GOTO 3000 !stop already ENDIF ELSE IF (TFLAG.EQ.1) THEN RVALS(1) = IX IF (NAXIS.GE.2) RVALS(2) = IY IF (NAXIS.GE.3) RVALS(3) = IZ CALL TBRWRR(TID,KCOUNT,NCOLS,TCOLS, + RVALS,STAT) ELSE IF (TFLAG.EQ.0) THEN DUMBUF(NDUM) = A(N) NDUM = NDUM + 1 IF (NDUM.GT.DUMLIM) THEN CALL STFPUT(TID,NDOFF,DUMLIM,DUMBUF,STAT) NDOFF = NDOFF + DUMLIM NDUM = 1 ENDIF ELSE IF (NAXIS.EQ.1) THEN WRITE(CBUF,10050) IX,RVALS(NCOLS) ELSE IF (NAXIS.EQ.2) THEN WRITE(CBUF,10051) IX,IY,RVALS(NCOLS) ELSE WRITE(CBUF,10052) IX,IY,IZ,RVALS(NCOLS) ENDIF CALL STTPUT(CBUF,STAT) ENDIF ENDIF ENDIF N = N + 1 100 CONTINUE N = XSAV + NPIX(1) IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IY.NE.EPIX(2))) THEN !if not last line of plane KOFF = ((IZ-1)*KPLANE) + (IY*NPIX(1)) + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 150 CONTINUE N = YSAV + KPLANE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IZ.NE.EPIX(3))) THEN !if not last plane KOFF = (IZ*KPLANE) + YOFF + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 200 CONTINUE C ELSE IF (XRFLAG.EQ.1) THEN ! <,XH DO 400, IZ=SPIX(3),EPIX(3) YSAV = N DO 350, IY=SPIX(2),EPIX(2) XSAV = N DO 300, IX=SPIX(1),EPIX(1) IF (A(N).LE.XH) THEN KCOUNT = KCOUNT + 1 RVALS(NCOLS) = A(N) IF (KCOUNT.GT.RMAX) THEN IF (MMAX.LT.0) THEN KCOUNT = RMAX GOTO 3000 !stop already ENDIF ELSE IF (TFLAG.EQ.1) THEN RVALS(1) = IX IF (NAXIS.GE.2) RVALS(2) = IY IF (NAXIS.GE.3) RVALS(3) = IZ CALL TBRWRR(TID,KCOUNT,NCOLS,TCOLS, + RVALS,STAT) ELSE IF (TFLAG.EQ.0) THEN DUMBUF(NDUM) = A(N) NDUM = NDUM + 1 IF (NDUM.GT.DUMLIM) THEN CALL STFPUT(TID,NDOFF,DUMLIM,DUMBUF,STAT) NDOFF = NDOFF + DUMLIM NDUM = 1 ENDIF ELSE IF (NAXIS.EQ.1) THEN WRITE(CBUF,10050) IX,RVALS(NCOLS) ELSE IF (NAXIS.EQ.2) THEN WRITE(CBUF,10051) IX,IY,RVALS(NCOLS) ELSE WRITE(CBUF,10052) IX,IY,IZ,RVALS(NCOLS) ENDIF CALL STTPUT(CBUF,STAT) ENDIF ENDIF ENDIF N = N + 1 300 CONTINUE N = XSAV + NPIX(1) IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IY.NE.EPIX(2))) THEN !if not last line of plane KOFF = ((IZ-1)*KPLANE) + (IY*NPIX(1)) + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 350 CONTINUE N = YSAV + KPLANE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IZ.NE.EPIX(3))) THEN !if not last plane KOFF = (IZ*KPLANE) + YOFF + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 400 CONTINUE C ELSE IF (XRFLAG.EQ.2) THEN ! XL,> DO 600, IZ=SPIX(3),EPIX(3) YSAV = N DO 550, IY=SPIX(2),EPIX(2) XSAV = N DO 500, IX=SPIX(1),EPIX(1) IF (A(N).GE.XL) THEN KCOUNT = KCOUNT + 1 RVALS(NCOLS) = A(N) IF (KCOUNT.GT.RMAX) THEN IF (MMAX.LT.0) THEN KCOUNT = RMAX GOTO 3000 !stop already ENDIF ELSE IF (TFLAG.EQ.1) THEN RVALS(1) = IX IF (NAXIS.GE.2) RVALS(2) = IY IF (NAXIS.GE.3) RVALS(3) = IZ CALL TBRWRR(TID,KCOUNT,NCOLS,TCOLS, + RVALS,STAT) ELSE IF (TFLAG.EQ.0) THEN DUMBUF(NDUM) = A(N) NDUM = NDUM + 1 IF (NDUM.GT.DUMLIM) THEN CALL STFPUT(TID,NDOFF,DUMLIM,DUMBUF,STAT) NDOFF = NDOFF + DUMLIM NDUM = 1 ENDIF ELSE IF (NAXIS.EQ.1) THEN WRITE(CBUF,10050) IX,RVALS(NCOLS) ELSE IF (NAXIS.EQ.2) THEN WRITE(CBUF,10051) IX,IY,RVALS(NCOLS) ELSE WRITE(CBUF,10052) IX,IY,IZ,RVALS(NCOLS) ENDIF CALL STTPUT(CBUF,STAT) ENDIF ENDIF ENDIF N = N + 1 500 CONTINUE N = XSAV + NPIX(1) IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IY.NE.EPIX(2))) THEN !if not last line of plane KOFF = ((IZ-1)*KPLANE) + (IY*NPIX(1)) + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 550 CONTINUE N = YSAV + KPLANE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IZ.NE.EPIX(3))) THEN !if not last plane KOFF = (IZ*KPLANE) + YOFF + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 600 CONTINUE ENDIF C GOTO 3000 C C find 1. pixel outside [XL,XH] C 2000 SIDE(1:) = 'outside ' C IF (XRFLAG.EQ.0) THEN DO 2200, IZ=SPIX(3),EPIX(3) YSAV = N DO 2150, IY=SPIX(2),EPIX(2) XSAV = N DO 2100, IX=SPIX(1),EPIX(1) IF ((A(N).LT.XL).OR.(A(N).GT.XH)) THEN KCOUNT = KCOUNT + 1 RVALS(NCOLS) = A(N) IF (KCOUNT.GT.RMAX) THEN IF (MMAX.LT.0) THEN KCOUNT = RMAX GOTO 3000 !stop already ENDIF ELSE IF (TFLAG.EQ.1) THEN RVALS(1) = IX IF (NAXIS.GE.2) RVALS(2) = IY IF (NAXIS.GE.3) RVALS(3) = IZ CALL TBRWRR(TID,KCOUNT,NCOLS,TCOLS, + RVALS,STAT) ELSE IF (TFLAG.EQ.0) THEN DUMBUF(NDUM) = A(N) NDUM = NDUM + 1 IF (NDUM.GT.DUMLIM) THEN CALL STFPUT(TID,NDOFF,DUMLIM,DUMBUF,STAT) NDOFF = NDOFF + DUMLIM NDUM = 1 ENDIF ELSE IF (NAXIS.EQ.1) THEN WRITE(CBUF,10050) IX,RVALS(NCOLS) ELSE IF (NAXIS.EQ.2) THEN WRITE(CBUF,10051) IX,IY,RVALS(NCOLS) ELSE WRITE(CBUF,10052) IX,IY,IZ,RVALS(NCOLS) ENDIF CALL STTPUT(CBUF,STAT) ENDIF ENDIF ENDIF N = N + 1 2100 CONTINUE N = XSAV + NPIX(1) IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IY.NE.EPIX(2))) THEN !if not last line of plane KOFF = ((IZ-1)*KPLANE) + (IY*NPIX(1)) + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 2150 CONTINUE N = YSAV + KPLANE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IZ.NE.EPIX(3))) THEN !if not last plane KOFF = (IZ*KPLANE) + YOFF + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 2200 CONTINUE C ELSE IF (XRFLAG.EQ.1) THEN ! <,XH DO 2400, IZ=SPIX(3),EPIX(3) YSAV = N DO 2350, IY=SPIX(2),EPIX(2) XSAV = N DO 2300, IX=SPIX(1),EPIX(1) IF (A(N).GT.XH) THEN KCOUNT = KCOUNT + 1 RVALS(NCOLS) = A(N) IF (KCOUNT.GT.RMAX) THEN IF (MMAX.LT.0) THEN KCOUNT = RMAX GOTO 3000 !stop already ENDIF ELSE IF (TFLAG.EQ.1) THEN RVALS(1) = IX IF (NAXIS.GE.2) RVALS(2) = IY IF (NAXIS.GE.3) RVALS(3) = IZ CALL TBRWRR(TID,KCOUNT,NCOLS,TCOLS, + RVALS,STAT) ELSE IF (TFLAG.EQ.0) THEN DUMBUF(NDUM) = A(N) NDUM = NDUM + 1 IF (NDUM.GT.DUMLIM) THEN CALL STFPUT(TID,NDOFF,DUMLIM,DUMBUF,STAT) NDOFF = NDOFF + DUMLIM NDUM = 1 ENDIF ELSE IF (NAXIS.EQ.1) THEN WRITE(CBUF,10050) IX,RVALS(NCOLS) ELSE IF (NAXIS.EQ.2) THEN WRITE(CBUF,10051) IX,IY,RVALS(NCOLS) ELSE WRITE(CBUF,10052) IX,IY,IZ,RVALS(NCOLS) ENDIF CALL STTPUT(CBUF,STAT) ENDIF ENDIF ENDIF N = N + 1 2300 CONTINUE N = XSAV + NPIX(1) IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IY.NE.EPIX(2))) THEN !if not last line of plane KOFF = ((IZ-1)*KPLANE) + (IY*NPIX(1)) + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 2350 CONTINUE N = YSAV + KPLANE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IZ.NE.EPIX(3))) THEN !if not last plane KOFF = (IZ*KPLANE) + YOFF + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 2400 CONTINUE C ELSE IF (XRFLAG.EQ.2) THEN ! XL,> DO 2600, IZ=SPIX(3),EPIX(3) YSAV = N DO 2550, IY=SPIX(2),EPIX(2) XSAV = N DO 2500, IX=SPIX(1),EPIX(1) IF (A(N).LT.XL) THEN KCOUNT = KCOUNT + 1 RVALS(NCOLS) = A(N) IF (KCOUNT.GT.RMAX) THEN IF (MMAX.LT.0) THEN KCOUNT = RMAX GOTO 3000 !stop already ENDIF ELSE IF (TFLAG.EQ.1) THEN RVALS(1) = IX IF (NAXIS.GE.2) RVALS(2) = IY IF (NAXIS.GE.3) RVALS(3) = IZ CALL TBRWRR(TID,KCOUNT,NCOLS,TCOLS, + RVALS,STAT) ELSE IF (TFLAG.EQ.0) THEN DUMBUF(NDUM) = A(N) NDUM = NDUM + 1 IF (NDUM.GT.DUMLIM) THEN CALL STFPUT(TID,NDOFF,DUMLIM,DUMBUF,STAT) NDOFF = NDOFF + DUMLIM NDUM = 1 ENDIF ELSE IF (NAXIS.EQ.1) THEN WRITE(CBUF,10050) IX,RVALS(NCOLS) ELSE IF (NAXIS.EQ.2) THEN WRITE(CBUF,10051) IX,IY,RVALS(NCOLS) ELSE WRITE(CBUF,10052) IX,IY,IZ,RVALS(NCOLS) ENDIF CALL STTPUT(CBUF,STAT) ENDIF ENDIF ENDIF N = N + 1 2500 CONTINUE N = XSAV + NPIX(1) IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IY.NE.EPIX(2))) THEN !if not last line of plane KOFF = ((IZ-1)*KPLANE) + (IY*NPIX(1)) + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 2550 CONTINUE N = YSAV + KPLANE IF ((N .GT. CHUNK) .AND. !get next chunk of data + (IZ.NE.EPIX(3))) THEN !if not last plane KOFF = (IZ*KPLANE) + YOFF + 1 CALL STFGET(IMNO,KOFF,CHUNK,IAV,A,STAT) KOFF = KOFF + IAV N = 1 + XOFF ENDIF 2600 CONTINUE ENDIF C C test, if no pixel found 3000 N = INDEX(CXB,' ') - 1 IF (N.LT.1) N = LEN(CXB) IF (KCOUNT.LE.0) THEN WRITE(CBUF,10000) SIDE,CXB(1:N) ELSE IF (KCOUNT.EQ.1) THEN WRITE(CBUF,20001) SIDE,CXB(1:N) ELSE WRITE(CBUF,20000) KCOUNT,SIDE,CXB(1:N) ENDIF ENDIF CALL STTPUT(CBUF,STAT) CALL STKWRI('OUTPUTI',KCOUNT,1,1,UNIT,STAT) INFO(2) = KCOUNT C C for image output we have to check more IF ((TFLAG.EQ.0) .AND. (NDUM.GT.1)) THEN DUMLIM = NDUM - 1 CALL STFPUT(TID,NDOFF,DUMLIM,DUMBUF,STAT) ENDIF RETURN C 10000 FORMAT('No pixel ',A,'[',A,'] ...') 10050 FORMAT('frame pixel (',I6,'), value = ',G12.5) 10051 FORMAT('frame pixel (',I6,',',I6,'), value = ',G12.5) 10052 FORMAT('frame pixel (',I6,',',I6,',',I6,'), value = ',G12.5) 20000 FORMAT(I6,' pixels ',A,'[',A,']') 20001 FORMAT('1 pixel ',A,'[',A,']') C END