FUNCTION DIcal_linearDN, in, fitsHdr, FN=fn

;;-----------------------------------------------------------------------------
;; PURPOSE:
;;	Used in the DI Calibration Pipeline.
;;	Linearizes the DN values to produce the expected output if the 
;;	instrument response were really linear with input radiance. 
;;
;; CALLING SEQUENCE:
;;	out = DIcal_linearDN(in, fitsHdr, FN=fn)
;;
;; 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:
;;	FN - Filename of the file that contains the coeficients for the 
;;		polynomial to divide by. The file should be an IDL saved file
;;		which contains an array of coeffiencent named C. c[0] is the 
;;		coeficient for the term of lowest degree.
;;
;; EXAMPLE:
;;      IDL> imgOut = DIcal_linearDN(imgIn, fitsHdr)
;;
;; PROCEDURES USED (i.e. called directly!):
;;
;; MODIFICATION HISTORY:
;;   2004-05-24  M. Desnoyer    Created
;;
;;-----------------------------------------------------------------------------

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

IF error NE 0 THEN BEGIN
	CATCH, /CANCEL
	message, 'Linearize DN Values - ' + !ERROR_STATE.MSG, /NONAME
ENDIF

;; Get the date of the file
date = sxpar(fitsHdr, 'OBSDATE', COUNT=c1)

;; Get the mode
mode = sxpar(fitsHdr, 'IMGMODE', COUNT=c2)

;; Make sure the header contains all needed info
IF (c1 EQ 0) OR (C2 EQ 0) THEN message, 'Invalid FITS Header', /NONAME

out = double(in)

;; Get the filename
IF not keyword_set(fn) THEN BEGIN
	serv = getSQLserver()
	db = 'di_calib'
	tbl = 'LIN_DN'
	sel = ['Filepath']
	cond = 'Date <= "'+date+'" AND Mode='+strtrim(mode,2)+' order by Date desc, Version desc limit 1'

	webfn = di_sqlquery(serv, db, tbl, sel, cond, dim=dim)

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

	fn1 = getLocFn(webFn, subdir='LIN_DN')
ENDIF ELSE fn1 = fn

;; Restore the coefficients
c = readfits(fn1[0], /silent)
IF n_elements(c) LT 4 THEN message, 'Invalid linearization polynomial',/noname
cdim = size(c)

;; Calculate the dividing factor for this integration time
dim = size(out,/dimensions)
IF cdim[0] EQ 3 THEN BEGIN
	;; This is a coefficient cube so we're getting the factors pixel by pixel
	factor = polymap(out,c)
	out = out / factor

	;; Fill in the FITS keywords
	fn1 = strupcase(strmid(fn1[0], strpos(fn1[0], path_sep(), /reverse_search)+1))
	FOR i=0, cdim[3]-1 DO BEGIN $
		fxaddpar, fitsHdr, 'LINAC'+strtrim(i,2), fn1
		fxaddpar, fitsHdr, 'LINBC'+strtrim(i,2), fn1
	ENDFOR
ENDIF ELSE BEGIN
	;; Factors are done on
	afactor = poly(out[0:dim[0]/2-1,*],c[*,0])
	IF cdim[0] EQ 2 THEN BEGIN
		;; Polynomial for each quadrant is different
		bfactor = poly(out[dim[0]/2:*,*],c[*,1])
		FOR i=0,cdim[1]-1 DO $
			fxaddpar, fitsHdr, 'LINBC'+strtrim(i,2), c[i,1]
	ENDIF ELSE BEGIN
		bfactor = poly(out[dim[0]/2:*,*],c[*,0])
		FOR i=0,cdim[1]-1 DO $
			fxaddpar, fitsHdr, 'LINBC'+strtrim(i,2), c[i]
	ENDELSE
	FOR i=0,cdim[1]-1 DO $
		fxaddpar, fitsHdr, 'LINAC'+strtrim(i,2), c[i]

	;; Perform the linearization
	out[0,0] = out[0:dim[0]/2-1,*] / afactor
	out[dim[0]/2,0] = out[dim[0]/2:*,*] / bfactor
ENDELSE

;; Update the FITS header
fxaddpar, fitsHdr, 'LINEARIZ', 'T'

RETURN, out

END