;pro examine_leaky_pixels
;
;PURPOSE:
;	To produce some plots of the bad types of pixels in dark and flat 
;	frames.
;
;INPUTS:
;	x - The center x-coordinate of a 3x3 region
;       y - The center y-coordinate of a 3x3 region
;KEYWORDS:
;	FROMMASK
;	    If this keyword is set, the Leaky and Hot Pixel masks will be 
;	    read in and the coordinates for the corresponding pixels will 
;	    be held in an array to be used.
;	AVGS - 
;	   0 - Plot individual pixels as we go through the for loop
;	   1 - Plot histograms at the end of the file
;***********************************************************************
pro Examine_Leaky_And_Hot_Pixels, x, y, ThisFilter=ThisFilter, FromMask = FromMask,$
	TvFlag = TvFlag, Avgs = Avgs, TestSlopeFit = TestSlopeFit

;Get the keywords from the structure file
@KeywordStruct.pro
cd, RootDir+Date

;Get plotting configuration
@PlotSettings.pro

;Search the given directory for the darks and the flats
DarkFiles = File_Search('Dark*'+strtrim(NReads,2)+'_Reads*.fits')
FlatFiles = File_Search('Flat*2*Reads*.fits')
DarkFile = DarkFiles(1)
DarkFileKey = strmid(DarkFile, stregex(DarkFile,'Reads_')+6)
print,DarkFile
FlatFile = FlatFiles(3)
FlatFileKey = strmid(FlatFile, stregex(FlatFile,'Reads_')+6)
print,FlatFile

If TvFlag eq 2 then begin 
   If Avgs then begin
     PSFileName = LeakyDir + DarkFileKey+'LeakyAndHotPixelsHistogram.ps'
   EndIf Else begin
     PSFileName = LeakyDir + DakrFileKey+'LeakyAndHotPixelsEx.ps'
   EndElse
     device, filename = PSFilename
EndIf

;If FromMask = 1, take the mask files and get all of the indices
If FromMask then begin
 ;Get the corresponding mask files 
 HotMaskFile  = File_Search(HotDir, 'HotPixelsMask*_'+ThisFilter+'_*.fits')
 LeakyMaskFile = File_Search(LeakyDir, 'LeakyPixelMask*_'+ThisFilter+'_*.fits')

 ;Read them in
 Fits_Read, HotMaskFile, HotMask, HotHeader
 Fits_Read, LeakyMaskFile, LeakyMask, LeakyHeader

 print, max(HotMask)
 print, max(LeakyMask)

 HotPixInd   = where(HotMask eq Max(HotMask))
 HotPixX     = HotPixInd mod Naxis1
 HotPixY     = HotPixInd/Naxis2

 LeakyPixInd = where(LeakyMask eq Max(LeakyMask))
 LeakyPixX   = LeakyPixInd mod Naxis1
 LeakyPixY   = LeakyPixInd / Naxis2

EndIf

;Form arrays for some averages
NLeaky = N_Elements(LeakyPixInd)
NNFracLeakyD = fltarr(NLeaky)
NNFracLeakyF = fltarr(NLeaky)

NHot = N_Elements(HotPixInd)
NNFracHotD = fltarr(NHot)
NNFracHotF = fltarr(NHot)

;Cycle through hot pixels first and then leaky pixels second
MaxPixInd = NHot
For PixType = 0, 1 do begin
  
  If PixType eq 1 then MaxPixInd = NLeaky
  ;Now look at a few cases in the dark and flat images
  For PixInd = 0L, MaxPixInd-1 do begin

    ;Determine what type of pixels are being examined here
    If PixType eq 0 then begin 
      X = HotPixX(PixInd)
      Y = HotPixY(PixInd)
      CommentStr = strtrim(KeyStr.UpperThresh,2) + ' counts above' 
    EndIf else begin 
      X = LeakyPixX(PixInd)
      Y = LeakyPixY(PixInd)
      CommentStr = strtrim(abs(KeyStr.LowerThresh),2) + ' counts below'
    EndElse

    ;Read the Region of the Dark Cube - Dark Diff(DD), Flat Diff(FD)
    DarkIm = Fits_Read_DataCube_Region(DarkFiles(0), DarkFCB, DarkHeader,$
	 	x-1, x+1, y-1, y+1, 0, NReads-1)
    DD12 = float(DarkIm(*,*,1)) - float(DarkIm(*,*,0))
    DDNNAvg = mean([DD12(0,1),DD12(1,0),DD12(2,1),DD12(1,2)])
   
    FlatIm = Fits_Read_DataCube_Region(FlatFiles(2), FlatFCB, FlatHeader, $
		x-1, x+1, y-1, y+1, 0, 1)
    FD12 = float(FlatIm(*,*,1)) - float(FlatIm(*,*,0))
    FDNNAvg = mean([FD12(0,1),FD12(1,0),FD12(2,1),FD12(1,2)])

    ;Add the ratio of the middle pixel to the average of nearest neighbors
    If PixType eq 0 then begin 
       NNFracHotD(PixInd) = float(DD12(1,1))/DDNNAvg
       NNFracHotF(PixInd) = float(FD12(1,1))/FDNNAvg
    Endif else if PixType eq 1 then begin
       NNFracLeakyD(PixInd) = float(DD12(1,1))/DDNNAvg
       NNFracLeakyF(PixInd) = float(FD12(1,1))/FDNNAvg
    EndIf

    If TvFlag ne 0 and Avgs eq 0 then begin 
      window, xsize = 700, ysize=700 
      !p.noerase = 0
      plot, DarkIm(1,1,*), yrange = [0, 6e4], $
	Title = 'Pixel Response in Dark w/key '+DarkFileKey ,$
	ytitle = 'Counts', xtitle = 'Read', $
	position = [0.1, 0.075, 0.8, 0.45], /normal, $
        background = white, color = black
      oplot, DarkIm(0,0,*), color = Red
      oplot, DarkIm(2,2,*), color = Blue
      oplot, DarkIm(0,1,*), color = Green
      xyouts, 0.025, 0.95, 'These Pixels were flagged because they were ' +$
          CommentStr + ' their '+$
	  'perpendicular neighbors in ' + ThisFilter +' flat images', $
          color=black, /normal
 
      !p.noerase = 1
      plot, FlatIm(1,1,*), yrange = [0, 6e4], $
	Title = 'Pixel Response in Flat w/key '+ FlatFileKey, $
	ytitle = 'Counts', $
        position = [0.1, 0.525, 0.8, 0.9], /normal, $
	color = black
      oplot, FlatIm(0,0,*), color = Red
      oplot, FlatIm(2,2,*), color = Blue
      oplot, FlatIm(0,1,*), color = Green
  
      xyouts, 0.85, 0.9, 'Pixel Coordinates', color   = Black, /normal
      xyouts, 0.85, 0.85, 'X : '+strtrim(X, 2), color = Black, /normal
      xyouts, 0.85, 0.80, 'Y : '+strtrim(Y, 2), color = Black, /normal 
      xyouts, 0.85, 0.75, 'Center', color = Black, /normal
      xyouts, 0.85, 0.70, 'Lower Left', color = Red, /normal
      xyouts, 0.85, 0.65, 'Upper Right', color = Blue, /normal
      xyouts, 0.85, 0.60, 'Left', color = Green, /normal

      xyouts, 0.85, 0.20, 'Flat Change : ', color = black, /normal 
      xyouts, 0.85, 0.175, $
	strtrim(long(FlatIm(1,1,1))-long(FlatIm(1,1,0)),2), color= Black,/normal
      xyouts, 0.925, 0.175, $
        strtrim(long(FlatIm(0,0,1))-long(FlatIm(0,0,0)),2), color= Red , /normal
      xyouts, 0.85, 0.150, $
        strtrim(long(FlatIm(2,2,1))-long(FlatIm(2,2,0)),2), color= Blue ,/normal

      xyouts, 0.85, 0.10, 'Dark Change : ', color = black, /normal
      xyouts, 0.85, 0.075, $
        strtrim(long(DarkIm(1,1,1))-long(DarkIm(1,1,0)),2), color= Black,/normal
      xyouts, 0.925, 0.075, $
        strtrim(long(DarkIm(0,0,1))-long(DarkIm(0,0,0)),2), color= Red , /normal
      xyouts, 0.85, 0.050, $
        strtrim(long(DarkIm(2,2,1))-long(DarkIm(2,2,0)),2), color= Blue ,/normal
      If TvFlag eq 1 then stop
      erase
    EndIf
    
  EndFor

EndFor
  
DarkClose = Fits_Close_DataCube(DarkFCB)
FlatClose = Fits_Close_DataCube(FlatFCB)

;Plot Histograms if TvFlag eq 2 and Avgs = 1
If TvFlag eq 2 and Avgs eq 1 then begin 

  !p.noerase = 1
  !p.charsize = 0.5
  ;Plot the histogram of the leaky ratio
  LeakyDHist = histogram(NNFracLeakyD, locations = LeakyDLocs, $
	nbins = 2000, min = -2)
  Plot, LeakyDLocs, LeakyDHist, psym = 10, xrange = [-2, 5], $
	Title = 'Ratios for the Leaky Pixels in Dark Frame ' $
	        +DarkFileKey, $
	Position = [0.05, 0.05, 0.45, 0.45], /normal

  LeakyFHist = histogram(NNFracLeakyF, locations = LeakyFLocs, $
	nbins = 2000)
  Plot, LeakyFLocs, LeakyFHist, psym = 10, xrange = [ 0, .4], $
	Title = 'Ratios for the Leaky Pixels in Flat Frame ' $
		+FlatFileKey, $
	position = [0.50, 0.05, 0.95, 0.45], /normal
  ;Plot the histograms of the hot ratios
  HotDHist = histogram(NNFracHotD, locations = HotDlocs, $
	nbins = 2000, min = -5)
  Plot, HotDLocs, HotDHist, psym = 10, xrange = [-5, 40], $
	Title = 'Ratios for the Hot Pixels in Dark Frame ' $
		+DarkFileKey, $
	position = [0.05, 0.55, 0.45, 0.95], /normal

  HotFHist = histogram(NNFracHotF, locations = HotFLocs, $
	nbins = 2000, min = -2)
  Plot, HotFLocs, HotFHist, psym = 10, xrange = [-.5, 2], $
	Title = 'Ratios for the Hot Pixels in Flat Frame ' $
		+FlatFileKey, $
	position = [0.50, 0.55, 0.95, 0.95], /normal

EndIf

If TvFlag eq 2 then device,/close

set_plot,'x'
stop

end