c  this file contains subroutines inputdeck, ordermap, inform, and 
c  output
c  
c:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
c  
c  
c  

      subroutine         inputdeck(inimage,outimage,xmaxfile,error)
c  
c  
c  
c.......................................................................
c  
c  author: Christopher W. Churchill
c
c  history: 20-Jan94   created, CWC
c           05-Mar94   added standard output warnings, error traps, CWC
c
c  purpose: load the user input deck, run simple checks, issue warnings
c
c  input:   the input deck file, which has list directed format:
c
c           inimage, outimage
c           xmaxfile, norders
c           xorder, yorder, gamma
c           rdnoise, gain
c           reject, toler
c           minflag, corrflag
c
c           inimage  = IRAF image to be background corrected
c           outimage = IRAF background surface
c           xmaxfile = the file containing the output of APMAP
c           norders  = the number of orders in inimage
c           xorder   = the polynomial order of the x terms
c           yorder   = the polynomial order of the y terms
c           gamma    = the power law scaling
c           rdnoise  = the system read noise [e-]
c           gain     = the system gain [DN/e-]
c           reject   = the sigma rejection factor 
c           toler    = the convergence tolerence for the fit
c           minflag  = (logical) flag for grid optimization
c           corrflag = (logical) flag for local term correction
c
c  output:  error    = (logical) input error flag
c
c  I/O: input file and standard output on error
c
c  options: none
c
c  routines called: none
c
c  description: this routine reads in the input deck and performs
c               several consistency and stupidity checks
c  
c.......................................................................
c  
      implicit undefined (a-z)
      logical            error
      character*80       deckname,inimage,outimage,xmaxfile
      include            'hamscatt.par'
      include            'hamscatt.com'
c
c
c
c
c  
c  query for the input deck
 11   write(stdout,10) 
      read(stdin,*,err=11) deckname
 10   format(1x,'HAMSCATT: Enter Input Deck Filename : ',$)
c
c  open the input deck, read and close
      open(unit=1,file=deckname,status='old',err=1001)
      read(1,*,err=1002) inimage, outimage
      read(1,*,err=1003) xmaxfile, norders
      read(1,*,err=1004) xorder, yorder, gamma
      read(1,*,err=1005) rdnoise, gain
      read(1,*,err=1006) reject, toler
      read(1,*,err=1007) minflag, corrflag
      close(unit=1)
c  
c  communicate input deck
      write(stdout,*) ' ==============================================='
      write(stdout,600) inimage,outimage,xmaxfile,norders,xorder,
     @                  yorder,gamma,rdnoise,gain,reject,toler
c  
c  useful checks on the input; error trapping
c  xorder & yorder > 0 ; gamma>0 ; readnoise>=0 ; gain >=1
      if (xorder.lt.1) then
       write(stdout,*) ' ERROR: bad input parameter - XORDER >=1'
       error = .true.
      end if
      if (yorder.lt.1) then
       write(stdout,*) ' ERROR: bad input parameter - YORDER >=1'
       error = .true.
      end if
      if (gamma.lt.0.d0) then
       write(stdout,*) ' ERROR: bad input parameter - GAMMA > 0'
       error = .true.
      end if
      if (rdnoise.lt.0.d0) then
       write(stdout,*) ' ERROR: bad input parameter - RDNOISE >= 0'
       error = .true.
      end if
      if (gain.lt.1.0d0) then
       write(stdout,*) ' ERROR: bad input parameter - GAIN >= 1'
       error = .true.
      end if
      if (reject.lt.0.) then
       write(stdout,*) ' ERROR: bad input parameter - REJECT > 0'
       error = .true.
      end if
c
c  issue warnings
      if (rdnoise.lt.gain) then
       write(stdout,*) ' WARNING: parameter GAIN > RDNOISE'
      end if
      if (rdnoise.eq.0.) then
       write(stdout,*) ' WARNING: parameter RDNOISE = 0.0'
      end if
      write(stdout,*) ' -----------------------------------------------'
c  
c  communicate option settings
c  set corrflag low if gamma=0 as a consistency measure
      if (minflag) then
       write(stdout,601)
      else
       write(stdout,602)
      end if
      if (gamma.ne.0) then
       write(stdout,603)
      else
       corrflag = .false.
       write(stdout,604)
      end if
      if (corrflag) then
       write(stdout,605)
      else
       write(stdout,606)
      end if
      write(stdout,*) ' ==============================================='
c  
      error = .false.
      return
c  
c  error=.true. ; communicate error and bail
 1001 write(stdout,*) ' ERROR: cannot open input deck'
      error = .true.
      goto 11
 1002 write(stdout,*) ' ERROR: input or output image names'
      write(stdout,*) '        incorrectly formatted'
      error = .true.
      return
 1003 write(stdout,*) ' ERROR: aperture database file or NORDERS'
      write(stdout,*) '        parameter incorrectly formatted'
      error = .true.
      return
 1004 write(stdout,*) ' ERROR: XORDER, YORDER, or GAMMA prameters'
      write(stdout,*) '        incorrectly formatted'
      error = .true.
      return
 1005 write(stdout,*) ' ERROR: RDNOISE or GAIN parameters'
      write(stdout,*) '        incorrectly formatted'
      error = .true.
      return
 1006 write(stdout,*) ' ERROR: sigma rejection factor or fitting '
      write(stdout,*) '        tolerance TOLER parameter'
      write(stdout,*) '        incorrectly formatted'
      error = .true.
      return
 1007 write(stdout,*) ' ERROR: Option flags (logical)'
      write(stdout,*) '        incorrectly formatted'
      error = .true.
      return
c  
c  formats
 600  format(1x,'Input Image  = ',a35,/,
     @       1x,'Output Image = ',a35,/,
     @       1x,'Order Map    = ',a35,/,
     @       1x,'norders      = ',i4,/,
     @       1x,'xorder (XD)  = ',i3,/,
     @       1x,'yorder       = ',i3,/,
     @       1x,'gamma        = ',f6.2,/,
     @       1x,'read-noise   = ',f6.2,/,
     @       1x,'gain         = ',f6.2,/,
     @       1x,'sigma reject = ',f6.2,/,
     @       1x,'tolerance    = ',1pe10.2)
 601  format(1x,'Option: interorder fitting     =  ON')
 602  format(1x,'Option: optimized fitting      =  OFF')
 603  format(1x,'Option: local term fitting     =  ON')
 604  format(1x,'Option: local term fitting     =  OFF')
 605  format(1x,'Option: local term correction  =  ON')
 606  format(1x,'Option: local term correction  =  OFF')
c  
c  
      end
c  
c  
c  
c.......................................................................
c  
c:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
c
c
c
      subroutine ordermap(xmaxfile,error)
c  
c
c
c.......................................................................
c  
c  author: Christopher W. Churchill
c
c  history: 02-Oct94   moved from main routine, CWC
c
c  purpose: load the order map file, check for incompatibilities with
c           physical sizes. returns ndata, the number of interorder
c           troughs ; and m, the number of data points to fit (prior
c           to sigma rejection)
c
c  input:   xmaxfile  = the file containing the output of APMAP
c
c  output:  error  = (logical) error flag 
c           xmax   = order map [pixel locations] (APEXTRACT tracings)
c
c  I/O: input of aperture database file (output of apmap) and
c       standard output on error
c
c  options: none
c
c  routines called: none
c
c  description: see subroutine comments ; this routine returns to vital
c               numbers, m and ndata.
c  
c.......................................................................
c
      implicit undefined (a-z)
      logical            error
      integer            i,j,k,col
      character*(*)      xmaxfile      
      include            'hamscatt.par'
      include            'hamscatt.com'      
c
c
c
c
c
c  open the file, error trap
      open(unit=1,file=xmaxfile,status='old',err=10)
c
c  loop over maxorders to check consistency with input norders
      do 01 j=1,maxorders
       do 02 k=1,ncols
        read(1,*,end=3,err=11) i,col,xmax(j,k)
 02    continue
 01   continue
      close(unit=1)
c 
c  enforce consistency between input deck and order map
 03   if (i.ne.norders) then
       write(stdout,*) 'WARNING: first column of order map corrupt?'
       write(stdout,*) 'number of orders does not agree with input.' 
      end if
      if ((j-1).ne.norders) then
       write(stdout,*) 'ERROR: parameter NORDERS does not agree'
       write(stdout,*) 'with number of orders tabulated in the'
       write(stdout,*) 'order map file'
       error = .true.
       return
      end if

c  all is well ; normal return
      return
c
c  error traps
 10   write(stdout,*) ' ERROR: cannot open order map file'
      error = .true.
      return
 11   write(stdout,*) ' ERROR: order map file corrupt?'
      error = .true.
      return
c
      end
c
c.......................................................................
c
c:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
c  
c  
c  
      subroutine         inform(info,error)
c  
c  
c  
c.......................................................................
c  
c  author: Christopher W. Churchill
c
c  history: 28-Aug94   created, CWC
c
c  purpose: interpret and communicate DNLS1 convergence check parameter
c           INFO
c
c  input:   info   = value of convergence check flag
c
c  output:  error  = (logical) error flag if info=0
c
c  I/O: standard output 
c
c  options: none
c
c  routines called: none
c
c  description: The accuracy of DNLS1 is controlled by the convergence 
c               parameter TOL.  This parameter is used in tests which 
c               make three types of comparisons between the current
c               solution and the last solution.  DNLS1 terminates when 
c               any of one of the tests is satisfied.  If TOL is less 
c               than the machine precision (as defined by the function 
c               R1MACH(4)), then DNLS1 attempts only to satisfy the test
c               to the accuracy of machine precision.  Further progress
c               is not usually possible.  This routine gives a brief
c               message as to which test has been satisfied.  Further
c               descriptions of these tests are in the dnls1.f comment
c               "cards" (ha ha!)
c  
c.......................................................................
c  
c  declarations
      implicit undefined (a-z)
      logical            error
      integer            info
      include            'hamscatt.par'
c  
c
c
c
c
c  the following messages are paraphrased from the dnl1se instruction
c  header, commented in the dnl1se.f source code file
c  
      if (info.eq.0) then
       write(stdout,*) ' improper input parameters.'
       error = .true.
      end if
c  
      if (info.eq.1) then
       write(stdout,*) ' **** 1st convergence condition satisfied ****'
       write(stdout,*) ' - algorithm estimates that the relative error'
       write(stdout,*) '   in the sum of squares is at most the input'
       write(stdout,*) '   tolerance'
      end if
c  
      if (info.eq.2) then
       write(stdout,*) ' **** 2nd convergence condition satisfied ****'
       write(stdout,*) ' - algorithm estimates that the relative error'
       write(stdout,*) '   between the last and current solutions is at'
       write(stdout,*) '   most the input tolerance'
      end if
c  
      if (info.eq.3) then
       write(stdout,*) ' * 1st & 2nd convergence condition satisfied *'
       write(stdout,*) ' - algorithm estimates that the relative error'
       write(stdout,*) '   in the sum of squares is at most the input'
       write(stdout,*) '   tolerance'
       write(stdout,*) ' - algorithm estimates that the relative error'
       write(stdout,*) '   between the last and current solutions is at'
       write(stdout,*) '   most the input tolerance'
      end if
c  
      if (info.eq.4) then
       write(stdout,*) ' **** 3rd convergence condition satisfied ****'
       write(stdout,*) ' - model vector is orthogonal to the columns'
       write(stdout,*) '   of the Jacobian to machine precision'
       write(stdout,*) ' - solution should be viewed with CAUTION'
      end if
c  
      if (info.eq.5) then
       write(stdout,*) ' *** convergence conditions not satisfied ***'
       write(stdout,*) ' maximum allowable calls to evaluate the model'
       write(stdout,*) ' vector exceeded'
      end if
c  
      if (info.eq.6) then
       write(stdout,*) ' *** convergence conditions not satisfied ***' 
       write(stdout,*) ' - input tolerance is too small, no further'
       write(stdout,*) '   reduction in the sum of squares is possible'
      end if
c  
      if (info.eq.7) then
       write(stdout,*) ' *** convergence conditions not satisfied ***' 
       write(stdout,*) ' - input tolerance is too small, no further'
       write(stdout,*) '   improvement in the reduction of the solution'
       write(stdout,*) '   is possible'
      end if
c  
      write(stdout,*) ' ==============================================='
c  
c  
      return
      end
c  
c.......................................................................
c
c:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
c
c
c
      subroutine     output(chisq,nu,chisqnu,varfit,vardata,m,a,n,nrej,
     @                      outimage)
c
c
c
c.......................................................................
c
c  author: Christopher W. Churchill
c
c  history: 04-Oct94   moved from other subroutines, CWC
c
c  purpose:  standard output of fitting stats and two output files
c            (1) a "DQF" file with rejected pixel locations
c            (2) a ".fit" file with the parameter results
c
c input: you figure it out
c
c output: ditto
c
c I/O: lots of it
c
c options: none
c 
c routines called: none
c
c description: who cares
c
c.......................................................................
c
      implicit undefined (a-z)
      integer            i,j,unit,nrej,icol,m,n 
      double precision   chisq,nu,chisqnu,varfit,vardata,a(1)
      character*(*)      outimage
      character*25       maskfile,fitfile
      include            'hamscatt.par'
      include            'hamscatt.com'

c
c
c
c
c
c  make mask and fitting output filenames
      do 17 j=1,len(outimage)
       if (outimage(j:j).ne.' ') i=j
 17   continue
      maskfile = outimage(1:i)//'.dqf' 
      fitfile  = outimage(1:i)//'.fit'  
c
c
c  fitting statistics to standard output
      write(stdout,602) chisqnu,varfit,vardata,chisq,nu,nrej
      write(stdout,*) ' ==============================================='
c
c
c  write fitting statistics and resulting fit coefficients to file
c  this file could be incorporated as an initial guess for a similar 
c  image ; HAMSCATT can be easily modified 
      unit = 1
      open(unit=unit,file=fitfile,status='unknown')
      write(unit,202)
      write(unit,203) varfit,chisq,nu,chisqnu,vardata
      do 04 j=1,n
       write(unit,*) j,a(j)
 04   continue
      close(unit=unit)
c
 202  format(1x,t5,'varfit',t15,'chisq',t25,'nu',t35,'chisqnu',
     @       t45,'vardata')
 203  format(1x,1p5e10.2,/)
 602  format(1x,'Reduced Chi-Squared = ',1pe11.3,/,
     @       1x,'Variance (FIT)      = ',1pe11.3,/,
     @       1x,'Variance (DATA)     = ',1pe11.3,/,
     @       1x,'Chi-Squared         = ',1pe11.3,/,
     @       1x,'Degrees of Freedom  = ',1pe11.3,/,
     @       1x,'Rejected Pixels     = ',i5)
c  
c
c  output ASCII mask ; this may be modified to be an IRAF DQF file
      unit=2
      open(unit=unit,file=maskfile,status='unknown')
      write(unit,99) reject
      if (m.eq.(ncols*ndata)) then
       write(unit,*) ' NO REJECTED PIXELS'
      else
       do 13 j=1,ndata
       do 15 icol=1,ncols
        if (sig(j,icol).eq.0) write(unit,100) icol,int(xmin(j,icol))
 15    continue
 13    continue
      end if
      close(unit=unit)
c
 99   format(1x,'     Pixel Rejection Mask',/,
     @       1x,'   sigma rejection = ',f5.2,/,/,
     @       1x,'  col',2x,'  row')
 100  format(1x,i5,2x,i5)
c
c
      return
      end
c  
c.......................................................................
c..eof