FORWARD_FUNCTION DIcal_RADtoIF, getHlambda
;;-----------------------------------------------------------------------------
;; PURPOSE:
;;	Used in the DI Calibration Pipeline.
;;	Converts a IR image from radience units to I/F units
;;
;; CALLING SEQUENCE:
;;	x = DIcal_RADtoIF(in, fitsHdr, SPECFN=specfn, GEOMFN=geomfn)
;;
;; REQUIRED INPUTS:
;;	in - The image to which this pipeline element is going to be applied
;;	fitsHdr - The FITS header for the image
;;
;; OUTPUTS:
;;	RETURN - the image after going through this calibration step
;;
;; OPTIONAL INPUT KEYWORDS:
;;	SPECFN - Filename of the spectral maps
;;	GEOMFN - Filename of the data for geometric calibration
;;	HLUTFN - Filename of the lookup table for H lambda
;;
;; EXAMPLE:
;;      IDL> imgIF = DIcal_RADtoIF(imgRad, fitsHdr)
;;
;; PROCEDURES USED (i.e. called directly!):
;;	SXPAR - retrieve parameters from the FITS header
;;	FXADDPAR - change items in the header
;;	getHLambda - Determines Hlambda
;;
;; MODIFICATION HISTORY:
;;   2004-07-22  M. Desnoyer    Created
;;   2005-04-12  M. Desnoyer    Reordered functions
;;
;;-----------------------------------------------------------------------------

;;-----------------------------------------------------------------------------
;; PURPOSE:
;;	Determines the solar spectral irradiance (W/[m^2 um]) at 1 AU for every
;;	pixel in the image. 
;;
;; CALLING SEQUENCE:
;;	H = getHlambda(fitsHdr, SPECFN=specfn, GEOMFN=geomfn)
;;
;; REQUIRED INPUTS:
;;	fitsHdr - The FITS header of the image you want Hlambda for
;;
;; OUTPUTS:
;;	RETURN - an image of the value of Hlambda at every pixel
;;
;; OPTIONAL INPUT KEYWORDS:
;;	SPECFN - Filename of the spectral maps
;;	GEOMFN - Filename of the data for geometric calibration
;;	LUTFN - Filename of the Hlambda LUT
;;
;; EXAMPLE:
;;      IDL> H = getHlambda(hdr)
;;
;; PROCEDURES USED (i.e. called directly!):
;;	SXPAR - Get items from the FITS header
;;	getIRSpectralMap - Gets the spectral map for the given temperature
;;
;; MODIFICATION HISTORY:
;;   2004-07-22  M. Desnoyer    Created
;;
;;-----------------------------------------------------------------------------
FUNCTION getHlambda, hdr, SPECFN=specfn, GEOMFN=geomfn, LUTFN=lutfn

;; Error handling
ON_ERROR, 2

;; Get the date
date = SXPAR(hdr, 'OBSDATE', COUNT=c1)

;; Detmine if the image went through geometric calibration
geom = SXPAR(hdr, 'GEOMCAL', count = c2)

IF (c1 EQ 0) OR (c2 EQ 0) THEN $
	message, 'Invalid FITS header', /NONAME

;; Get the spectral map in um
lambda = (getIRSpectralMap(hdr, DOGEOM=geom, GEOMFN=geomfn, $
	FN=specfn))[*,*,0]

;; Get the H lambda lookup table
IF NOT keyword_set(lutfn) THEN BEGIN
	db = 'di_calib'
	tbl = 'SOLAR_SPECTRUM'
	sel = ['Filepath']
	cond = 'Date <="'+date+'" order by Date desc, Version desc limit 1'

	fn1 = di_sqlquery(getSQLServer(),db,tbl,sel,cond, DIM=dim)

	IF min(dim) EQ 0 THEN message, 'Could not find Hlambda LUT', /noname

	fn1 = getLocfn(fn1, subdir='SOLAR_SPECTRUM')
ENDIF ELSE fn1 = lutfn

lut = readFITS(fn1[0], lhdr, /silent)
IF !ERROR_STATE.CODE LT 0 THEN message, !ERROR_STATE.MSG, /noname

;; The LUT is two rows. The first row is the wavelength (um) and the second row
;; is the irradiance values. So, interpolate for every pixel using the wavelength
;; map
out = interpol(lut[*,1],lut[*,0],lambda,/spline)

;; Update the FITS header
fxaddpar, hdr, 'HLUTFN', $
	strupcase(strmid(fn1[0], strpos(fn1[0], path_sep(), /reverse_search)+1))

RETURN, out

END

;;-----------------------------------------------------------------------------
;; Main function
;;-----------------------------------------------------------------------------
FUNCTION DIcal_RADtoIF, in, fitsHdr, SPECFN=specfn, GEOMFN=geomfn, HLUTFN=hlutfn

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

IF error NE 0 THEN BEGIN
	CATCH, /CANCEL
	message, 'Convert from Rad to I/F - ' + !ERROR_STATE.MSG, /NONAME
ENDIF

;; Get the time of observation
time = SXPAR(fitsHdr, 'FNISODOY', COUNT=c1)

;; Get the current units used
units = SXPAR(fitsHdr, 'BUNIT', COUNT=c2)

;; Make sure that all of the items were in the header
IF (c1 EQ 0) OR (c2 EQ 0) THEN $
	message, 'Invalid FITS header', /NONAME

;; Check the units of the file to see if this can be applied
IF strtrim(units,2) EQ 'I/F' THEN RETURN, in ;; This has already been applied
IF strtrim(units,2) NE 'W/[m^2 sr um]' THEN $
	message, 'Image must be in radiance units to convert to I/F'

D = getSunTgtDist(fitsHdr)
H = getHlambda(fitsHdr, GEOMFN=geomfn, SPECFN=specfn, LUTFN=hlutfn)

;; Change the data unit in the header
fxaddpar, fitsHdr, 'BUNIT', 'I/F'
fxaddpar, fitsHdr, 'IOFCALD', D

RETURN, double(in)*(!dPI*D^2)/H
END