;^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^Lance Simms, Stanford University 2009
;pro Dark_Current_From_Slope_xRG
;
;PURPOSE:
;	To find the hot pixels as well as the volcano regions in the HxRG detectors.
;       This script will also create histograms of dark current.
;
;INPUTS:
;   FITSFILENAME: str
;       - The full path to a slopefitted dark exposure
;
;OUTPUTS:
;   FILE 1: KeyStr.HotDir+'HotMask.fits'
;     - This file contains a mask of the high dark current pixels.  If a pixel 
;       in the mask has a value of 0, its dark current falls below DARKMAX.  
;       Otherwise, the pixel will have its slope in (ADU/sec) as its value.
;
;KEYWORDS:
;       TVFLAG: int
;         0 - Don't plot
;         1 - Plot to the XWindow
;         2 - Plot to .eps  
;	UPDATEMASK: int
;	  0 - Add the hot pixels to the hotpixel mask
;	  1 - Don't add the hot pixels to the mask
; 
;**************************************************************************

pro Dark_Current_From_Slope_xRG, FitsFileName, $ 
	 TvFlag = TvFlag, $
	 UpdateMask = UpdateMask

If N_Elements(UpdateMask)    eq 0 then UpdateMask     = 0
If N_Elements(TvFlag)        eq 0 then TvFlag        = 0

;Include KeywordStruct.pro file with all the defaults for keywords
@KeywordStruct_xRG.pro
@PlotSettings_xRG.pro

;Indices with skipping of reference pixels and channel 0
NAxis1 = KeyStr.NAxis1
NAxis2 = KeyStr.NAxis2
Ref    = KeyStr.Ref

;If H4RG is used, take first channel as obsolete
If DetStr eq 'H4RG-10-007' then begin
  Ch0W   = KeyStr.Ch0W
EndIf Else Begin
  Ch0W   = Ref/2
EndElse
 
;Get a Few Kewyords
FileCheck = Return_File_Params_xRG(FitsFileName)

;Read the slopefitted file in - The zeroth frame is the slope fit
Fits_Read, FitsFileName, SlopeCube, SlopeH
SlopeIm =fltarr(KeyStr.NAxis1, KeyStr.NAxis2)
SlopeIm(*,*) = SlopeCube(*,*,0)

;Histogram the pixels if TVFLAG==1
if TvFlag gt 0 then begin
  SlopeImRealValues = where(SlopeIm ne 0, NumReal)
  if KeyStr.DetStr eq 'H1RG-022' then begin
    DCHist=histogram(SlopeIm(SlopeImRealValues),locations=loc_pix,nbins=2000,min=-1,max=20)
    CGain = 5.47 ;e-/ADU
  endif else begin
    DCHist=histogram(SlopeIm(SlopeImRealValues),locations=loc_pix,nbins=2000,min=-1,max=20)
  endelse

  MaxHist = max(DCHist)
  MaxPlot = MaxHist*1.1
  XTitle = 'ADU/s'
  YTitle = 'Number'
  XMinADU   = 0   &  XMinElec = XMinADU*CGain
  XMaxADU   = .5  &  XMaxElec = XMaxADU*CGain
  MeanDC    = mean(SlopeIm)*CGain
  MedDC     = median(SlopeIm)*CGain
  ModeDC    = loc_pix(where(DCHist eq max(DCHist)))*CGain
  !P.FONT=0

  If TvFlag eq 2 then begin
    device, filename = KeyStr.HotDir+KeyStr.RawObjectKey+'_Histogram.eps', $
            bits_per_pixel=8, /encapsul, ysize=5.75, xsize=8.0,/inches
    device, set_font='Times-Roman'
  EndIf
  plot,loc_pix,DCHist,psym=10,xrange=[XMinADU,XMaxADU], yrange=[0,MaxHist], $
       background=fsc_color('white'), color =fsc_color('black'), $
       XTitle=XTitle, YTitle=YTitle, position=[0.1,0.1,0.9,0.9]
  axis,xaxis=1, xrange=[XMinElec,XMaxElec], /nodata, $
       xtitle= 'e!U-!N/s', color = fsc_color('black'), $
       xtickv=findgen(12)*0.2*(XMaxElec-XMinElec), xticks=10, xminor=5

  xyouts, (!X.CRANGE(1)-!X.CRANGE(0))*0.60+!X.CRANGE(0),  $
          (!Y.CRANGE(1)-!Y.CRANGE(0))*0.80+!Y.CRANGE(0), $
          'Detector        :  '+ KeyStr.DetStr, color=fsc_color('black')
  xyouts, (!X.CRANGE(1)-!X.CRANGE(0))*0.60+!X.CRANGE(0),  $
          (!Y.CRANGE(1)-!Y.CRANGE(0))*0.75+!Y.CRANGE(0), $
          'Temperature : 160 K', color=fsc_color('black')
  xyouts, (!X.CRANGE(1)-!X.CRANGE(0))*0.60+!X.CRANGE(0),  $
          (!Y.CRANGE(1)-!Y.CRANGE(0))*0.70+!Y.CRANGE(0), $
          'Mean            :  '+ strtrim(string(MeanDC,format='(F10.4)'),2) + ' e!U-!N/s', $
           color=fsc_color('black')
  xyouts, (!X.CRANGE(1)-!X.CRANGE(0))*0.60+!X.CRANGE(0),  $
          (!Y.CRANGE(1)-!Y.CRANGE(0))*0.65+!Y.CRANGE(0), $
          'Median         :  '+ strtrim(string(MedDC,format='(F10.4)'),2) + ' e!U-!N/s', $
           color=fsc_color('black')
  xyouts, (!X.CRANGE(1)-!X.CRANGE(0))*0.60+!X.CRANGE(0),  $
          (!Y.CRANGE(1)-!Y.CRANGE(0))*0.60+!Y.CRANGE(0), $
          'Mode            :  '+ strtrim(string(ModeDC,format='(F10.4)'),2) +' e!U-!N/s', $
           color=fsc_color('black')

  if TvFlag eq 2 then device, /close

endif

;Dark Threshold
MaxDark = 1.0 ;In ADU/second

;Make a hot mask
HotPix  = where(SlopeIm gt MaxDark, NumHot)
HotMask = bytarr(KeyStr.NAxis1, KeyStr.NAxis2)
HotMask(HotPix)=1

;Check to see if a Hot Mask already exists.  If it does not, make one
HotMaskExists = File_Search(KeyStr.HotDir+'HotMask.fits')

If HotMaskExists eq '' then begin
  ;Prepare a Header
  MkHdr, HotHeader, HotMask
  FxAddPar, HotHeader, 'NUMHOT' , NumHot, 'Number of hot pixels on sca'
  FxAddPar, HotHeader, 'NUMEXP' , 1, 'Number of Exposures Used'
  FxAddPar, HotHeader, 'FILE1'  , KeyStr.RawObjectKey, 'Filename'

  ;Write the mask to the appropriate directory
  Fits_Write, KeyStr.HotDir+'HotMask.fits', HotMask, HotHeader
EndIf Else Begin
  Fits_Read, KeyStr.HotDir+'HotMask.fits', HotMaskOld, HotHeader
  NumHotOld = Long(SxPar(HotHeader, 'NUMHOT'))
  NumExp    = Long(SxPar(HotHeader, 'NUMEXP'))
  HotMask   = HotMask+HotMaskOld
  FxAddPar, HotHeader, 'NUMEXP', NumExp+1, 'Number of Exposures Used'  
  FxAddPar, HotHeader, 'FILE'+Strtrim(NumExp+1), KeyStr.RawObjectKey, 'Filename'
  Fits_Write, KeyStr.HotDir+'HotMask.fits', HotMaskOld, HotHeader

EndElse

;There is a wierd kink in the histogram
MinKink  = 1.
MaxKink  = 3.
KinkDark = where(SlopeIm gt MinKink and SlopeIm lt MaxKink) 
KinkMask = bytarr(KeyStr.NAxis1, KeyStr.NAxis2)
KinkMask(KinkDark)=1
Fits_Write, KeyStr.HotDir+'KinkInDarkDistributionMask.fits', KinkMask

stop

end