pro DIcal, inFn, DNFN=dnfn, RADFN=radfn, REVRADFN=revradfn, IFFN=iffn, REVIFFN=reviffn,$
	VERBOSE=verbose, SILENT=silent, FORCE=force, MODULES=modules, $
	CHOSEN_ONLY=chosen_only, BADPIXS=badpixs, DESPIKE=despike, $
	FILLGAPS=fillgaps, DENOISE=denoise, GEOM=geom, MTF=mtf, $
	DESMEARALG=desmearalg, DARKFN=darkfn, FLATFN=flatfn, $
	MAXGAPSIZE=maxgapsize, USEDARKMODEL=usedarkmodel, $
	COMPRESSFN=compressfn, ADCFN=adcfn, BADPIXSFN=badpixsfn, $
	GEOMFN=geomfn, MTFFN=mtffn, VISCONSTFN=visconstfn, LINDNFN=lindnfn, $
	SPECFN=specfn, CONSTMAPFN=constmapfn, GAINFN=gainfn, SPIKESIG=spikesig,$
	SPIKEITER=spikeIter, SPIKEBOX=spikebox, SPIKEMED=spikemed, SPIKEALG=spikealg,$
	BADPIXSINTERP=badpixsinterp, MISSINGINTERP=missinginterp, $
	HLUTFN=hlutfn, COMPRESSMETH=compressmeth, MTFALG=mtfalg, MTFPARAM=mtfparam

;;-----------------------------------------------------------------------------
;; PURPOSE:
;;	Performs the calibration on a Deep Impact FITS image.
;;	1. Temperatures and Voltages calibrated
;;	2. Decompression
;;	3. Flag saturated pixels
;;	4. Correct for uneven ADC bit weighting
;;	5. Linearize DN values (IR only)
;;	6. Remove electrical crosstalk
;;	7. Subtract the dark frame (dark found using interpolation)
;;	8. Normalize quadrant gains
;;	9. Divide by the normalized flat field
;;	10. Flag bad pixels 
;;	11. Transfer smear removed (VIS only)
;;	12. Image converted from DN units to radiance and/or I/F units
;;	13. Bad Pixels and Gaps are Reclaimed
;;      14. Cosmic Ray Removal
;;	15. Denoising
;;	16. Deconvolution
;;	17. Geometric calibration
;;
;;	By default, modules 2,3,5,6,7,9,10,11,12 and 13 are applied
;;
;;	For those modules that are applied by default and have associated
;;	keywords, if the keyword is explicitly set to 0, the module is not
;;	applied.
;;	
;;	The output FITS files consists of:
;;	Primary Data Unit - Calibrated Image
;;	1st Extension - compressed version of the flags image. To extract, 
;;		use the di_readfits_flagmap routine.
;;	2nd Extension - Wavelength map (IR only)
;;	3rd Extension - Delta Wavelength map (IR only)
;;	4th Extension - A SNR map giving an estimate of the SNR for that pixel
;;
;; CALLING SEQUENCE:
;;	DIcal, inFn, RADFN=radfn, IFFN=iffn, /VERBOSE, /SILENT, /FORCE, $
;;		MODULES=modules, /CHOSEN_ONLY, /BADPIXS, /DESPIKE, /FILLGAPS,$
;;		/DENOISE, /GEOM, /MTF, DESMEARALG=desmearalg, /RAD, /IOVERF, $
;;		DARKFN=darkfn, FLATFN=flatfn, MAXGAPSIZE=maxgapsize, $
;;		/USEDARKMODEL
;;
;; REQUIRED INPUTS:
;;	inFn - The filename of the FITS image to calibrate
;;
;; OUTPUTS:
;;
;; OPTIONAL INPUT KEYWORDS:
;;	DNFN - Filename of the calibrate file in DN
;;	RADFN - Filename of the calibrated file in radience units to output
;;	REVRADFN - Filename of the reversibly calibrated file in radiance units
;;		to output
;;	IFFN - Filename of the calibrated file in I/F units to output. RADFN
;;		and/or IFFN must be set
;;	REVIFFN - Filename of the reversibly calibrated file in I/F units to
;;		output
;;	VERBOSE - If applied then status messages will be printed
;;	SILENT - If flagged, then error message won't be displayed
;;	FORCE - If flagged then an error in a calibration step won't abort
;;		the calibration
;;	MODULES - A activeModules structure specifying which modules to use.
;;		This allows complete control of what modules are applied 
;;		(Any module whose structure element is flagged is applied)
;;	CHOSEN_ONLY - Only apply those modules that are explicity chosen are
;;		applied. This can be used to apply badPixs, despike, fillGaps,
;;		denoise, geom, mtf, desmear, dark frame and flat field
;;		operations only.
;;	BADPIXS - Applies the module to interpolate over bad pixels.  
;;	DESPIKE - Applies the module to despike the image.
;;	FILLGAPS - Causes gaps in data whose size is less than MAXGAPSIZE 
;;		(default=9999) to be interpolated over.
;;	DENOISE - Applies the module that tries to remove random gaussian noise
;;	GEOM - Applies a geometric rubber sheet correction
;;	MTF - Applies a convolution filter to compensate for MTF and scattered
;;		light
;;	DESMEARALG - Determines which desmearing algorithm to use (VIS only)
;;		0 - No desmearing occurs
;;		1 - Remove smear using the POC values (default for modes 1-6,9)
;;		2 - Remove smear using column averaging (default for modes 7&8)
;;	DARKFN -  Filename of the dark FITS file to sutract from the image.
;;		If this is not set, then a dark frame is found based on the 
;;		USEDARKMODEL keyword. Setting this keyword is equivilant to 
;;		explicity choosing the the dark subtraction operation.
;;	FLATFN - Filename of the flat FITS file to divide the image by. If
;;		this is not set, the default flat field is used. Setting this
;;		keyword is equivilant to explicity choosing the flat field
;;		compensation module
;;	MAXGAPSIZE - The maximum size of a data gap to fill if the Fill gap
;;		module is run.
;;	USEDARKMODEL - If flagged then the dark frame to subtract is found
;;		by using a pre-calculated dark model at every pixel. 
;;		Explicity setting this keyword to any value is equivilant of
;;		choosing the dark subtraction module.
;;	COMPRESSFN - filename of the decompression LUT to use
;;	ADCFN - filename of the bit weighting correction LUT to use
;;	BADPIXSFN - filename of the bad pixel map to use
;;	GEOMFN - filename of the geometric tie points to use
;;	LINDNFN - filename of the linearization polynomial to use
;;	MTFFN - filename of the convolution filter to use
;;	SPECFN - filename of the spectral map to use
;;	CONSTMAPFN - filename of the ST calibration constant map to use
;;	VISCONSTFN - filename of the calibration constant to use
;;	GAINFN - filename of the gain maps to use
;;	SPIKESIG - threshold for the despiking algorithm in distance from mean
;;		in standard deviations
;;	SPIKEITER - maximum number of iterations for the despking routine
;;	SPIKEBOX - How big a box to use to calculate statistics for despiking
;;	SPIKEMED - Should despiking use Median? Otherwise it will use mean
;;	SPIKEALG - Despiking routine algorithm
;;		1 = Sigma filtering
;;		2 = imgclean (VIS only)
;;	BADPIXSINTERP - Should the bad pixels be interpolated over.
;;	MISSINGINTERP - Should missing data be interpolated over.
;;	HLUTFN - The filename of the H lambda lookup table for IR images
;;	COMPRESSMETH - The method used to distribute DN values in decompression
;;		1 - average value
;;		2 - gaussian distribution
;;		3 - random uniform distribution
;;	MTFALG - Algorithm used for deconvolution 
;;		1 - constrained least squares
;;		2 - Richardson-Lucy with wavelet noise suppression
;;	MTFPARAM - Algorithm specific parameter for the deconvolution algorithm
;;
;; EXAMPLE:
;;      IDL> DIcal,'in.fit',radFn='rad.fit',/VERBOSE
;;
;; PROCEDURES USED (i.e. called directly!):
;;	READFITS - Reads a FITS file
;;	WRITEFITS - Writes a FITS file
;;	SXPAR - Reads a FITS parameter
;;	DIcal_VIS - Calibrates a VIS image
;;	DIcal_IR - Calibrates an IR image
;;
;; MODIFICATION HISTORY:
;;   2004-05-25  M. Desnoyer    Created
;;   2005-04-05  B. Carcich     Modified to do check_FITS,...,/SILENT
;;   2005-04-12  M. Desnoyer  	Removed call to compile routines
;;   2005-05-26  M. Desnoyer    Added new options to despiking on deconvolution
;;
;;-----------------------------------------------------------------------------

@dical_flags

doverbose = keyword_set(VERBOSE)
IF not doverbose then begin
  ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
  ;;;BTCarcich2005-04-14
  ;;; Use saved quiet state
  savequiet=!QUIET
  !QUIET = 1
ENDIF

; Do error handling 
;CATCH, error
error=0
curStep = 0

IF error NE 0 THEN BEGIN
	;; Throw an error if an error causes things to stop
	CATCH, /CANCEL
	IF keyword_set(SILENT) THEN $
		message, /NOPRINT
	message, inFn+': '+!ERROR_STATE.MSG
ENDIF

;; Make sure an output file is specified
dn = keyword_set(dnfn)
ioverf = keyword_set(iffn)
revIoverf = keyword_set(revIffn)
rad = keyword_set(radfn)
revRad = keyword_set(revRadfn)
IF (rad EQ 0) AND (ioverf EQ 0) AND (dn EQ 0) AND (revIoverf EQ 0) AND (revRad EQ 0) THEN $
	message, 'Output filename must be specified', /NONAME

;; Open the FITS file
imgIn = readFITS(inFn,hdr,/silent)

;; Make sure we got a valid file
IF size(imgIn,/n_dimensions) NE 2 THEN $
	message, 'Invalid FITS file', /NONAME $
ELSE $
	IF doverbose THEN print, 'File opened'

;; Get the flag map
catch, valid
IF valid EQ 0 THEN BEGIN
	flags = di_readfits_flagmap(inFn)
ENDIF ELSE BEGIN
	sz = size(imgIn,/dimensions)
	flags = intarr(sz[0],sz[1])
	message, /reset
ENDELSE
catch, /cancel

;; Get the detector
dtctr = sxpar(hdr,'IMGH021', COUNT=c1)

;; Make sure the items needed were in the FITS header
IF (c1 EQ 0) THEN message, 'Invalid FITS header', /NONAME

;; Determine which modules to apply
IF keyword_set(modules) EQ 0 THEN BEGIN
	modules = {activeModules}
	
	;; Modules that are defaulted as on unless the chosen_only flag is set
	NotChosen_only = keyword_set(CHOSEN_ONLY) EQ 0
	modules.tempvolt = 0
	modules.flagSat = NotChosen_only
	modules.bitWeight =  keyword_set(adcfn)
	modules.linDN = (NotChosen_only OR keyword_set(lindnfn)) $
		AND (dtctr EQ 1)
	modules.xTalk = NotChosen_only
	
	;; Modules that require more complex logic
	modules.dark = NotChosen_only OR (keyword_set(DARKFN)) OR $
		(n_elements(USEDARKMODEL) GT 0)
	modules.gain =  keyword_set(GAINFN)
	modules.flat = (NotChosen_only OR (keyword_set(FLATFN)))
	IF n_elements(DESMEARALG) GT 0 THEN $
		modules.desmear = DESMEARALG NE 0 $
	ELSE modules.desmear = NotChosen_only AND (dtctr EQ 0)
	IF n_elements(BADPIXS) GT 0 THEN $
		modules.badPixs = BADPIXS $
	ELSE modules.badPixs = NotChosen_only
	IF n_elements(DESPIKE) GT 0 THEN $
		modules.despike = DESPIKE $
	ELSE modules.despike = 0
	IF n_elements(FILLGAPS) GT 0 THEN $
		modules.fillGaps = FILLGAPS $
	ELSE modules.fillGaps = NotChosen_only
	IF n_elements(DENOISE) GT 0 THEN $
		modules.denoise = DENOISE $
	ELSE modules.denoise = 0
	IF n_elements(GEOM) GT 0 THEN $
		modules.geom = GEOM $
	ELSE modules.geom = 0
	IF n_elements(MTF) GT 0 THEN $
		modules.mtf = MTF $
	ELSE modules.mtf =0
ENDIF

;; Insert all of the calibration FITS keywords with default values
dical_insertCalKeywords, hdr

;; For adding to headers later on
stars='***********************************************************************'

;; Determine which calibration pipeline to enter based on dectector
IF dtctr EQ 0 THEN BEGIN
	IF doverbose THEN print, 'Calibrating VIS image'
	imgOut = DIcal_VIS(imgIn,hdr,modules, verbose=verbose, silent=silent,$
		force=force, maxgapsize=maxgapsize, desmearalg=desmearalg, $
		dn=dn, rad=rad, revrad=revrad, ioverf=ioverf, revioverf=revioverf,$
		flatfn=flatfn, darkfn=darkfn, usedarkmodel=usedarkmodel, $
		compressfn=compressfn, adcfn=adcfn, badpixsfn=badpixsfn, $
		geomfn=geomfn, mtffn=mtffn, visconstfn=visconstfn, flags=flags,$
		gainfn=gainfn, spikesig=spikesig, spikeiter=spikeiter, spikealg=spikealg, $
		spikebox=spikebox, spikemed=spikemed, badpixsinterp=badpixsinterp, $
		compressmeth=compressmeth, missinginterp=missinginterp, mtfalg=mtfalg,$
		mtfparam=mtfparam, snr=snr)
ENDIF ELSE IF dtctr EQ 1 THEN BEGIN
	IF doverbose THEN print, 'Calibrating IR image'
	imgOut = DIcal_IR(imgIn,hdr,modules, verbose=verbose, silent=silent,$
		force=force, maxgapsize=maxgapsize,  rad=rad, revrad=revrad,$
		ioverf=ioverf, revioverf=revioverf, dn=dn, darkfn=darkfn, $
		usedarkmodel=usedarkmodel, compressfn=compressfn, adcfn=adcfn,$
		badpixsfn=badpixsfn, mtffn=mtffn, gainfn=gainfn,$
		lindnfn=lindnfn, specfn=specfn, constmapfn=constmapfn, $
		spikesig=spikesig, spikeiter=spikeiter, spikebox=spikebox, $
		spikemed=spikemed, spikealg=spikealg, flags=flags, specMap=specMap, $
		badpixsinterp=badpixsinterp, hlutfn=hlutfn, flatfn=flatfn, $
		compressmeth=compressmeth, missinginterp=missinginterp, mtfalg=mtfalg,$
		mtfparam=mtfparam, snr=snr)
	
	;; Make the headers for the extensions
	mkhdr, specHdr0, float(specMap[*,*,0]), /image
	sxaddpar, specHdr0, 'BSCALE', 1.0, ' Data scaling factor'
	sxaddpar, specHdr0, 'BZERO', 0.0, ' Data offset value'
	sxaddpar, specHdr0, 'BUNIT', 'MICROMETER', ' Physical Units of Data'
	sxaddpar, specHdr0, 'COMMENT', stars
	specHdr1 = specHdr0
	sxaddpar, specHdr0, 'COMMENT', '** This extension is a wavelength map for the'
	sxaddpar, specHdr0, 'COMMENT', '** primary image array'
	sxaddpar, specHdr0, 'COMMENT', stars
	sxaddpar, specHdr1, 'COMMENT', '** This extension is a spectral bandwidth map for the'
	sxaddpar, specHdr1, 'COMMENT', '** primary image array'
	sxaddpar, specHdr1, 'COMMENT', stars
	sxaddpar, specHdr0, 'EXTNAME', 'WAVELENGTH', $
		' Unique name for this image extension'
	sxaddpar, specHdr1, 'EXTNAME', 'DELTA_WAVELENGTH', $
		' Unique name for this image extension'
ENDIF ELSE message, 'Invalid detector ID', /NONAME

;; Create the header for the SNR map
mkhdr, snrHdr, float(snr), /image
sxaddpar, snrHdr, 'COMMENT', stars
sxaddpar, snrHdr, 'COMMENT', '** This extension is a signal-to-noise ratio map for the'
sxaddpar, snrHdr, 'COMMENT', '** primary image array'
sxaddpar, snrHdr, 'COMMENT', stars
sxaddpar, snrHdr, 'BSCALE', 1.0, ' Data scaling factor'
sxaddpar, snrHdr, 'BZERO', 0.0, ' Data offset value'
sxaddpar, snrHdr, 'BUNIT', 'N/A', ' Data is unitless'
sxaddpar, snrHdr, 'EXTNAME', 'SNR', ' Unique name for this image extension'

;; Create a list of the files to create
fns=['','']
IF dn THEN        fns = [[fns], [dnfn,     'DN']]
IF rad THEN       fns = [[fns], [radfn,    'RAD']]
IF revrad THEN    fns = [[fns], [revradfn, 'RADREV']]
IF ioverf THEN    fns = [[fns], [iffn,     'IF']]
IF revioverf THEN fns = [[fns], [reviffn,  'IFREV']]
n = n_elements(fns[0,*]) - 1
;; Write the output FITS file(s)
FOR i=1L, n DO BEGIN

	idata = i-1L
	fn = fns[0,i]
	caltype = fns[1,i]

	lclhdr = hdr[*,idata]

        ;; Make modifications to the header
	fxaddpar, lclhdr, 'CALTYPE', caltype $
                , ' REVersible, DN=data#, RADiance, IF=I/F' $
                , after='bunit'
	tmp = getBiasCols(lclhdr, complement=imgCols)
	tmp = getPOCRows(lclhdr, complement=imgRows)
	subWindow = imgOut[imgCols[0]+1:imgCols[1]-1, imgRows[0]+1:imgRows[1]-1 ,idata]
	fxaddpar, lclhdr, 'DATAMIN', min(subWindow), ' Minimum pixel in image array'
	fxaddpar, lclhdr, 'DATAMAX', max(subWindow), ' Maximum pixel in image array'
	fxaddpar, lclhdr, 'MEDPVAL', median(subWindow), ' Median pixel in image array'
	fxaddpar, lclhdr, 'STDPVAL', stddev(subWindow), ' Standard Deviation of pixel values'
	fxaddpar, lclhdr, 'PSATVAL', -999
	fxaddpar, lclhdr, 'PSATNUM', -999

	keyNames = ['BADPXCT','MISSPXCT','DESPIKCT','INTERPCT','PSATPXCT','SATPXCT','ASATPXCT','ULTCMPCT']
	for j=0,7 do begin
		w = where((flags[*,*,idata] AND ishft(1,j)) NE 0,flgCnt)
		fxaddpar, lclhdr, keyNames[j],flgCnt
	endfor

        imgOutTmp = float(imgOut[*,*,idata])
	check_FITS, imgOutTmp, lclhdr, /UPDATE, /FITS, /SILENT
	writeFITS, fn, imgOutTmp, lclhdr
	di_writefits_flagmap, fn, flags[*,*,idata]
	IF dtctr EQ 1 THEN BEGIN
		writeFITS, fn, float(specMap[*,*,0]), specHdr0, /append
		writeFITS, fn, float(specMap[*,*,1]), specHdr1, /append
	ENDIF
	writeFITS,fn,float(snr),snrHdr, /append
	IF doverbose THEN print, 'Calibrated file written: '+fn
ENDFOR
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;BTCarcich2005-04-14
;;; Use saved quiet state
;;;!QUIET = 0
IF not doverbose THEN BEGIN
  !QUIET = savequiet
ENDIF
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

RETURN
END