;pro dark_median
;
;PURPOSE:
;	To take a collection of dark images and get the median value of the 
;	pixels and then write them to a fits file
;
;INPUTS:
;
;OUTPUTS:
;
function Dark_Median, rootdir=rootdir, date=date,  reduceddir = reduceddir, $
		blocklength = blocklength, $ 
		ObjectName=ObjectName, NReads = NReads, ThisFilter = ThisFilter,$
		OverRide = OverRide

  If not(keyword_set(date)) then Date='07Apr30'
  If not(keyword_set(ReducedDir)) then ReducedDir=$
	   '~/Desktop/LSST/KPNO/Images/Reduced/'+Date+'/'
  If not(keyword_set(RootDir)) then RootDir = '/Volumes/SEA_DISC/'+Date+'/'
  If not(keyword_set(ObjectName)) then ObjectName = 'Dark'
  If not(keyword_set(NReads)) then NReads = 20
  If not(keyword_set(ThisFilter)) then ThisFilter = 'G'
  If not(keyword_set(OverRide)) then OverRide = 0
  ;cd,rootdir+Date+'/'
	
  ;Block Length for Offset
  if not(keyword_set(BlockLength)) then BlockLength = 0

  ;***************************************************************************************
  ;Check to see if the Median Dark already exists and only write a fits file if it doesn't
  DarkMedianFile = file_search(ReducedDir+'MedianDark*'+strtrim(NReads,2)+'Reads*.fits')

  ;If Dark does not exist, form the median and write the image
  If DarkMedianFile(0) eq '' or OverRide then begin
    print, 'Forming Median Dark Image for ' + strtrim(NReads,2) + ' Read Ramp' 

    ;Get all of the Dark Files for this Number of Reads
    Dark_Files = file_search(ObjectName+'*'+strtrim(NReads,2)+'_Reads*.fits')
    Dummy = 0 
    Valid_Indices = 0

    ;Get a few kewyords
    Fits_Read,  Dark_Files(0), NoData, FitsHeader_Dark, /header_only       
    Naxis1     = long(sxpar(FitsHeader_Dark, 'NAXIS1'))          
    Naxis2     = long(sxpar(FitsHeader_Dark, 'NAXIS2'))         
    TotalReads = sxpar(FitsHeader_Dark ,'NAXIS3')      
    ArrSize = long(Naxis1*Naxis2)
	
    ;Only use the images that are certain to have Blank as the filter
    for i = 0, n_elements(Dark_Files) - 1 do begin
        Fits_Read, Dark_Files(i), NoData, FitsHeader, /header_only
        Filter = sxpar(FitsHeader ,'FILTER')
	print, i, ' ' , Filter
        If stregex(Filter, 'Blank' , /boolean) and $
	   not stregex(Dark_Files(i), '42_44', /boolean) then $
           valid_indices = [valid_indices, i]
        endfor

	;Form an Array big enough for all of the darks of a given Naxis3
	NumDark_Files = N_elements(Valid_Indices)
	FullIms = LonArr(Naxis1, Naxis2, NumDark_Files-1)
	
	;Form Cube of Darks and then Take Median at the End 
	for DarkNum = 1, NumDark_Files - 1 do begin
	   index = Valid_Indices(DarkNum)
	   Fits_Read,  Dark_Files(index), NoData, FitsHeader, /header_only
           Filter = sxpar(FitsHeader ,'FILTER')
	   print, index, ' ', Filter, ' ', Dark_Files(index) 

           ;Get Difference of Dark Frames
	   FullIms(*,*, DarkNum - 1) = return_cube_diff(Dark_Files(index), $
		   ArrSize, Naxis1, Naxis2, 1, TotalReads)

	endfor
	
	;Write a fits file with the Median Dark
	medim = median(FullIms, dimension = 3)
        writefits,ReducedDir+'MedianDark_'+Strtrim(NReads,2)+'Reads.fits',medim
	print, 'Median Dark Image written with filename: '
	print, ReducedDir+'MedianDark_'+Strtrim(NReads,2)+'Reads.fits'

	return, 1

  EndIf else begin
        
	;The File already exists
	print, 'Median Dark Image: ' +  ReducedDir+'MedianDark_'+Strtrim(NReads,2)+'Reads.fits' + ' already exists...continuing'
	return, 2
 
  EndElse

stop

end