! $Id$
! $Name$
!
! AUTHOR: C. Guirao
!
! NAME: baches_calib.prg
!   Semi-automatic wavelength procedure for echelle spectra BACHES.
!   The calibration is based in the ECHELLE package of MIDAS which assumes
!   raw-images oriented so that orders appears along the rows, withe low orders 
!   up and wavelength increases from left to right, i.e. red is up and left.
!   START and STEP descriptors of the image must be positive.
!
! INPUT:
!   - 2-dimensional FITS/BDF image with the spectrum of flat field lamp.  
!   - 2-dimensional FITS/BDF image with the spectrum of a calibration lamp.
!   - MIDAS table of the calibration lamp in FITS format. Default "thar.fit"
!   - Number of orders to be calibrated. Default 26
!   - Slit width. Default 10
!   - Final tolerance on RMS error. Default 0.3
!   - Final degree of polynomial defining dispersion. Default 3
!
! OUTPUT:
!   - ECHELLE session files for BACHES. Result of the "SAVE/ECHE baches"
!     generates: bachesLINE.tbl and bachesORDE.tbl
!     
!
! SUFIX CONVENTION:
!   b with bias subtraction
!   d with dark-current subtraction
!   f with flat-field division
!   1 averaged to one-dimensional image
!   w calibrated in wavelength
!   r rebined in wavelength
!   s calibrated in flux with a standard star
!   m all orders merged 
!
! USAGE: @@ baches_calib <ff> <lamp> [table] [orders] [width] [tol] [poly]
!
! EXAMPLES:
!      @@ baches_calib ff.fit th               (26 orders, 10 width, 0.2 tol, 3 poly)
!      @@ baches_calib ff.fit th ? 25 5 0.5 4  (25 orders, 5 width, 0.5 tol, 5 poly)
!      @@ baches_calib ff.fit th neon 5 5 0.5 4 (neon table, 5 orders, 5 width, 0.5 tol, 4 poly)
!
!*******************************************************************************
! ECHO/FULL
!
! INPUT PARAMS:
! p1 REQUIRED: FITS or BDF 2-dimension of image of spectrum of a flat field lamp
! p2 REQUIRED: FITS or BDF 2-dimension of image of spectrum of a calibration lamp
! p3 OPTIONAL: Table containing lines of calib. lamp.
! p4 OPTIONAL: Orders to be calibrated (starting with lower orders, red)
! p5 OPTIONAL: Slit width in pixels 
! p6 OPTIONAL: Final tolerance on RMS error
! p7 OPTIONAL: Final degree of polynomial defining dispersion
!
! Edit here to change defaults:
!
 DEFINE/PARA p1 ? I "Enter 2-dimensional image with the spectrum of a flat field lamp? "
 DEFINE/PARA p2 ? I "Enter 2-dimensional image with the spectrum of a calibration lamp? "
 DEFINE/PARA p3 thar.fit  TAB
 DEFINE/PARA p4 26  NUM
 DEFINE/PARA p5 10  NUM
 DEFINE/PARA p6 0.3 NUM
 DEFINE/PARA p7 3   NUM
!
! Convert paramaters into variables:
!
 DEFINE/LOCA cmd/c/1/20         "baches_calib"
 DEFINE/LOCA ff_frame/c/1/40    {p1}
 DEFINE/LOCA lamp_frame/c/1/40  {p2}
 DEFINE/LOCA calib_table/c/1/8  {p3}
 DEFINE/LOCA orders/i/1/1       {p4}
 DEFINE/LOCA width/i/1/1        {p5}
 DEFINE/LOCA tolerance/r/1/1    {p6}
 DEFINE/LOCA poly/i/1/1         {p7}
!
! Other local definitions
!
 DEFINE/LOCA frame/c/1/40
 DEFINE/LOCA basename/c/1/40    "" ? +lower_levels
 DEFINE/LOCA ext/c/1/10         "" ? +lower_levels
 DEFINE/LOCA session/c/1/40	"baches"
 DEFINE/LOCA session_file1/c/1/40
 DEFINE/LOCA session_file2/c/1/40
 DEFINE/LOCA tmpname/c/1/40
 DEFINE/LOCA total_first/i/1/1   0        
 DEFINE/LOCA total_last/i/1/1    0    
 DEFINE/LOCA iteration/i/1/1     0    

!
! MAIN, main:
!
 WRITE/OUT
 WRITE/OUT "PARAMETERS FOR THIS CALIBRATION:"
 WRITE/OUT "==============================="
 WRITE/OUT "Flat field  = {ff_frame}"
 WRITE/OUT "Calibration lamp = {lamp_frame}"
 WRITE/OUT "Calibration table = {calib_table}"
 WRITE/OUT "Num. of orders  = {orders}"
 WRITE/OUT "Slit width = {width}"
 WRITE/OUT "Tolerance on RMS = {tolerance}"
 WRITE/OUT "Polinomyal degree  = {poly}"
 WRITE/OUT
!
! Ask before any remove
!
! WRITE/OUT "{cmd}: cleaning previous calibration files"
WRITE/KEYW reply/c/1/1 y
INQUIRE/KEYW reply "{cmd}: Do you want to continue [yn] (y)?"
IF reply(:1) .eq. "y" THEN
  GOTO continue
ELSE
  RETURN/EXIT
ENDIF
!
! Clean
!
continue:
 @@ clean.prg

!
! Set context ECHELLE, if not yet already
!
 SHOW/CONT echelle > NULL
 IF  OUTPUTI(5) .eq. 0 THEN
   SET/CONTEXT echelle >NULL
 ENDIF
 INITIA/ECHE 
 RESET/DIS
!
! Check if an order identification exist from a previous session
!
 session_file1 = session + "IORDE.fit"
 tmpname = session + "ORDE"
 IF m$exist(session_file1) .eq. 1 THEN
   WRITE/OUT "{cmd}: An order identification already done in a previus session"
   WRITE/KEYW reply/c/1/1 y
   INQUIRE/KEYW reply "{cmd}: Do you want to use it [yn] (y)?"
   IF reply(:1) .eq. "y" THEN
     INDISK/FITS {session_file1} {tmpname}
     INITIA/ECHE {session}
     GOTO continue_with_line_identification
   ENDIF
 ENDIF

!
! Start creating a display
!
 CREATE/DISP 0 590,390,0,25
 LOAD/LUT heat

!
! Gets basename and extension of input file
!
 @@ file_ext {ff_frame}
 frame = basename
!
! Convert FITS format into MIDAS BDF format
!
 IF ext .ne. ".bdf" THEN
   @@ compute_fits {ff_frame} {frame}
 ENDIF
 CUTS/IMAG {frame} =sigma
 DISPLA/ECHE {frame}
 SET/ECHE widthi = {width}
 SET/ECHE nbordi = {orders}
 DEFINE/HOUG {frame}
 SAVE/ECHE {session} 
 session_file1 = session + "ORDE.tbl"
 tmpname = session + "IORDE.fit"
 OUTDISK/FITS {session_file1} {tmpname}
!
! Show current ECHELLE parameters
! before asking to continue
!
 SHOW/ECHE w
 WRITE/KEYW reply/c/1/1 y
 WRITE/OUT
 WRITE/OUT "{cmd}: Order identification finished"
 INQUIRE/KEYW reply "{cmd}: Do you want to continue [yn] (y)?"
 IF reply(:1) .eq. "y" THEN
   GOTO continue_with_line_identification
 ELSE
   RETURN/EXIT
 ENDIF

continue_with_line_identification:
!
! Create two displays: 1 for reference lamp, 0 for input lamp
!
 CLEAR/CHAN OVER
 CREATE/DISP 0 590,390,0,25
 LOAD/LUT heat
 CREATE/GRAP 0 590,440,590,25
 CREATE/DISP 1 590,390,590,25
 INDISK/FITS th_ref.fit th_ref.bdf
 INDISK/FITS th_ref_ident.fit th_ref_ident.tbl
 LOAD/LUT heat
 CUTS/IMAG th_ref =sigma
 DISPLA/ECHE th_ref
 CLEAR/CHAN OVER
 LOAD/TABLE th_ref_ident :X :Y :IDENT 0 6
!
 ASSIGN/DISP 0
!
! Gets basename and extension of calibration lamp frame
!
 @@ file_ext {lamp_frame}
 frame = basename
!
! Convert FITS format into MIDAS BDF format
!
 IF ext .ne. ".bdf" THEN
   @@ compute_fits {lamp_frame} {frame}
 ENDIF
 LOAD/LUT heat
! DISPLA/ECHE th
 CUTS/IMAG {frame} =sigma
 DISPLA/ECHE {frame}
!
 SET/ECHE wlc = {frame}
 SET/ECHE slit = {width}
 EXTRAC/ECHE {wlc} extwlc
 thres2 = 10.0
 SEARCH/ECHE extwlc 
!
 LOAD/SEAR
 SAVE/ECHE baches
 CLEAR/CHAN OVER
!
!
! Gets basename and extension of calibration lamp frame
!
 @@ file_ext {calib_table}
 frame = basename
!
! Convert FITS format into MIDAS BDF format
!
 IF ext .ne. ".bdf" THEN
  INDISK/FITS {calib_table} {frame}
 ENDIF
 LOAD/LUT heat
 SET/ECHE lincat = {frame}.tbl
!
! Start with tol=2 and dc=3 (poly)
!
 SET/ECHE tol = 2
 wlciter(2) = 0.2
 SET/ECHE dc = 3

!
! First identification using reference file
!
 IDENTI/ECH 
 SET/ECHE guess=baches
 SAVE/ECHE {session}
 SET/ECHE wlcmtd=guess
 DELETE/DISPL 1
!
! Set the tolerance 0.2 higher than requested
! Set the dc +1 higher than requested
!
 tol = {tolerance} + 0.2 
 dc = {dc} + 1

begin_loop:
 iteration = {iteration} + 1
 $rm -f IDENTI.res
 IDENTI/ECHE >>IDENTI.res
 SHOW/ECHE w
 $tail -13 IDENTI.res
 $grep TOTAL IDENTI.res | gawk '{print $7}' | write/key total_last
 IF total_last .gt. total_first THEN
   SAVE/ECHE {session}
   total_first = total_last
   GOTO begin_loop
 ELSE 
   IF tol .gt. tolerance then
     tol = {tol} - 0.1
     total_first = 0
     GOTO begin_loop
   ELSE
     IF dc .gt. poly then
       dc = {dc} - 1
       total_first = 0
       GOTO begin_loop
     ELSE
       GOTO end_loop
     ENDIF
   ENDIF
 ENDIF 

end_loop:

!
! Save session in fits
! 
 WRITE/OUT "{cmd}: Identification completed after {iteration} iterations."
 OUTDISK/FITS bachesORDE.tbl bachesORDE.fit
 OUTDISK/FITS bachesLINE.tbl bachesLINE.fit

!RESET/DISPL
 WRITE/KEYW reply/c/1/1 y
 INQUIRE/KEYW reply "{cmd}: Do you want to clean temporary files [yn] (y)?"
 IF reply(:1) .eq. "y" THEN
   GOTO clean
 ELSE
   GOTO noclean
 ENDIF
!
! Clean/NoClean
!
clean:
 @@ clean.prg
!
noclean:
 RESET/DISP
 ECHO/OFF
 RETURN/EXIT