diff --git a/Tracker/TrackerRecAlgs/PairVertex/scripts/NtupleReader.py b/Tracker/TrackerRecAlgs/PairVertex/scripts/NtupleReader.py new file mode 100644 index 0000000000000000000000000000000000000000..590370734773a812764c99d1ee205afaf4beae9f --- /dev/null +++ b/Tracker/TrackerRecAlgs/PairVertex/scripts/NtupleReader.py @@ -0,0 +1,94 @@ +#!/usr/bin/env python +from PlotHelper import * +from Utilities import * + +#CUTS +REALCUT={ + 'longTracks':[2,None,None,None],} +# 'Track_nLayers':[6,None,None,None], +# 'Track_nDoF':[0,None,None,None], +# 'Track_Chi2':[0,None,None,None], +# 'Preshower0_charge':[0,None,None,None], +# 'Preshower1_charge':[0,None,None,None] +# } + +REALLAMBDAS=[ + lambda event: event.Timing0_charge+Timing1_charge+Timing2_charge+Timing3_charge>0, + lambda event: event.Calo0_charge+Calo1_charge+Calo2_charge+Calo3_charge>0, + lambda event: all([event.Track_nLayers[i]>6 for i in event.longTracks]) +] + +def GeneralPlotter(arrDict,dcuts=[],lcuts=[],showplot=True,savefig=False): + nans=0 + + for event in t: + if ApplyAllCuts(event,dcuts,lcuts): + for name in arrDict.keys(): + value=nan + if 'lkwargs' in list(arrDict[name].keys()): + value = arrDict[name]['lambda'](event,**arrDict[name]['lkwargs']) + else: value = arrDict[name]['lambda'](event) + checkReal=True + if type(value)==list: checkReal = all([v != inf and v !=nan and v !=-1*inf for v in value]) + else: checkReal = value != inf and value !=nan and value !=-1*inf + + if checkReal and type(value)==list: + if arrDict[name]['plotterfunc'] in [configBasicHist]: arrDict[name]['array']+=value + else: + for iv in range(len(value)): + arrDict[name]['array'][iv]+=[value[iv]] + elif checkReal: arrDict[name]['array']+=[value] + else: nans+=1 + + print(f"# of NaNs:{nans}") + for name in arrDict.keys(): + isPrintout= arrDict[name]['plotterfunc'] in [PrintBoolPercent] + print(name) + print('LEN:',len(arrDict[name]['array'])) + arrDict[name]['plotterfunc']( arrDict[name]['array'], **arrDict[name]['plotkwargs']) + if not isPrintout: + plt.title(user_input+name) + if savefig: plt.savefig(user_input+name) + if showplot: plt.show() +COUNT=0 +def ReadoutTest(event): + print(len(event.Track_nDoF)) + global COUNT + COUNT+=1 + if COUNT>100: raise KeyboardInterrupt + +# +MainArrDict={ + 'CutTest':{ + 'array':[], + 'lambda': ApplyAllCuts, + 'lkwargs':{'dcuts':[REALCUT],'lcuts':[]}, + 'plotterfunc': PrintBoolPercent, + 'plotkwargs' : {}, + 'descr':'Check how much a cut cuts.' + }, + 'BasicHist':{ + 'array':[], + 'lambda': ReadoutTest, + 'plotterfunc': configBasicHist, + 'plotkwargs' : {}, + 'descr':'Single variable tests.' + },} + +TestArrDict=dict() +testname='BasicHist' #Any key from MainArrDict +TestArrDict[testname]=MainArrDict[testname] + + + +#EXECUTE +if __name__ == "__main__": + DP_fnames=['Faser-Physics-008301-00000-00049-r0011-PHYS'] + for f in DP_fnames: + user_input=f + t = TChain("nt") + nfiles = 0 + nfiles += t.Add("/eos/home-s/sshively/CALYPZONE/custpython/"+user_input+".root") + print(f) + + GeneralPlotter(TestArrDict,dcuts=[{}],lcuts=[]) \ No newline at end of file diff --git a/Tracker/TrackerRecAlgs/PairVertex/scripts/Utilities.py b/Tracker/TrackerRecAlgs/PairVertex/scripts/Utilities.py index ff711357eed1ed411749dcf9ed2e800546f358c3..43eba0369ec626cf87400a803d4ca6d0655dbdce 100644 --- a/Tracker/TrackerRecAlgs/PairVertex/scripts/Utilities.py +++ b/Tracker/TrackerRecAlgs/PairVertex/scripts/Utilities.py @@ -21,6 +21,9 @@ BASICCUT={ 'TrackCount':[0,None,None,None] } + + + def BoolDaughterCut(event): if event.TrackCount>=2 and len(event.tracksTruthBarcode)>=2: return 2 in event.tracksTruthBarcode and 3 in event.tracksTruthBarcode diff --git a/Tracker/TrackerRecAlgs/PairVertex/src/PairVertexAlg.cxx b/Tracker/TrackerRecAlgs/PairVertex/src/PairVertexAlg.cxx index be6c6fcb17102fb3631c5283a36e155a89603e2b..f923da18ece0745adf49977796b1f95b1ed94d5a 100644 --- a/Tracker/TrackerRecAlgs/PairVertex/src/PairVertexAlg.cxx +++ b/Tracker/TrackerRecAlgs/PairVertex/src/PairVertexAlg.cxx @@ -445,6 +445,8 @@ StatusCode PairVertexAlg::execute(const EventContext &ctx) const int clustercount=0; //loop over all clusters + if (clusterContainer->size()>0) + { for (auto collection : *clusterContainer) { Identifier id = collection->identify(); @@ -463,7 +465,7 @@ StatusCode PairVertexAlg::execute(const EventContext &ctx) const //loop over all RDOs in the cluster for (auto rdoID : clusterRDOList) { - if (simDataCollection->count(idRDO) > 0) + if (h_collectionMap->count(rdoID) > 0) { @@ -475,26 +477,35 @@ StatusCode PairVertexAlg::execute(const EventContext &ctx) const { barcode = depositPair.first->barcode(); - // int pdgID = 0//depositPair.first->pdg_id(); - // float edeposit = 0.0//depositPair.second; + int pdgID = depositPair.first->pdg_id(); + float edeposit = depositPair.second; - // TotalEnergyDeposited+=edeposit; //sum all deposited energies in RDO for each cluster + TotalEnergyDeposited+=edeposit; //sum all deposited energies in RDO for each cluster - // if (barcode==2 || barcode==3) //its a daughter, do something. - // { - // isDaughter=true; - // if (pdgID>0){TotalEnergyDeposited_daughterNeg+=edeposit;} - // if (pdgID<0){TotalEnergyDeposited_daughterPos+=edeposit;} - // } + if (barcode==2 || barcode==3) //its a daughter, do something. + { + isDaughter=true; + if (pdgID>0){TotalEnergyDeposited_daughterNeg+=edeposit;} + if (pdgID<0){TotalEnergyDeposited_daughterPos+=edeposit;} + } } - // m_clusterETotPos.push_back(TotalEnergyDeposited_daughterPos); - // m_clusterETotNeg.push_back(TotalEnergyDeposited_daughterNeg); - // m_clusterETot.push_back(TotalEnergyDeposited); - // m_clusterEFracPos.push_back(TotalEnergyDeposited_daughterPos/TotalEnergyDeposited); - // m_clusterEFracNeg.push_back(TotalEnergyDeposited_daughterNeg/TotalEnergyDeposited); - // m_clusterLocation.push_back(clusterLocation); - // m_clusterIsDaughter.push_back(isDaughter); + m_clusterETotPos.push_back(TotalEnergyDeposited_daughterPos); + m_clusterETotNeg.push_back(TotalEnergyDeposited_daughterNeg); + m_clusterETot.push_back(TotalEnergyDeposited); + m_clusterLocation.push_back(clusterLocation); + m_clusterIsDaughter.push_back(isDaughter); + if (TotalEnergyDeposited>0) + { + m_clusterEFracPos.push_back(TotalEnergyDeposited_daughterPos/TotalEnergyDeposited); + m_clusterEFracNeg.push_back(TotalEnergyDeposited_daughterNeg/TotalEnergyDeposited); + } + else + { + m_clusterEFracPos.push_back(NaN); + m_clusterEFracNeg.push_back(NaN); + } + } } clustercount++; @@ -502,6 +513,7 @@ StatusCode PairVertexAlg::execute(const EventContext &ctx) const m_clusterCount=clustercount; } + } m_tree->Fill();