SUBROUTINE PARAPZ
C
C PARAPREDRIZ V0.62 - prepare a script for drizzling an input EIS
C image onto all affected output superimage tiles.
C This version is designed for running in parallel.
C
C It is an IRAF application written using the F77/VOS.
C
C Originally based on a test application called taste.
C
C Richard Hook, July 1997
C
C History:
C
C Modified for single image mode, Richard Hook, August 1997
C Modified for output as .hhh format, September 1997
C Modified to include context map, November 1997
C Modified to include "flprc" before the drizzle, November 1997
C Modified to also set "contab" for the global context table, 24/11/97
C Modified to re-order and re-name the output, 26/11/97
C Improved reading of double precision quantities and version, 17/12/97
C Changed fillval to indef, 27/1/98
C Increased buffer size for DATA, 8/6/98
C Rebuilt for new astrometric headers and added PROJ=COE or TAN, 11/2/99
C PARAPREDRIZ created for parallel output streams to different files, 15/2/99
C Modified for new convention with context information in header, 24/3/99
C Modified to make context images optional and output some useful
C numbers, 30/3/99
C Minor modifications for WFI pipeline data testing, 30.6.99
C
C Version 0.4, 3rd August 1999 - new "bestfit" option and slightly
C                              different use of STEM.
C
C Version 0.5, 20th April 2000 - correct minor problem with "single"
C                                mode to avoid .hhh when one asks for .fits
C
C Version 0.6, 10th August 2000 - put expkey="exptime" so that it works
C                                 with drizzle 0.85 and uses the exposure
C                                 time information.
C
C Version 0.61, 11th September 2001 - cleaner use of time functions to
C                                     get things working on HPUX 11.
C
C Version 0.62, 3rd May 2002 - correct bug when closing files, pointed
C                              out by Thomas Erben
C
      IMPLICIT NONE
      INTEGER NC,ISTAT,ID,MPIX,DIMS(7)
      INTEGER DATTYP,NDIMS,NIM
      INTEGER NXST,NXEN,NYST,NYEN,J,NX,NY,K,N,XZ,YZ
      INTEGER I,IA,JA,IMLD,I1,I2,LD,IP,IS,USEOF,OI1,OI2
      INTEGER II,JJ,KK,IOLD,JOLD,IITIME,NPARA,NP,NSEC
      CHARACTER*200 CHARS
      CHARACTER*48 VERS
      CHARACTER*80 IMAGES,MASK,DIR,OUTDAT
      CHARACTER*11 PNAME
      CHARACTER*3 PROJ
      CHARACTER*5 EXT
      INTEGER MAXIM
      PARAMETER (MAXIM=1000)
      INTEGER XMIN(MAXIM),XMAX(MAXIM),YMIN(MAXIM),YMAX(MAXIM)
      INTEGER XI,XA,YI,YA,NXP,NYP,XPIXOF,YPIXOF
      CHARACTER*80 IMAGE(MAXIM)
      CHARACTER*80 DATA(MAXIM)
      INTEGER XS(MAXIM),YS(MAXIM)
      INTEGER IMAX,IMIN,JMAX,JMIN
      CHARACTER*4 STEM
      CHARACTER*8 MODE
      CHARACTER*1 XSIGN,YSIGN
      CHARACTER*24 ICTIME
      EXTERNAL IITIME,ICTIME
      DOUBLE PRECISION CRVAL1,CRPIX1,CRVAL2,CRPIX2
      DOUBLE PRECISION CDELT1,CDELT2,XCO(100),YCO(100)
      DOUBLE PRECISION RAOFF,DECOFF,FLXSCL,PIBY,SCALE
      DOUBLE PRECISION RAREF,DECREF,RPP,XPOFF,YPOFF
      DOUBLE PRECISION OUTRA,OUTDEC,XOUT,YOUT
      DOUBLE PRECISION OXMIN,OXMAX,OYMAX,OYMIN
      LOGICAL TOP
      LOGICAL CON
      LOGICAL OVER
      LOGICAL SINGLE

      PARAMETER(PIBY=3.141592653589/180.0)
      PARAMETER (USEOF=-2)

C Initializations
      XPOFF=0.0
      YPOFF=0.0
      XPIXOF=0
      YPIXOF=0

      TOP=.FALSE.
      DO I=1,MAXIM
         IMAGE(I)=' '
      ENDDO

      IMAGES=' '
      IMAX=-999
      JMAX=-999
      IMIN=999
      JMIN=999
      NIM=1

C Get mode
      CALL UCLGST('outmode',MODE,ISTAT)

C Set version 
      VERS='PARAPREDRIZ Version 0.62 (3rd May 2002)'
      IF(MODE.NE.'single') THEN
         CALL UMSPUT('+ '//VERS(1:48),1,0,ISTAT)
         CALL UMSPUT(' ',1,0,ISTAT)
      ENDIF

C Get the output projection (this is just passed on to Drizzle)
      CALL UCLGST('proj',PROJ,ISTAT)

C Check whether we want context images or not
      CALL UCLGSB('context',CON,ISTAT)

C Get the stem for the output file names
      CALL UCLGST('stem',STEM,ISTAT)

C Pad the name
      DO I=1,4
         IF(STEM(I:I).EQ.' ') STEM(I:I)='-'
      ENDDO

C Get the number of parallel streams
      CALL UCLGSI('npara',NPARA,ISTAT)
      IF(NPARA.LT.1) THEN
         CALL UMSPUT('! Invalid number of parallel streams',
     :               1,0,ISTAT)
         GO TO 999
      ENDIF

C If the mode is single then we need to get the name 
C for the output as well as the position of its centre
      IF(MODE.EQ.'single') THEN
         CALL UCLGST('output',OUTDAT,ISTAT)
         CALL UCLGSD('outra',OUTRA,ISTAT)
         OUTRA=OUTRA*PIBY
         CALL UCLGSD('outdec',OUTDEC,ISTAT)
         OUTDEC=OUTDEC*PIBY
         SINGLE=.TRUE.

C Check for multiple streams in this mode
         IF(NPARA.NE.1) THEN
            CALL UMSPUT(
     :   '! Parallel output not possible for single mode',
     :      1,0,ISTAT)
            GO TO 9999
         ENDIF
      ELSE
         SINGLE=.FALSE.
      ENDIF

C Get the size of the output image
      CALL UCLGSI('outnx',NX,ISTAT)
      CALL UCLGSI('outny',NY,ISTAT)
 
C Get the reference pixel positions
      CALL UCLGSD('raref',RAREF,ISTAT)
      RAREF=RAREF*PIBY
      CALL UCLGSD('decref',DECREF,ISTAT)
      DECREF=DECREF*PIBY

C Get the output scale (in arcsecs/pix)
      CALL UCLGSD('outscl',SCALE,ISTAT)

C Get the coordinates of the central section (this is
C just an offset to avoid negative section numbers)
      CALL UCLGSI('xcen',XZ,ISTAT)
      CALL UCLGSI('ycen',YZ,ISTAT)

C Convert from arcsecs/pix to radians/pixel
      RPP=SCALE/(3600.0/PIBY)

C Get the name of the input directory
      CALL UCLGST('direct',DIR,ISTAT)
      IF(DIR.EQ.' ') THEN
	 LD=0
      ELSE
         CALL LENSTR(DIR,I1,I2)
         LD=I2-I1+1
         IF(I1.EQ.0 .AND. I2.EQ.0) LD=0
      ENDIF

C Check that the directory ends in either $ or /
C and if not put a / on
      IF(LD.GT.0 .AND. 
     :   DIR(I2:I2).NE.'/' .AND. DIR(I2:I2).NE.'$') THEN
         DIR(I2+1:I2+1)='/'
         I2=I2+1
         LD=LD+1
      ENDIF

C Get the name of the input images and initialise template processing
      CALL UCLGST('image',IMAGES,ISTAT)

C Get the name of the mask, this can be "default" or blank
      CALL UCLGST('weight',MASK,ISTAT)

      IF(LD.GT.0) THEN
         CALL TIMOTP(DIR(I1:I2)//IMAGES,IMLD,ISTAT)
      ELSE
         CALL TIMOTP(IMAGES,IMLD,ISTAT)
      ENDIF

      IF(ISTAT.NE.0) THEN
         GO TO 9999
      ELSE
         TOP=.TRUE.
      ENDIF

C If we are in single mode we need to find out the range
C of super-image pixel coordinates covered by the output image
      IF(SINGLE) THEN
         IF(PROJ.EQ.'COE') THEN
            CALL COEPRO(OUTRA,OUTDEC,1,RAREF,DECREF,RPP,
     :               XOUT,YOUT)
         ELSE
            CALL TANPRO(OUTRA,OUTDEC,1,RAREF,DECREF,RPP,
     :               XOUT,YOUT)
         ENDIF

         OXMIN=XOUT-(NX/2)
         OYMIN=YOUT-(NY/2)
         OXMAX=XOUT+(NX/2)
         OYMAX=YOUT+(NY/2)
      ENDIF

C Before we start looping we can write out all the things
C which never change - only in the single case
      IF(MODE.EQ.'single') THEN
         CALL UMSPUT('# '//VERS(1:48),1,0,ISTAT)
         WRITE(CHARS(1:38),'(''# Starting at '',A24)') 
     :         ICTIME(IITIME())
         CALL UMSPUT(CHARS(1:38),1,0,ISTAT)
         CALL UMSPUT('#',1,0,ISTAT)
         WRITE(CHARS,'(''drizzle.outnx='',I5)') NX
         CALL UMSPUT(CHARS,1,0,ISTAT)
         WRITE(CHARS,'(''drizzle.outny='',I5)') NY
         CALL UMSPUT(CHARS,1,0,ISTAT)
         WRITE(CHARS,'(''drizzle.raref='',F13.7)') RAREF/PIBY
         CALL UMSPUT(CHARS,1,0,ISTAT)
         WRITE(CHARS,'(''drizzle.decref='',F13.7)') DECREF/PIBY
         CALL UMSPUT(CHARS,1,0,ISTAT)
         WRITE(CHARS,'(''drizzle.outscl='',F13.7)') SCALE
         CALL UMSPUT(CHARS,1,0,ISTAT)
         CALL UMSPUT('drizzle.pixfrac=1.0',1,0,ISTAT)
         CALL UMSPUT('drizzle.out_un="counts"',1,0,ISTAT)
         CALL UMSPUT('drizzle.expkey="exptime"',1,0,ISTAT)
         CALL UMSPUT('drizzle.exp_sc=no',1,0,ISTAT)
         CALL UMSPUT('drizzle.wt_scl="header"',1,0,ISTAT)
         CALL UMSPUT('drizzle.fillval="indef"',1,0,ISTAT)
         CALL UMSPUT('drizzle.proj="'//PROJ//'"',1,0,ISTAT)
      ENDIF

C Get the names of the input data images
C Start of main loop
      DO WHILE(.TRUE.)
       CALL TIMXTP(IMLD,IMAGE(NIM),ISTAT)

       IF(ISTAT.EQ.USEOF) GO TO 888

       CALL UIMOPN(IMAGE(NIM),1,ID,ISTAT)
       IF(ISTAT.NE.0) THEN
         CALL UMSPUT('! Unable to access image '//IMAGE(NIM),
     :               1,0,IS)
         GO TO 999
       ENDIF

C Get size and shape
       CALL UIMGID(ID,DATTYP,NDIMS,DIMS,ISTAT)

C Get the EIS astrometric header information
       CALL EIGTCO(ID,CRVAL1,CRPIX1,CRVAL2,CRPIX2,
     :                  CDELT1,CDELT2,MPIX,FLXSCL,
     :                  XCO,YCO,RAOFF,DECOFF,NC,ISTAT)
       IF(ISTAT.NE.0) THEN
         CALL UMSPUT('! Unable to get EIS header information',
     :               1,0,ISTAT)
         GO TO 999
       ENDIF

C Close the image
       CALL UIMCLO(ID,ISTAT)

C Now work out where all 4 corners go in terms of pixel
C offsets from the centre of the superimage
      CALL CORNER(DIMS(1),DIMS(2),
     :   CRVAL1,CRPIX1,CRVAL2,CRPIX2,
     :   CDELT1,CDELT2,MPIX,
     :   XCO,YCO,RAOFF,DECOFF,NC,
     :   RAREF,DECREF,XPOFF,YPOFF,RPP,PROJ,
     :   XMIN(NIM),XMAX(NIM),YMIN(NIM),YMAX(NIM),ISTAT)

      IF(ISTAT.NE.0) THEN
         CALL UMSPUT('! Error calculating corner positions',
     :    1,0,ISTAT)
         GO TO 999
      ENDIF

      IF(SINGLE) THEN
         IF(XMAX(NIM).GT.OXMIN .AND. YMAX(NIM).GT.OYMIN .AND.
     :      XMIN(NIM).LT.OXMAX .AND. YMIN(NIM).LT.OYMAX) THEN
            OVER=.TRUE.
         ELSE
            OVER=.FALSE.
            CALL UMSPUT('#',1,0,ISTAT)
            CALL UMSPUT(
     :  '#-No overlap for input image: '//IMAGE(NIM),
     :   1,0,ISTAT)
            CALL UMSPUT('#',1,0,ISTAT)

            GO TO 777
         ENDIF
      ELSE
         OVER=.TRUE.
      ENDIF

C Now we have to handle the single and multiple cases
C differently. We can only handle the single option here
C and do the multiple output image options afterwards when
C we know all about the images
      IF(SINGLE) THEN
         CALL UMSPUT('#',1,0,ISTAT)
         CALL UMSPUT(
     :  '#---------Input image: '//IMAGE(NIM),
     :   1,0,ISTAT)

         IF(MASK.EQ.'default') THEN
            CALL LENSTR(IMAGE(NIM),I1,I2)

C Check for standard extensions (.fits and .hhh)
            IF(IMAGE(NIM)(I2-4:I2).eq.'.fits') THEN
               CALL UMSPUT(
     : 'drizzle.in_mask="'//IMAGE(NIM)(I1:I2-5)//'.weight.fits"',
     :               1,0,ISTAT)
            ELSE IF(IMAGE(NIM)(I2-3:I2).eq.'.hhh') THEN
               CALL UMSPUT(
     : 'drizzle.in_mask="'//IMAGE(NIM)(I1:I2-4)//'.weight.hhh"',
     :               1,0,ISTAT)
            ELSE
               CALL UMSPUT(
     : 'drizzle.in_mask="'//IMAGE(NIM)(I1:I2)//'.weight"',
     :               1,0,ISTAT)
            ENDIF

         ELSE IF(MASK.EQ.' ') THEN
            CALL UMSPUT('drizzle.in_mask=""',1,0,ISTAT)
         ELSE
            CALL LENSTR(MASK,I1,I2)
            CALL UMSPUT('drizzle.in_mask="'//MASK(I1:I2)//'"',
     :        1,0,ISTAT)
         ENDIF

         CALL UMSPUT('#',1,0,ISTAT)
         WRITE(CHARS,'(''drizzle.xpoff='',I6)') NINT(OXMIN)
         CALL UMSPUT(CHARS,1,0,ISTAT)
         WRITE(CHARS,'(''drizzle.ypoff='',I6)') NINT(OYMIN)
         CALL UMSPUT(CHARS,1,0,ISTAT)

         CALL LENSTR(OUTDAT,OI1,OI2)
         IF(OUTDAT(OI2-4:OI2).EQ.'.fits') THEN
            OI2=OI2-5
            EXT='.fits'
         ELSE IF(OUTDAT(OI2-3:OI2).EQ.'.hhh') THEN
            OI2=OI2-4
            EXT='.hhh '
         ELSE
            EXT='.fits'
         ENDIF

         CALL LENSTR(IMAGE(NIM),I1,I2)

C Flush the process cache to avoid resource problems
         CALL UMSPUT('flprc',1,0,ISTAT)

         IF(CON) THEN
           CALL UMSPUT('drizzle '//IMAGE(NIM)(I1:I2)//' '
     :       //OUTDAT(OI1:OI2)//EXT//' outweig='
     :       //OUTDAT(OI1:OI2)//'.weight'//EXT//' outcont='
     :       //OUTDAT(OI1:OI2)//'.context'//EXT,
     :               1,0,ISTAT)
         ELSE
           CALL UMSPUT('drizzle '//IMAGE(NIM)(I1:I2)//' '
     :       //OUTDAT(OI1:OI2)//EXT//' outweig='
     :       //OUTDAT(OI1:OI2)//'.weight'//EXT//' outcont=""',
     :               1,0,ISTAT)
         ENDIF
       ENDIF

777    CONTINUE
       NIM=NIM+1
      ENDDO

888   CONTINUE

      NIM=NIM-1

C Check on the total number of images
      IF(NIM.EQ.0) THEN
         CALL UMSPUT('! No input images were found',1,0,ISTAT)
         ISTAT=1
	 GO TO 999
      ELSE
         WRITE(CHARS,'(''# '',I6,'' input images found'')') NIM
	 CALL UMSPUT(CHARS,1,0,ISTAT)
         CALL UMSPUT(' ',1,0,ISTAT)
      ENDIF

      IF(.NOT.SINGLE) THEN

C There are now two options - tiles and bestfit
       IF(MODE.EQ.'tiles') THEN

        N=1
        DO K=1,NIM

C We now try to work out which images to coadd onto.
         NXST=NINT(FLOAT(XMIN(K))/NX)
         NXEN=NINT(FLOAT(XMAX(K))/NX)
         NYST=NINT(FLOAT(YMIN(K))/NY)
         NYEN=NINT(FLOAT(YMAX(K))/NY)

C Loop over all the images which the data will touch
         DO J=NYST,NYEN
           DO I=NXST,NXEN

              DATA(N)=IMAGE(K)
              XS(N)=I
              YS(N)=J

C Find range
              IF(I.GT.IMAX) IMAX=I
              IF(I.LT.IMIN) IMIN=I
              IF(J.GT.JMAX) JMAX=J
              IF(J.LT.JMIN) JMIN=J

              N=N+1
           ENDDO
         ENDDO
        ENDDO
      ELSE

C This is the "bestfit" option, new in v0.4

C First find the extreme values on the superimage
        XI=XMIN(1)
        XA=XMAX(1)
        YI=YMIN(1)
        YA=YMAX(1)
        DO K=1,NIM
           IF(XMIN(K).LT.XI) XI=XMIN(K)
           IF(XMAX(K).GT.XA) XA=XMAX(K)
           IF(YMIN(K).LT.YI) YI=YMIN(K)
           IF(YMAX(K).GT.YA) YA=YMAX(K)
        ENDDO

C Now work out the appropriate offsets so that the output images
C all fit neatly into the smallest possible set of output images
        IF(MOD((XA-XI+1),NX).EQ.0) THEN
           NXP=(XA-XI+1)/NX
        ELSE
           NXP=(XA-XI+1)/NX+1
        ENDIF
        IF(MOD((YA-YI+1),NY).EQ.0) THEN
           NYP=(YA-YI+1)/NY
        ELSE
           NYP=(YA-YI+1)/NY+1
        ENDIF

C Now we work out an offset in pixels from the reference pixel
C in the superimage to the centre of the central section in the
C bestfit area
        IF(MOD(NXP,2).EQ.0) THEN
           XPIXOF=(XA+XI)/2-NX/2
        ELSE
           XPIXOF=(XA+XI)/2
        ENDIF

        IF(MOD(NYP,2).EQ.0) THEN
           YPIXOF=(YA+YI)/2-NY/2
        ELSE
           YPIXOF=(YA+YI)/2
        ENDIF

C We now try to work out which images to coadd onto.
        N=1
        DO K=1,NIM
         NXST=NINT(FLOAT(XMIN(K)-XPIXOF)/NX)
         NXEN=NINT(FLOAT(XMAX(K)-XPIXOF)/NX)
         NYST=NINT(FLOAT(YMIN(K)-YPIXOF)/NY)
         NYEN=NINT(FLOAT(YMAX(K)-YPIXOF)/NY)

C Loop over all the images which the data will touch
         DO J=NYST,NYEN
           DO I=NXST,NXEN

              DATA(N)=IMAGE(K)
              XS(N)=I
              YS(N)=J

C Find range
              IF(I.GT.IMAX) IMAX=I
              IF(I.LT.IMIN) IMIN=I
              IF(J.GT.JMAX) JMAX=J
              IF(J.LT.JMIN) JMIN=J

              N=N+1
           ENDDO
         ENDDO
        ENDDO
      ENDIF

      N=N-1

C Now we work through the values again and write
C out what is needed, this is the same for "bestfit" and "tiles"
      IOLD=-999
      JOLD=-999
      NSEC=0
      NP=0

      DO J=JMIN,JMAX
       DO I=IMIN,IMAX
        DO K=1,N
         IF(XS(K).EQ.I .AND. YS(K).EQ.J) THEN

C Apply the offsets
          II=I+XZ
          JJ=J+YZ

          IF(II.LT.0) THEN
            XSIGN='-'
          ELSE
            XSIGN='+'
          ENDIF

          IF(JJ.LT.0) THEN
            YSIGN='-'
          ELSE
            YSIGN='+'
          ENDIF

          IA=ABS(II)
          JA=ABS(JJ)

          CALL LENSTR(DATA(K),I1,I2)

C If we have moved to a new section let everyone know
          IF(IOLD.NE.II .OR. JOLD.NE.JJ) THEN

C Increment section counter
             NSEC=NSEC+1

C Increment parallel stream counter and check for wrap-around
             NP=NP+1
             IF(NP.GT.NPARA) NP=NP-NPARA

C If we are still in the first loop we need to open the file
             IF(NP.EQ.NSEC) THEN
                WRITE(PNAME,'(A4,''-'',I3,''.cl'')') STEM,NP
                DO KK=6,8
                   IF(PNAME(KK:KK).EQ.' ') PNAME(KK:KK)='0'
                ENDDO

                OPEN(10+NP,FILE=PNAME,IOSTAT=ISTAT)
                IF(ISTAT.NE.0) THEN
                   CALL UMSPUT('! Unable to open output file '
     :                         //PNAME,1,0,ISTAT)
                   GO TO 999
                ELSE
                   CALL UMSPUT(
     :        '# Writing drizzle commands to script file '
     :                         //PNAME,1,0,ISTAT)
                ENDIF

C Write a header appropriately
C Before we start looping we can write out all the things
C which never change
      WRITE(10+NP,'(''# '',A48)') VERS
      WRITE(10+NP,'(''# Starting at '',A24)') ICTIME(IITIME())
      WRITE(10+NP,'(A)')'#'
      WRITE(CHARS,'(''drizzle.outnx='',I5)') NX
      WRITE(10+NP,'(A)')CHARS(1:19)
      WRITE(CHARS,'(''drizzle.outny='',I5)') NY
      WRITE(10+NP,'(A)')CHARS(1:19)
      WRITE(CHARS,'(''drizzle.raref='',F13.7)') RAREF/PIBY
      WRITE(10+NP,'(A)')CHARS(1:28)
      WRITE(CHARS,'(''drizzle.decref='',F13.7)') DECREF/PIBY
      WRITE(10+NP,'(A)')CHARS(1:28)
      WRITE(CHARS,'(''drizzle.outscl='',F13.7)') SCALE
      WRITE(10+NP,'(A)')CHARS(1:28)
      WRITE(10+NP,'(A)')'drizzle.pixfrac=1.0'
      WRITE(10+NP,'(A)')'drizzle.out_un="counts"'
      WRITE(10+NP,'(A)')'drizzle.expkey="exptime"'
      WRITE(10+NP,'(A)')'drizzle.exp_sc=no'
      WRITE(10+NP,'(A)')'drizzle.wt_scl="header"'
      WRITE(10+NP,'(A)')'drizzle.fillval="indef"'
      WRITE(10+NP,'(A)')'drizzle.proj="'//PROJ//'"'

             ENDIF

             WRITE(10+NP,
     :  '(''#--------------------------------------'')')
           WRITE(CHARS,'(''print "Starting new coadded section: '''//
     :          ',A4,''CA'',A1,I2,A1,I2,''"'')')
     :          STEM,XSIGN,IA,YSIGN,JA
             IF(CHARS(45:45).EQ.' ') CHARS(45:45)='0'
             IF(CHARS(48:48).EQ.' ') CHARS(48:48)='0'
             WRITE(10+NP,'(A)') CHARS(1:80)   
             IOLD=II
             JOLD=JJ
          ENDIF

C Write the input weight image name
          IF(MASK.EQ.'default') THEN

C Check for standard extensions (.fits and .hhh)
            IF(DATA(K)(I2-4:I2).EQ.'.fits') THEN
               WRITE(10+NP,'(A)') 
     : 'drizzle.in_mask="'//DATA(K)(I1:I2-5)//'.weight.fits"'
            ELSE IF(IMAGE(NIM)(I2-3:I2).EQ.'.hhh') THEN
               WRITE(10+NP,'(A)')
     : 'drizzle.in_mask="'//DATA(K)(I1:I2-4)//'.weight.hhh"'
            ELSE
               WRITE(10+NP,'(A)')
     : 'drizzle.in_mask="'//DATA(K)(I1:I2)//'.weight"'
            ENDIF

          ELSE IF(MASK.EQ.' ') THEN
               WRITE(10+NP,'(A)') 'drizzle.in_mask=""'
          ELSE
            CALL LENSTR(MASK,I1,I2)
            WRITE(10+NP,'(A)') 'drizzle.in_mask="'//MASK(I1:I2)//'"'
          ENDIF

          WRITE(CHARS,'(''drizzle.xpoff='',I6)') I*NX+XPIXOF-NX/2
          WRITE(10+NP,'(A)') CHARS(1:20)
          WRITE(CHARS,'(''drizzle.ypoff='',I6)') J*NY+YPIXOF-NY/2
          WRITE(10+NP,'(A)') CHARS(1:20)

C Flush process cache
          WRITE(10+NP,'(''flprc'')')

          CHARS(1:9+I2-I1)='drizzle '//DATA(K)(I1:I2)
          IF(CON) THEN
CCC           WRITE(CHARS,'(''drizzle '',A<I2-I1+1>,1X,A4,''CA'',
           WRITE(CHARS(10+I2-I1:),'(1X,A4,''CA'',
     :                     A1,I2,A1,I2,''.hhh'',
     :          '' outweig='',A4,''CA'',A1,I2,A1,I2,
     :          ''.weight.hhh'',
     :          '' outcont='',A4,''CA'',A1,I2,A1,I2,
     :          ''.context.hhh'')')
CCC     :          DATA(K)(I1:I2),STEM,XSIGN,IA,YSIGN,JA,
     :          STEM,XSIGN,IA,YSIGN,JA,
     :          STEM,XSIGN,IA,YSIGN,JA,
     :          STEM,XSIGN,IA,YSIGN,JA
          ELSE
CCC           WRITE(CHARS,'(''drizzle '',A<I2-I1+1>,1X,A4,''CA'',
           WRITE(CHARS(10+I2-I1:),'(1X,A4,''CA'',
     :                     A1,I2,A1,I2,''.hhh'',
     :          '' outweig='',A4,''CA'',A1,I2,A1,I2,
     :          ''.weight.hhh outcont=""'')')
     :          STEM,XSIGN,IA,YSIGN,JA,
CCC     :          DATA(K)(I1:I2),STEM,XSIGN,IA,YSIGN,JA,
     :          STEM,XSIGN,IA,YSIGN,JA
          ENDIF

C Fill gaps
          IP=I2-I1+18
          IF(CHARS(IP:IP).EQ.' ') CHARS(IP:IP)='0'
          IF(CHARS(IP+3:IP+3).EQ.' ') CHARS(IP+3:IP+3)='0'
          IF(CHARS(IP+25:IP+25).EQ.' ') CHARS(IP+25:IP+25)='0'
          IF(CHARS(IP+28:IP+28).EQ.' ') CHARS(IP+28:IP+28)='0'
          IF(CON) THEN
             IF(CHARS(IP+57:IP+57).EQ.' ') CHARS(IP+57:IP+57)='0'
             IF(CHARS(IP+60:IP+60).EQ.' ') CHARS(IP+60:IP+60)='0'
          ENDIF

          WRITE(10+NP,'(A)') CHARS(1:IP+73)
          WRITE(10+NP,'(''#---'')')
          CHARS(1:I2-I1+1)=DATA(K)(I1:I2)
          WRITE(CHARS(I2-I1+2:),'('' -> '',
     :                  A4,''CA'',A1,I2,A1,I2,'' [Stream '',
     :                  I2,'']'')')
     :         STEM,XSIGN,IA,YSIGN,JA,NP
               IF(CHARS(IP-5:IP-5).EQ.' ') CHARS(IP-5:IP-5)='0'
               IF(CHARS(IP-2:IP-2).EQ.' ') CHARS(IP-2:IP-2)='0'
          ENDIF
         ENDDO
        ENDDO
       ENDDO
      ENDIF

C Close the image template matching
999   CONTINUE

C This line had a mistake before May 2002
      DO I=11,10+NPARA
         CLOSE(I)
      ENDDO

      IF(TOP) CALL TIMCTP(IMLD,ISTAT)

C Finally list out the basic numbers and sizes of what is going to 
C happen - only in tiles mode
      IF(.NOT.SINGLE) THEN
         CALL UMSPUT(' ',1,0,ISTAT)
         WRITE(CHARS,'(''# Number of input images is '',I5)') NIM
         CALL UMSPUT(CHARS,1,0,ISTAT)
         WRITE(CHARS,'(''# Number of drizzle operations is '',I5)') N
         CALL UMSPUT(CHARS,1,0,ISTAT)
         WRITE(CHARS,'(''# Number of parallel streams is '',I5)') NPARA
         CALL UMSPUT(CHARS,1,0,ISTAT)
         WRITE(CHARS,'(''# Number of output sections is '',I5)') NSEC
         CALL UMSPUT(CHARS,1,0,ISTAT)
         IF(CON) THEN
            WRITE(CHARS,
     : '(''# Total space required for output sections is '',F10.3,
     :   '' (Mb)'')') 
     :    FLOAT(NSEC*NX*NY)*4.0*2.5/FLOAT(1024*1024)
         ELSE
            WRITE(CHARS,
     :    '(''# Total space required for output sections is '',F10.3,
     :      '' (Mb)'')')
     :    FLOAT(NSEC*NX*NY)*4.0*2.0/FLOAT(1024*1024)
         ENDIF
         CALL UMSPUT(CHARS,1,0,ISTAT)

         CALL UMSPUT(' ',1,0,ISTAT)
      ENDIF

9999  CONTINUE
      RETURN
      END