;+
; PROJECT:
;	GOES-12 Solar X-ray Imager (SXI)
;
; NAME:
;	SXIG12_PREP
;
; PURPOSE:
;	Apply corrections and calibrations to raw images.
;
; EXPLANATION:
;	TBS
;
; CALLING SEQUENCE:
;   sxig12_prep,input_path,input_filename,output_path,output_filename,image_type,  $
;
; INPUTS:
;	input_filename	=	name of level-0 FITS file to read (string)
;	input_path		=	name of input directory (string)
;
; OUTPUTS:
;	output_filename	=	name of level-1 FITS file to write
;	output_path		=	name of output directory (string)
;
; OPTIONAL INPUT KEYWORDS:
;	N/A
;
; EXTERNAL CALLS:
;	FAKE_IMAGES,filename,fake,replace_file
;   HEAD2STC,header (function)
;   IMAGE_STATISTICS,image1 (function)
;   MAKE_GIFS,image,index,image_type,image_level,gif,stamp
;	READFITS,input_filename,header,exten_no=exten_no (function)
;   SXI_021_HDR,input_header,master_header,statistics,image_type
;   SXIM_FILL,image,counter1,counter2 (function)
;	SXIM_DESPIKE,image0	(function)
;   WRITEFITS,output_filename,image,header
;
; CALLED BY:
;	SXIM_BATCH or operational Python code
;
; METHOD:
;	TBS
;
; HISTORY:
;   003.000 Modified by SMH 8/6/2001
;           Added IDL error reset if backgrounds aren't found by
;           IDL code. This special handling simplifies python driver code.
;   003.001 Modified by VJP 8/13/2001
;           added background correction logic: sxim_bkgd
;           changed linear interpolation in time between short and long bkgd
;           put in new flat (with correct y-coord)
;           eliminated registration step for non-sun image PNGs
;   003.002 Modified by SMH 8/14/2001
;           Fixed typo in code to get long dark file
;           Updated version code written to L-1 FITS header
;   003.003 Modified by SMH 8/24/2001
;           Applies despike threshold of 5 and updates relevant keywords.
;   003.004 Modified by SMH 8/27/2001
;           Added handling for negative or zero MCP voltages.
;   003.005 Modified by SMH 8/28/2001
;           Swapped XCEN and YCEN and changed sign for CROTA correct
;           to correct image registration in level-1 PNG
;   003.006 Modified by SMH 9/04/2001
;           Added logic to do sqrt scaling (gamma=0.5) for linearly
;           scaled non-radshld images
;   003.007 Modified by SMH 9/05/2001
;           Added count for saturated pixels (SAT_PIX FITS keyword)
;   003.008 Modified by SMH 9/19/2001
;           make_browse now does constant scaling
;   003.009 Modified by SMH 9/19/2001
;           Changed sign on CROTA derotation code below
;   003.010 Modified by SMH 9/21/2001
;           Moved image statistics out of non-dark frame logic-block
;           so stats will be computed for dark frames.
;   003.011 Modified by SMH 9/21/2001
;           Changed gamma levels from 0.5 to 1.0 for linear images
;           Goes along with changes to make_browse.pro
;   003.012 Modified by SMH 10/3/2001
;           Fixed problem of loading dark file names into L1 FITS
;           header.
;   003.100	Modified by SMH 10/29/2001
;           Many updates: broke out background processing into
;           separate routine, cleaned up extraneous code, improved
;           comments.
;   003.101 Modified by SMH 10/30/2001
;           Replaced flat001.bin with flat002_pc.bin. Earlier flat
;           was improperly oriented.
;   003.102 Modified by SMH 11/6/2001
;           Modified sxim_interp.pro to interpolate over zeroes in
;           addition to -1's.  Our main problem is proton events,
;           not data drops.
;   003.103 Modified by SMH 11/19/2001
;           Changed make_browse level-1 scaling displayed on PNGs
;   003.200 Modified by SMH 02/01/2002
;           Added FITS extensions for pixels affected by saturation,
;           leaky reference pixels, TLM dropouts, median filtering.
;   003.201 Modified by SMH 02/06/2002
;           Changed name from SXIM_PREP to SXIG12_PREP to conform to
;           new naming standard.
;           Included usage of SXIG12_READ to handle file errors, e.g.,
;           files with no image array etc.
;           Corrected an error with extensions header defaults for
;           NAXIS1
;           Extracted PNG product generation into SXIG12_browsegen
;   003.202 Modified by SMH 05/16/2002
;           Removed history prior to version 3.000
;           Removed call to SXIM_BKGD for dark frame processing
;           Set BSCALE to 1 to represent DN/S as output units.  This was
;           done because of the high uncertainties in absolute calibration
;   003.203 Modified by SMH 05/21/2001
;           Changed output_filename so that it now includes the '.FTS' extension
;   003.204 Modified by SMH 05/29/2002
;           Changed so that an error in Image_Type results in 'NA' being put
;           in the level-1 FITS header
;   003.205 Modified by SMH 6/11/2002
;           Updated sxig12_get_bg call, set default scaleoption on sxig12_read to
;           'LIN' so that LOG images are converted to a linear scale
;           Changed code so that LOG images have average of model background
;           subtracted.  This helps beat down the digitization noise on the
;           backgrounds since they were taken with the linear amplifier.
;   003.206 Modified by SMH 6/11/2002
;           Changed sxig12_make_browse to include backup for level-1 image
;           code determination
;   003.207 Modified by SMH 6/12/2002
;           No change to this file, but changed date-time format in browse
;           images (sxi_make_browse) to be more compatible with SQL and
;           CCSDS formats.
;   003.208 Modified by SMH 6/13/2002
;           No change to this file, but changed sxig12_read so that old FITS keywords,
;           e.g., CROTA, are updated in both the header parameter and in the
;           index keyword paramter
;   003.209 Modified by SMH 6/23/2002
;           No change to this file, but fixed coordinate and sign errors with
;           HASS corrections in sxig12_read.pro
;   003.210 Modified by SMH 6/23/2002
;           No change to this file, but in sxig12_browsegen:
;           Added missing=0 keyword to ROT IDL procedure to avoid
;           extrapolation beyond CCD edge in PNG files
;   003.211 Modified by SMH 6/25/2002
;           Fixed error in computation of FIXPIX array that led to
;           failures, particularly for windowed images and those
;           with RFI problems.
;   003.212 Modified by SMH 6/27/2002
;           Fixed index1.FILENAME so that the *.FTS extension was
;           no longer included.
;   003.213 Modified by SMH 1/17/2003
;           Removed return for DARK FRAMES.  We want to process them
;           to have statistics in the meta-data. Updated Level-1 header
;           template, used to not fill 80 columns on END statement.
;           Added code to include DN values for saturated pixels taken
;           from nearest short exposure image in SAT_PIX extension.
;   003.214 Modified by SMH 1/22/2003
;           Ensured that dark frames are converted to floating point.
;           Added 'catch' error handler to replace ON_IOERROR handler
;           A little clean up on IO (e.g. /get_lun etc.)
;   003.215 Modified by SMH 1/23/2003
;           Made correction so that satpix is not undefined when
;           last_file = 'NULL'
;           --commented out change because of production string problem
;   003.300 Modified by SMH 1/23/2003
;           Fixed bugs with sat_pix extension.
;   003.310 Modified by SMH 4/22/2003
;           Upgrades external to this routine:
;             Better background model based on post PLT data
;             Background subtraction for 'OPEN' light leak
;             Improved error handling in SXI FITS read programs
;             Hardcoded BSCALE=1.0
;
;-
pro sxig12_prep,input_path,input_filename,output_path,output_filename,image_type

message,/reset

  CATCH, Error_status

; This statement begins the error handler:
  IF (Error_status NE 0) THEN BEGIN
    errstr = strtrim(string(Error_status),2) + ' '
    result_message = 'ERROR: ' + errstr + !ERROR_STATE.MSG
    close,/all

  ; detailed error logging
    openw, /append, u, 'sxiprep_logfile.txt',/get_lun
    printf,u,'sxiprep trapped an error'
    printf,u,'FILENAME: ',input_filename
    printf,u,'System Time: ',systime(0)
    printf,u,'error_state:'
    printf,u,'name = ', !error_state.name
    printf,u,'block = ', !error_state.block
    printf,u,'code = ', !error_state.code
    printf,u,'sys_code = ', !error_state.sys_code
    printf,u,'msg = ', !error_state.msg
    printf,u,'sys_msg = ', !error_state.sys_msg
    printf,u,'msg_prefix = ', !error_state.msg_prefix
    free_lun,u

    return
  ENDIF

;;;;;;;;;;;;;;;;;;;;;
;                   ;
; Read config. file ;
;                   ;
;;;;;;;;;;;;;;;;;;;;;

dark_path=''
fake=0
new_fake=0
replace_file=0
last_file=''

openr,u,'config.asc',/get_lun
readf,u,dark_path
readf,u,fake
readf,u,new_fake
readf,u,replace_file
readf,u,last_file
free_lun,u

;***************************************************************************
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;                                  ;
; Image generator - For Test Only! ;
;                                  ;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
if fake eq 1 then imagef=fake_images(input_path,input_filename,image_type,new_fake,replace_file)
;***************************************************************************
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;                           ;
; Read in level-0 FITS file ;
;                           ;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

filename=input_path+input_filename
image0=sxig12_read(filename,hdr0,pinfo,index=index0,fn_struct=in_fn_struct,scaleoption='LIN')

;*****Dangerous to return both the index and the header!!!!!

if !error_state.msg ne '' then begin
  output_filename = 'ERROR: '+!error_state.msg
  close,/all
  return
endif

;*********
  if fake eq 1 then image0=imagef  ;*******test only code
;*********

; Make file names and structures          ;
out_fn_struct=in_fn_struct
out_fn_struct.ProductCode = 'B'
out_fn_struct.VersionCode = 'B'
output_filename=sxi_make_name(out_fn_struct)
;print,input_filename
;print,in_fn_struct

; Get Level-1 primary header template and load basic parameters.
;  NOTE - NEVER EDIT HEADER TEMPLATE IN IDL - IT TRUNCATES TRAILING SPACES.
new_hdr=strarr(79)
level1_hdr_filename='level1_hdr.txt'

openr,u,level1_hdr_filename,/get_lun
readf,u,new_hdr
free_lun,u
hdr1 = sxig12_LoadLevel1Header(hdr0,new_hdr)
index1=fitshead2struct(hdr1)
index1.level0id=index0.FILENAME

index1.FILENAME=sxi_make_name(out_fn_struct,/kernel)

if strmid(image_type,0,5) ne 'ERROR' then begin
  index1.IMG_CODE=image_type
endif else begin
  index1.IMG_CODE='NA'
endelse

index1.VERLEVL1='003.310'

; Get Level-1 extension header template and load basic parameters
ext_hdr=strarr(12)
extension_hdr_filename='Extension_hdr.txt'
openr,u,extension_hdr_filename,/get_lun
readf,u,ext_hdr
free_lun,u
ext1_hdr1=ext_hdr
ext1_hdr1=sxi_sxmodpar(ext1_hdr1,'EXTNAME','SAT_PIX')
ext2_hdr1=ext_hdr
ext2_hdr1=sxi_sxmodpar(ext2_hdr1,'EXTNAME','ZERO_PIX')
ext3_hdr1=ext_hdr
ext3_hdr1=sxi_sxmodpar(ext3_hdr1,'EXTNAME','MISS_PIX')
ext4_hdr1=ext_hdr
ext4_hdr1=sxi_sxmodpar(ext4_hdr1,'EXTNAME','FIX_PIX')
ext5_hdr1=ext_hdr
ext5_hdr1=sxi_sxmodpar(ext5_hdr1,'EXTNAME','MED_PIX')

;;;;;;;;;;;;;;;;;;;;;;;;;;;
;                         ;
; Basic processing for    ;
; all images              ;
;                         ;
;;;;;;;;;;;;;;;;;;;;;;;;;;;

;Count and identify saturated pixels nearest neighbors
index1.sat_pix=sxig12_satpix(image0,ext1_hdr1,satarrtmp)

if size(satarrtmp,/type) ne 0 then satarr=satarrtmp

;Get short exposure pixels to append to list

if index1.sat_pix gt 0 and $
   last_file ne 'NULL' and $
   index1.EXPTIME ge 2.0 and $
   index1.WAVELNTH ne 'RDSH' then begin

  tmp=anytim2jd(index1.date_obs)
  jdsat=double(tmp.int)+tmp.frac

  prev_file=''

  openr,u,last_file,/get_lun
  readf,u,prev_file
  free_lun,u

  unsatfile=input_path+prev_file
  unsat=sxig12_read_one(unsatfile,unsathdr,pinfo,index=unsatindex)

  if size(unsat,/n_dimensions) eq 0 then goto, EXIT_SATPIX

  tmp=anytim2jd(unsatindex.date_obs)
  jdunsat=double(tmp.int)+tmp.frac

; Verify compatibility of images
  if abs(jdsat-jdunsat) le 90.0/86400.0 and $
     unsatindex.WAVELNTH eq index1.WAVELNTH and $
     unsatindex.MCP_GAIN eq index1.MCP_GAIN then begin

    unsat_reg=rot(unsat,0.0,1.0,(255.5+index1.XCEN/5.0)-unsatindex.XCEN/5.0, $
         (255.5+index1.YCEN/5.0)-unsatindex.YCEN/5.0,/interp)

    satpix=long(satarrtmp(0,*))+long(satarrtmp(1,*))*512
    unsatpixDNs=unsat_reg(satpix)

    ext1_hdr1=sxi_sxmodpar(ext1_hdr1,'NAXIS1',3)
    ext1_hdr1=sxi_sxmodpar(ext1_hdr1,'BITPIX',-32)
    ext1_hdr1=sxi_sxmodpar(ext1_hdr1,'SRC_FILE',unsatindex.filename,after='EXTVER')
    satarr=fltarr(3,index1.sat_pix)
    satarr(0,*)=satarrtmp(0,*)
    satarr(1,*)=satarrtmp(1,*)
    satarr(2,*)=unsatpixDNs

  endif else begin
    satarr=satarrtmp
  endelse

endif

EXIT_SATPIX: print,''
;EXIT_SATPIX: begin
;  close,/all
;  if size(satarrtmp,/n_dimensions) eq 2 then satarr=satarrtmp
;end


; Count and identify zero (0) pixels
index1.zero_pix=sxig12_zeropix(image0,ext2_hdr1,zeroarr)

; Count, indentify, and FIX missing (-1) pixels.
index1.fix_pix=sxi_interp(400,image0,image1,miss_pix)
index1.miss_pix=miss_pix

if index1.miss_pix gt 0 then begin
  ext3_hdr1=sxi_sxmodpar(ext3_hdr1,'NAXIS2',index1.miss_pix)
  dummy=where(image0 eq -1, temp_count)
  misspix=intarr(2,temp_count)
  misspix(0,*)=dummy mod 512
  misspix(1,*)=fix(dummy/512)

;  dummy_list=float(image0 eq -1)*image1     ;select all original missing pixels
;  dummy=where(dummy_list ne -1,temp_count)  ;select fixed pixels

  dummy=where(image0 ne image1,temp_count)   ;select fixed pixels

  if temp_count gt 0 then begin
    ext4_hdr1=sxi_sxmodpar(ext4_hdr1,'NAXIS2',index1.fix_pix)
    fixpix=intarr(2,temp_count)
    fixpix(0,*)=dummy mod 512
    fixpix(1,*)=fix(dummy/512)
  endif else begin
    ext4_hdr1=sxi_sxmodpar(ext4_hdr1,'NAXIS',0)
    sxdelpar,ext4_hdr1,'NAXIS1'
    sxdelpar,ext4_hdr1,'NAXIS2'
  endelse

endif else begin
  ext3_hdr1=sxi_sxmodpar(ext3_hdr1,'NAXIS',0)
  sxdelpar,ext3_hdr1,'NAXIS1'
  sxdelpar,ext3_hdr1,'NAXIS2'
  ext4_hdr1=sxi_sxmodpar(ext4_hdr1,'NAXIS',0)
  sxdelpar,ext4_hdr1,'NAXIS1'
  sxdelpar,ext4_hdr1,'NAXIS2'
endelse

; Apply Median filter for spike/track removal
index1.med_thrs=-1
image1=sxi_despike(image1,fixed_pixels,threshold=index1.med_thrs,spike_list=med_pixels)
index1.med_pix=fixed_pixels
ext5_hdr1=sxi_sxmodpar(ext5_hdr1,'NAXIS2',index1.med_pix)

if index1.med_pix eq 0 then begin
  ext5_hdr1=sxi_sxmodpar(ext5_hdr1,'NAXIS',0)
  sxdelpar,ext5_hdr1,'NAXIS1'
  sxdelpar,ext5_hdr1,'NAXIS2'
endif

;;;;;;;;;;;;;;;;;;;;;;;
;                     ;
; START: Logic Block  ;
; for Non-Dark Frames ;
;                     ;
;;;;;;;;;;;;;;;;;;;;;;;
browse_msg=''
;print,index0.object,index0.wavelnth
if (index0.object ne 'DARK' and index0.wavelnth ne 'RDSH') then begin
  browse_msg=sxig12_get_bg(dark_path,index1,background,shrtdark,longdark)
  index1.SHRTDARK=shrtdark
  index1.LONGDARK=longdark

  if strmid(browse_msg,0,5) eq 'ERROR' then begin
    if size(background,/n_dimensions) ne 2 then background=fltarr(512,512)
    ;WHAT IS ERROR HANDLING STRATEGY AND MESSAGING HERE!!!
  endif

  ; Subtract background
  if index0.LIN_LOG eq 'LIN' then image1=image1-background
  if index0.LIN_LOG eq 'LOG' then image1=image1-total(background)/n_elements(background)

  ; Divide by flat field
  pth=''
  flatname='flat002_pc.bin'
  openr,u,pth+flatname,/get_lun
  tmp=assoc(u,fltarr(512,512))
  flat=tmp(0)
  free_lun,u
  image1=image1/flat

  ; Remove scattering background

  ;;;;;;;;;;;;;;;;;;;;;;;
  ;                     ;
  ; Compute calibration ;                    480.2 is the average value of the factor
  ; keyword values      ;                    in Pizzo 20 Apr 1999 dectector report
  ;                     ;                    (considering 8% and 12% QE values).
  ;;;;;;;;;;;;;;;;;;;;;;;
  image1=image1/float(index1.exptime)       ;normalize for exposure (DN/sec)
  ;     Need to deal with log scaling also!

  index1.BZERO=0.0
  index1.BSCALE=1.0
;STUB for actual calibration, assumes TSYNC packet available
;  if index0.MCP1K_V ge 400.0 and index0.MCP1K_V le 1000.0 then begin
;    MCP_VOLTS = index0.MCP1K_V
;  endif else begin
;    MCP_VOLTS = index0.MCP_GAIN
;    browse_msg='ERROR: MCP GAIN'
;  endelse

;  if MCP_VOLTS ge 400.0 and MCP_VOLTS le 1000.0 then begin
;    index1.BSCALE=float(1.0/(480.2*(float(MCP_VOLTS)/1000.0)^14.2))   ;detected photons/DN
;  endif else begin
;    index1.BSCALE = 1.0
;    browse_msg='ERROR: MCP GAIN'
;  endelse

endif

;;;;;;;;;;;;;;;;;;;;;;;
;                     ;
; END: Logic Block  ;
; for Non-Dark Frames ;
;                     ;
;;;;;;;;;;;;;;;;;;;;;;;

;Compute image statistics
image_stats,image1,index1

;;;;;;;;;;;;;;;;;;;;;;;;;
;                       ;
; Final header keywords ;
; and write FITS file   ;
;                       ;
;;;;;;;;;;;;;;;;;;;;;;;;;

index1.DATE=long_time_str(BIN_DATE(SYSTIME(0)))
;hdr1=sxi_struct2fitshead(index1)  ;don't like this routine

Lev0_booleans=['EXTEND','LIN_DSBL','SAD_DSBL']
nt=n_tags(index1)
keyword=tag_names(index1)

for i=0,nt-3 do begin
  tmp=where(Lev0_booleans eq keyword(i),cnt)
  if cnt eq 0 then begin
    if i ne 24 then hdr1=sxi_sxmodpar(hdr1,keyword(i),index1.(i)) ; exception for DATE-OBS
  endif else begin
    case index1.(i) of
      0:  bool='F'
      1:  bool='T'
    endcase
    hdr1=sxi_sxmodpar(hdr1,keyword(i),bool)
  endelse
endfor

outfile=sxi_make_name(out_fn_struct)

writefits,output_path+outfile,float(image1),hdr1

writefits,output_path+outfile,satarr,ext1_hdr1,/append
writefits,output_path+outfile,zeroarr,ext2_hdr1,/append
writefits,output_path+outfile,misspix,ext3_hdr1,/append
writefits,output_path+outfile,fixpix,ext4_hdr1,/append
writefits,output_path+outfile,med_pixels,ext5_hdr1,/append

;Save most recent short exposure filename for 'compositing', i.e.
; including pixel values in SATPIX extension of NEXT image
; if it (the next image) includes saturated pixels.
;print,last_file
if last_file ne 'NULL' then begin
  if index1.WAVELNTH ne 'RDSH' and index1.EXPTIME lt 0.5 then begin
    openw,u,last_file,/get_lun
    printf,u,outfile
    free_lun,u
  endif
endif

close,/all  ; This is a kludge to fix an unclosed file unit somewhere in the code

end