C @(#)polfil.for	17.1.1.1 (ES0-DMD) 01/25/02 17:40:35
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 POLFIL(CASE,INIM,TSTIM,IN2IM,
     +			OUTIM,NPIX,WINDOW,THRESH,FILL)
C
C ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C	
C.IDENTIFICATION
C  subroutine POLFIL			version 1.00	880912
C  subroutine ISCAN			version 1.00	880912
C  subroutine OVRLAP			version 1.30	881028
C  subroutine PROGRS			version 2.20	880919
C  subroutine MEMDSK			version 1.00	891129
C  K. Banse				ESO - Garching
C 
C.KEYWORDS
C  polygon fill
C 
C.PURPOSE
C  fill pixels inside given polygon (including borderline pixels)
C 
C.ALGORITHM
C  get beginning + end of polygon per line in UP, DOWN
C  replace all pixels inside with constant or other pixel values
C	
C.IN/OUTPUT:
C call as   POLFIL(CASE,INIM,TSTIM,IN2IM,OUTIM,NPIX,WINDOW,THRESH,FILL)
C
C input par:
C CASE:		char. exp.	= C for constant, I for image
C INIM:		R*4		input image
C TSTIM:	R*4		test image
C IN2IM:	R*4		second input image
C OUTIM:	R*4		output image
C NPIX:         I*4 array       NPIX of images
C WINDOW:       I*4 array       lower left + upper right corner of window
C                               in images
C THRESH:       R*4 array       low,high threshold for polygon detection
C FILL:         R*4             fill value
C
C.VERSIONS
C  1.00         from version 1.00 as of 860707
C
C-----------------------------------------------------------------------
C
      IMPLICIT NONE
C
      CHARACTER*1    CASE
C
      INTEGER     WINDOW(4),NPIX(2)
      INTEGER     NOLINS,NOPIX,OFF,UP,DOWN,SW,KX,JX
      INTEGER     M,NN,N,KLO,KHI,MLO,MHI,JLO,JHI
      INTEGER     LOWUP
C      
      REAL        INIM(*),TSTIM(*),OUTIM(*),IN2IM(*)
      REAL        THRESH(2),FILL
C      
C  initalize pointer to lower left corner of window
      OFF = ( (WINDOW(2) - 1) * NPIX(1) ) + WINDOW(1) - 1
      NOLINS = WINDOW(4) - WINDOW(2) + 1
      NOPIX = WINDOW(3) - WINDOW(1) + 1
C
C  here we go
C      
      SW = 0
      IF (CASE.EQ.'C') THEN
         DO 1000 N=1,NOLINS
            IF (SW.NE.0) WRITE(*,12000) (N-1)
            SW = 0
            UP = 0
            DOWN = 0
            KX = 0
            JX = 0
            JLO = 0
            JHI = 0
            DO 400, NN=OFF+1,OFF+NOPIX
               M = NN - 1
               IF ( (TSTIM(NN).GE.THRESH(1)) .AND. 
     +              (TSTIM(NN).LE.THRESH(2)) ) THEN
                  IF (UP.NE.M) THEN
                     IF (DOWN.EQ.0) THEN
                        SW = 1 
                        KLO = LOWUP(TSTIM,THRESH,NPIX,1,NN)
                        KHI = LOWUP(TSTIM,THRESH,NPIX,2,NN)
                        KX = NN
                     ELSE
                        JX = NN
                        JLO = LOWUP(TSTIM,THRESH,NPIX,1,NN)
                        JHI = LOWUP(TSTIM,THRESH,NPIX,2,NN)
                     ENDIF
                  ENDIF
                  UP = NN
               ELSE
                  IF (UP.EQ.M) THEN
                     MLO = LOWUP(TSTIM,THRESH,NPIX,1,M)    !check line below
                     MHI = LOWUP(TSTIM,THRESH,NPIX,2,M)    !check line above
                     IF (DOWN.EQ.0) THEN
                        DOWN = 1 - DOWN
                        IF (KX.NE.M) THEN
                           IF ((MLO.EQ.1).AND.(MHI.EQ.1)) THEN
                              SW = 0
                              DOWN = 0
                           ELSE IF ((KLO.EQ.MLO) .AND. 
     +                              (KHI.EQ.MHI)) THEN
                              SW = 1 - SW 
                              DOWN = 0
                           ENDIF
                        ELSE IF (MHI.NE.MLO) THEN
                           IF ((KLO.EQ.MLO) .AND. (KHI.EQ.MHI)) THEN
                              SW = 1 - SW 
                              DOWN = 0
                           ENDIF
                        ENDIF
                     ELSE
                        IF ((JX.NE.M) .OR.(MHI.NE.MLO)) THEN
                           IF ((JLO.NE.MLO) .OR. (JHI.NE.MHI)) THEN
                              SW = 1 - SW
                              DOWN = 0
                           ENDIF
                        ELSE
                           SW = 1 - SW
                           DOWN = 0
                        ENDIF
                     ENDIF
                  ENDIF
               ENDIF
C 
               IF (SW.NE.0) OUTIM(NN) = FILL
400         CONTINUE
            OFF = OFF + NPIX(1)
1000     CONTINUE
C      
C  use second input image
      ELSE
         DO 2000 N=1,NOLINS
            IF (SW.NE.0) WRITE(*,12000) (N-1)
            SW = 0
            UP = 0
            DOWN = 0
            KX = 0
            JX = 0
            JLO = 0
            JHI = 0
            DO 1400, NN=OFF+1,OFF+NOPIX
               M = NN - 1
               IF ( (TSTIM(NN).GE.THRESH(1)) .AND.
     +              (TSTIM(NN).LE.THRESH(2)) ) THEN
                  IF (UP.NE.M) THEN
                     IF (DOWN.EQ.0) THEN
                        SW = 1
                        KLO = LOWUP(TSTIM,THRESH,NPIX,1,NN)
                        KHI = LOWUP(TSTIM,THRESH,NPIX,2,NN)
                        KX = NN
                     ELSE
                        JX = NN
                        JLO = LOWUP(TSTIM,THRESH,NPIX,1,NN)
                        JHI = LOWUP(TSTIM,THRESH,NPIX,2,NN)
                     ENDIF
                  ENDIF
                  UP = NN
               ELSE
                  IF (UP.EQ.M) THEN
                     MLO = LOWUP(TSTIM,THRESH,NPIX,1,M)    !check line below
                     MHI = LOWUP(TSTIM,THRESH,NPIX,2,M)    !check line above
                     IF (DOWN.EQ.0) THEN
                        DOWN = 1 - DOWN
                        IF (KX.NE.M) THEN
                           IF ((MLO.EQ.1).AND.(MHI.EQ.1)) THEN
                              SW = 0
                              DOWN = 0
                           ELSE IF ((KLO.EQ.MLO) .AND. 
     +                              (KHI.EQ.MHI)) THEN
                              SW = 1 - SW
                              DOWN = 0
                           ENDIF
                        ELSE IF (MHI.NE.MLO) THEN
                           IF ((KLO.EQ.MLO) .AND. (KHI.EQ.MHI)) THEN
                              SW = 1 - SW
                              DOWN = 0
                           ENDIF
                        ENDIF
                     ELSE
                        IF ((JX.NE.M) .OR.(MHI.NE.MLO)) THEN
                           IF ((JLO.NE.MLO) .OR. (JHI.NE.MHI)) THEN
                              SW = 1 - SW
                              DOWN = 0
                           ENDIF
                        ELSE
                           SW = 1 - SW
                           DOWN = 0
                        ENDIF
                     ENDIF
                  ENDIF
               ENDIF
C
               IF (SW.NE.0) OUTIM(NN) = IN2IM(NN)
1400        CONTINUE
            OFF = OFF + NPIX(1)
2000     CONTINUE
C
      ENDIF
C 
C  that's it folks
      RETURN
C 
12000 FORMAT(' line',I5,' has problems...')
C 
      END

      INTEGER FUNCTION LOWUP(A,THR,NPIX,LFLAG,PIXNO)
C      
      IMPLICIT NONE
C      
      REAL        A(*),THR(2)
C      
      INTEGER     PIXNO,NPIX(2),LFLAG
      INTEGER     N,K
C      
      LOWUP = 0
      N = PIXNO - 1
      IF (LFLAG.EQ.1) THEN                  !check lower line
         N = N - NPIX(1)
         IF (N.LT.1) RETURN
C 
         DO 200, K=1,3
            IF ((A(N).GE.THR(1)).AND.(A(N).LE.THR(2))) GOTO 1000
            N = N + 1
200      CONTINUE
      ELSE
         N = N + NPIX(1)
         IF (N.GT.(NPIX(1)*NPIX(2))) RETURN
C 
         DO 400, K=1,3
            IF ((A(N).GE.THR(1)).AND.(A(N).LE.THR(2))) GOTO 1000
            N = N + 1
400      CONTINUE
      ENDIF
      RETURN
C      
1000  LOWUP = 1
      RETURN
      END

      SUBROUTINE OVRLAP(STARTA,STEPA,NPIXA,STARTB,STEPB,NPIXB,
     +                  BEGIN,END,STAT)
C
C++++++++++++++++++++++++++++++++++++++++++++++++++
C
C.IDENTIFICATION
C  subroutine OVRLAP            version 1.30      881028
C  K. Banse                     ESO - Garching
C
C.KEYWORDS
C  conjunction
C
C.PURPOSE
C  determine overlapping region of two "lines" in real space
C
C.ALGORITHM
C  straight forward
C
C.INPUT/OUTPUT
C  call as   OVRLAP(STARTA,STEPA,NPIXA,STARTB,STEPB,NPIXB,BEGIN,END,STAT)
C
C  input par:
C  STARTA:	R*8	start of "line" A
C  STEPA:	R*8	stepsize of A
C  NPIXB:	I*4	no. of steps in A
C  STARTB:	R*8	start of "line" B
C  STEPB:	R*8	stepsize of B
C  NPIXB:	I*4	no. of steps in B
C
C  output par:
C  BEGIN:	R*8	begin of overlapping interval
C  END:		R*8	end of overlapping interval
C  STAT:	I*4	return status, should be = 0
C
C--------------------------------------------------
C
      IMPLICIT NONE
C
      INTEGER      NPIXA,NPIXB,STAT
C
      DOUBLE PRECISION STARTA,STEPA,STARTB,STEPB
      DOUBLE PRECISION EB,BEGIN,END
C
C  compute position of last pixel
      STAT = 0
      END = STARTA + (NPIXA-1)*STEPA
      EB = STARTB + (NPIXB-1)*STEPB
C
C  negative stepsize
      IF (STEPA.LT.0.D0) THEN
         IF (STARTA.LT.STARTB) THEN            !low begin
            BEGIN = STARTA
         ELSE
            BEGIN = STARTB
         ENDIF
         IF (END.LT.EB) END = EB            !high end
C
         IF (END.GT.BEGIN) STAT = 1            !no overlap...!
C
C  positive stepsize
      ELSE
         IF (STARTA.LT.STARTB) THEN            !high begin
            BEGIN = STARTB
         ELSE
            BEGIN = STARTA
         ENDIF
         IF (END.GT.EB) END = EB            !low end
C
         IF (END.LT.BEGIN) STAT = 1            !no overlap...!
      ENDIF
C
      RETURN
      END

      SUBROUTINE PROGRS(INC,CINC,PRCNT,CHECK)
C
C++++++++++++++++++++++++++++++++++++++++++++++++++
C
C.IDENTIFICATION
C  subroutine PROGRS		version 2.00	860609
C  K. Banse			ESO - Garching
C  					2.20	880919
C
C.KEYWORDS
C  time, progress of work
C
C.PURPOSE
C  display current time + percentage of work already done
C
C.ALGORITHM
C  straight forward 
C
C.INPUT/OUTPUT
C  call as PROGRS(INC,CINC,PRCNT,CHECK)
C
C  input par:
C  INC:		I*4		no. of percent to increase
C  CINC:	I*4		increment for check
C  
C  in/output par:
C  PRCNT:	I*4		percent of processing done
C  CHECK:	I*4		check value
C 
C.VERSIONS
C  2.10		do not write to logfile - only to terminal...
C  2.20		move to FORTRAN 77 + use system independent calls
C 
C--------------------------------------------------
C
      IMPLICIT NONE
C      
      INTEGER     INC,CINC,PRCNT,CHECK
      INTEGER     N,IAV,KLOG,NLOG,STAT
      INTEGER     UNIT(1),NULLO
C      
      CHARACTER   TIME*40,PBUF*60
C      
      IF (PRCNT.GT.100) RETURN
C
C  get current time, convert percentage to ASCII + display all that
      TIME = ' '
      CALL GENTIM(TIME)
      DO 100 N=40,2,-1
         IF (TIME(N:N).NE.' ') THEN
            IAV = N
            GOTO 200
         ENDIF
100   CONTINUE
      IAV = 1
C 
200   PBUF(1:) = ' '
      WRITE(PBUF(1:),10000) TIME(1:IAV),PRCNT
C      
C  set log flag to 0, so we do not write into the log file...
      CALL STKRDI('LOG',1,1,IAV,KLOG,UNIT,NULLO,STAT)
      NLOG = 0
      CALL STKWRI('LOG',NLOG,1,1,UNIT,STAT)
      CALL STTPUT(PBUF,STAT)
      CALL STKWRI('LOG',KLOG,1,1,UNIT,STAT)
C      
C  update
      PRCNT = PRCNT + INC
      CHECK = CHECK + CINC
      RETURN
C      
10000 FORMAT(A,I4,'% done ... ')
C      
      END

      SUBROUTINE MEMDSK(IMNO,FLAG,A,NPIX,STA,WSIZE,STAT)
C
C++++++++++++++++++++++++++++++++++++++++++++++++++
C
C.IDENTIFICATION
C  subroutine MEMDSK            version 1.00      891129
C  K. Banse                     ESO - Garching
C
C.KEYWORDS
C  memory, diskfile
C
C.PURPOSE
C  copy data stored in virtual memory to right place in disk file
C
C.ALGORITHM
C  use STFPUT interface
C
C.INPUT/OUTPUT
C  call as   MEMDSK(IMNO,FLAG,A,NPIX,STA,WSIZE,STAT)
C
C  input par:
C  IMNO:	Integer		file no.
C  FLAG:	Integer		= 0, buffer A holds full frame
C				= 1, buffer A holds data from window start on
C				= 2, buffer A holds only the data of the window
C  A:		Real		buffer with data
C  NPIX:	Integer		no. of pixels in A
C  STA:		Integer		start of window in total frame
C  WSIZE:	Integer		no. of pixels in window
C
C  output par:
C  STAT:	Integer		return status, should be = 0
C
C--------------------------------------------------
C
      IMPLICIT NONE
C 
      INTEGER     IMNO,FLAG,NPIX(*),STA(*),WSIZE(*),STAT
      INTEGER     FELEM,SIZE,N,IOFF
C 
      REAL        A(*)
C 
      FELEM = ((STA(2)-1) * NPIX(1)) + STA(1)
      IF (FLAG .EQ. 0) THEN
         SIZE = WSIZE(2) * NPIX(1)              !WSIZE(2) lines...
         CALL STFPUT(IMNO,FELEM,SIZE,A(FELEM),STAT)
      ELSE IF (FLAG .EQ. 1) THEN
         SIZE = WSIZE(2) * NPIX(1)              !WSIZE(2) lines...
         CALL STFPUT(IMNO,FELEM,SIZE,A(1),STAT)
      ELSE IF (FLAG .EQ. 2) THEN
         IOFF = 1
         DO 1000 N=1,WSIZE(2)                   !again WSIZE(2) lines ...
            CALL STFPUT(IMNO,FELEM,WSIZE(1),A(IOFF),STAT)
            IOFF = IOFF + WSIZE(1)
            FELEM = FELEM + NPIX(1)
1000     CONTINUE
      ENDIF
C 
      RETURN
      END