#PlotEnergyHistogram.py			Lance Simms, LLNL 2010
#
#PURPOSE:
#  To plot the energy histogram for a CZTCubeSat run.
#  The energy histogram will either represent the number of electrons 
#  counted or the actual energy reported by GEANT4.
#
#KEYWORDS:
#  SimFileName: str
#    The relative path to the file that contains the data to be histogrammed
#    Assumed to live in /Users/lances/Desktop/Software/GEANT4/work/CZTCubeSat/
#  ExpFileName: str
#    The relative path to the .dat file that shows Garrett's data
#    Assumed to live in /Users/lances/Desktop/PoloSat/Code/
#
#EXAMPLES:
#  run PlotEnergyHistogram.py --ExpFileName=SampleData_r20c12_def_NoCuts.dat
#
##########
#
#Get the filename and options from command line
SaveFig=1
import os
execfile(os.path.expanduser('~')+'/Desktop/Software/Python/CZTCubeSat/CZTCubeSatSetup.py')
execfile(os.path.expanduser('~')+'/Desktop/Software/Python/CZTCubeSat/BigFontPlotSettings.py')

from ReturnPlottableHistogram import *
from ReturnRebinnedHistogram import *
import getopt
import pdb

##PARSE KEYWORDS
keywords    = ['SimFileName=','ExpFileName=']

#DEFAULT VALUES FOR KEYWORDS
SimFileName           = 'EHistogram.dat'
ExpFileName           = 'Out.dat'

#NOW PARSE KEYWORDS
opts, extraparams = getopt.getopt(sys.argv[1:],'',keywords)
for o,p in opts:
  if o in ['--SimFileName']:
    SimFileName = p
  elif o in ['--ExpFileName']:
    ExpFileName = p

E,Num = readcol('/Users/lances/Desktop/Software/GEANT4/work/CZTCubeSat/'+SimFileName, twod=False)

#Make it into histogram values
EHist, NumHist = ReturnPlottableHistogram(E, Num)

#Allow to Rebin
RebinNum = 4
EHistRebinned, NumHistRebinned = ReturnRebinnedHistogram(E, Num, RebinNum)
EHistRebinned, NumHistRebinned = ReturnPlottableHistogram(EHistRebinned, NumHistRebinned)

#Now histogram the values
PlotNoHist=0
if PlotNoHist:
  mplot.figure(0)
  mplot.clf()
  mplot.plot(EHist, NumHist)
  mplot.xlabel('\textbf{Energy (keV)}', fontsize=16)
  mplot.ylabel('Count')
  mplot.xlim([0,62])

#Now histogram the values
mplot.figure(1)
mplot.clf()
mplot.plot(EHistRebinned, NumHistRebinned, label='GEANT4 Simulation')
mplot.xlabel(r'Energy (keV)', fontsize=16, fontweight='bold')
mplot.ylabel('Count')
mplot.xlim([0,62])

if SaveFig == 1:
  mplot.savefig("EnergyHistogram.png")

############################################
#Choose whether or not to plot Garrett's data
PlotRealData = 1
if PlotRealData == 1:
  LineNum, EC, EMax, EL, ER, EBL, EB, EBR, ETL, ET, ETR, ECut1, ECut2 = \
     readcol('/Users/lances/Desktop/PoloSat/Code/'+ExpFileName, twod=False)

#Multiply the data by 668.539 to get the numbers in eV (0.089 Volts ~ 59.5 keV)
ConversionFactor = 668.539
EC = ConversionFactor*EC

MinHist=0
MaxHist=70
N, Bins = histogram(EC, 120, range=[MinHist,MaxHist])

MinHist = 0
NBins   = 180
Bins    = 0.4*frange(NBins)+MinHist+.15
Bins, N = histOutline.histOutline(EC, binsIn=Bins)

mplot.plot(Bins,N, label='Real CZT Data')
mplot.xlim([0,80])
mplot.legend(loc=[.1,.8])
if SaveFig == 1:
  mplot.savefig("EnergyHistogramWRealData.png")

#Choose whether or not to plot Garrett's data
if PlotRealData == 1:
  LineNum, EC, EMax, EL, ER, EBL, EB, EBR, ETL, ET, ETR, ECut1, ECut2 = \
     readcol('/Users/lances/Desktop/PoloSat/Code/'+ExpFileName, twod=False)

#Try to remove offsets and plot sums of pixels
CenRowOffset = 0.022
#EC = EC-CenRowOffset
EL = EL-CenRowOffset
ER = ER-CenRowOffset

#Now remove the offsets from the upper and bottom rows
BotRowOffset = 0.094
EB  = EB-BotRowOffset
EBLOffset = 0.0935
EBL = EBL-EBLOffset
EBROffset = 0.094
EBR = EBR-EBROffset

TopRowOffset  = 0.0925
ET  = ET-TopRowOffset
ETLOffset = 0.0939
ETL = ETL-ETLOffset
ETROffset = 0.0915
ETR = ETR-ETROffset

ETotal = EC+EL+ER+EB+EBL+EBR+ET+ETL+ETR

##Go through and get the two max pixels
EMaxArr  = zeros(size(EC), dtype=float)
EMax2Arr = zeros(size(EC), dtype=float)
E2TotalArr = zeros(size(EC), dtype=float)
for i in arange(size(EC)):
  TheseEs = [EC[i], ET[i], EB[i], EL[i], ER[i], ETL[i], ETR[i], EBL[i], EBR[i]]
  TheseEsSorted = numpy.sort(TheseEs)
  EMaxArr[i]    = TheseEsSorted[8]
  EMax2Arr[i]   = TheseEsSorted[7]
  E2TotalArr[i] = EMaxArr[i]+EMax2Arr[i]

#Bin the Total Energy
ETotal = ETotal*ConversionFactor

MinHist = -30
NBins   = 180
Bins    = 0.5*frange(NBins)+MinHist
Bins, N = histOutline.histOutline(ETotal, binsIn=Bins)
 
mplot.figure(2)
mplot.clf()
mplot.plot(Bins, N, label='All 9 pixels')

#Bin the Top 2 pixels
E2TotalArr = E2TotalArr*ConversionFactor

MinHist = -10
NBins   = 180
Bins    = 0.5*frange(NBins)+MinHist
Bins, N = histOutline.histOutline(E2TotalArr, binsIn=Bins)

mplot.plot(Bins, N, label='Top 2 pixels')
mplot.legend(loc=[.1,.8])
if SaveFig == 1:
  mplot.savefig("EnergyHistogramNoCuts_DoubleAnd9Pixels.png")

#Choose whether or not to plot Garrett's data
if PlotRealData == 1:
  LineNum, EC, EMax, EL, ER, EBL, EB, EBR, ETL, ET, ETR, ECut1, ECut2 = \
     readcol('/Users/lances/Desktop/PoloSat/Code/'+ExpFileName, twod=False)

