5 import os, sys, ROOT, argparse
7 from DCubeHistograms
import MyHistoFiller
9 if __name__ ==
"__main__":
10 parser = argparse.ArgumentParser(prog=
'createDCubeSimHistograms', formatter_class=argparse.ArgumentDefaultsHelpFormatter)
11 parser.add_argument(
'-i',
'--inputFile', help=
'choose input ROOT file', default=
'NSWPRDValAlg.sim.ntuple.root', type=str)
12 parser.add_argument(
'-o',
'--outputFile', help=
'choose output ROOT file', default=
'NSWPRDValAlg.dcube.root', type=str)
13 parser.add_argument(
'--doCSC', help=
'turn off CSC if using Run4 input ROOT file', default=
False, action=
'store_true')
14 parser.add_argument(
'--doMM', help=
'turn off MM if using Run2 input ROOT file', default=
False, action=
'store_true')
15 parser.add_argument(
'--doSTGC', help=
'turn off STGC if using Run2 input ROOT file', default=
False, action=
'store_true')
16 parser.add_argument(
'--RPCsel', help=
'Choose eta_sector selections for RPC, e.g. positive_1 for positive eta and sector 1, None_None for no selection', default=
'None_None', type=str)
17 parser.add_argument(
'--MDTsel', help=
'Choose eta_sector selections for MDT, e.g. positive_1 for positive eta and sector 1, None_None for no selection', default=
'None_None', type=str)
18 parser.add_argument(
'--CSCsel', help=
'Choose eta_sector selections for CSC, e.g. positive_1 for positive eta and sector 1, None_None for no selection', default=
'None_None', type=str)
19 parser.add_argument(
'--TGCsel', help=
'Choose eta_sector selections for TGC, e.g. positive_1 for positive eta and sector 1, None_None for no selection', default=
'None_None', type=str)
20 parser.add_argument(
'--MMsel', help=
'Choose eta_sector selections for MM, e.g. positive_1 for positive eta and sector 1, None_None for no selection', default=
'None_None', type=str)
21 parser.add_argument(
'--sTGCsel', help=
'Choose eta_sector selections for sTGC, e.g. positive_1 for positive eta and sector 1, None_None for no selection', default=
'None_None', type=str)
22 Options = parser.parse_args()
24 ROOT.gROOT.SetBatch(
True)
26 if not os.path.exists(Options.inputFile):
27 print (
'ERROR: File %s does not exist'%Options.inputFile)
30 inputFile = ROOT.TFile(Options.inputFile,
"READ")
32 print (
'ERROR: Failed to open file %s'%Options.inputFile)
34 inputTree = inputFile.Get(
"NSWValTree")
36 print (
'ERROR: NSWValTree does not exist in file %s'%Options.inputFile)
39 nEntries = inputTree.GetEntries()
41 print (
'ERROR: NSWValTree of file %s has 0 entries'%Options.inputFile)
44 outputFile = ROOT.TFile(Options.outputFile,
"RECREATE")
46 print (
'ERROR: Failed to open file %s'%Options.outputFile)
50 outputFile.mkdir(
"simulation/")
51 ODir = outputFile.GetDirectory(
"simulation/")
54 RPCselections = Options.RPCsel.split(
"_")
55 RPC_eta = RPCselections[0]
56 if RPCselections[1] !=
"None":
57 RPC_sector = int (RPCselections[1])
59 RPC_sector = RPCselections[1]
61 MDTselections = Options.MDTsel.split(
"_")
62 MDT_eta = MDTselections[0]
63 if MDTselections[1] !=
"None":
64 MDT_sector = int (MDTselections[1])
66 MDT_sector = MDTselections[1]
68 CSCselections = Options.CSCsel.split(
"_")
69 CSC_eta = CSCselections[0]
70 if CSCselections[1] !=
"None":
71 CSC_sector = int (CSCselections[1])
73 CSC_sector = CSCselections[1]
75 TGCselections = Options.TGCsel.split(
"_")
76 TGC_eta = TGCselections[0]
77 if TGCselections[1] !=
"None":
78 TGC_sector = int (TGCselections[1])
80 TGC_sector = TGCselections[1]
82 MMselections = Options.MMsel.split(
"_")
83 MM_eta = MMselections[0]
84 if MMselections[1] !=
"None":
85 MM_sector = int (MMselections[1])
87 MM_sector = MMselections[1]
89 sTGCselections = Options.sTGCsel.split(
"_")
90 sTGC_eta = sTGCselections[0]
91 if sTGCselections[1] !=
"None":
92 sTGC_sector = int (sTGCselections[1])
94 sTGC_sector = sTGCselections[1]
97 for i
in range(inputTree.GetEntries()):
108 for ntruth
in range(0,len(inputTree.MuEntry_ParticlePt)):
109 truthhists += [
MyHistoFiller( chamber_name =
"TruthInfo", eta_sel =
None, sector_sel =
None )]
110 truthhists[ntruth].
fill(inputTree,ntruth)
114 if RPC_eta ==
"positive":
115 rpc_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.RPC_SIM_stationEta[nrpcHit])) >= 0
116 elif RPC_eta ==
"negative":
117 rpc_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.RPC_SIM_stationEta[nrpcHit])) < 0
119 rpc_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.RPC_SIM_stationEta[nrpcHit])) < 9
121 if RPC_sector ==
"None":
122 rpc_sector_sel =
lambda s: MyHistoFiller.Eta(ord(s.RPC_SIM_stationPhi[nrpcHit])) < 10
124 rpc_sector_sel =
lambda s: MyHistoFiller.Eta(ord(s.RPC_SIM_stationPhi[nrpcHit])) == RPC_sector
126 for nrpcHit
in range(0,len(inputTree.RPC_SIM_localPostionX)):
127 rpchists += [
MyHistoFiller( chamber_name =
"RPC_Sim", eta_sel = rpc_eta_sel, sector_sel = rpc_sector_sel )]
128 rpchists[nrpcHit].
fill(inputTree, nrpcHit)
131 if MDT_eta ==
"positive":
132 mdt_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.MDT_Sim_stationEta[nmdtHit])) >= 0
133 elif MDT_eta ==
"negative":
134 mdt_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.MDT_Sim_stationEta[nmdtHit])) < 0
136 mdt_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.MDT_Sim_stationEta[nmdtHit])) < 9
138 if MDT_sector ==
"None":
139 mdt_sector_sel =
lambda s: MyHistoFiller.Eta(ord(s.MDT_Sim_stationPhi[nmdtHit])) < 10
141 mdt_sector_sel =
lambda s: MyHistoFiller.Eta(ord(s.MDT_Sim_stationPhi[nmdtHit])) == MDT_sector
143 for nmdtHit
in range(0,len(inputTree.MDT_Sim_hitLocalPositionX)):
144 mdthists += [
MyHistoFiller( chamber_name =
"MDT_Sim", eta_sel = mdt_eta_sel, sector_sel = mdt_sector_sel )]
145 mdthists[nmdtHit].
fill(inputTree, nmdtHit)
148 if Options.doCSC ==
True:
149 if CSC_eta ==
"positive":
150 csc_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.CSC_Sim_stationEta[ncscHit])) >= 0
151 elif CSC_eta ==
"negative":
152 csc_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.CSC_Sim_stationEta[ncscHit])) < 0
154 csc_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.CSC_Sim_stationEta[ncscHit])) < 9
156 if CSC_sector ==
"None":
157 csc_sector_sel =
lambda s: MyHistoFiller.Eta(ord(s.CSC_Sim_stationPhi[ncscHit])) < 10
159 csc_sector_sel =
lambda s: MyHistoFiller.Eta(ord(s.CSC_Sim_stationPhi[ncscHit])) == CSC_sector
161 for ncscHit
in range(0,len(inputTree.CSC_Sim_hitGlobalPositionX)):
162 cschists += [
MyHistoFiller( chamber_name =
"CSC_Sim", eta_sel = csc_eta_sel, sector_sel = csc_sector_sel )]
163 cschists[ncscHit].
fill(inputTree, ncscHit)
166 if TGC_eta ==
"positive":
167 tgc_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.TGC_Sim_stationEta[ntgcHit])) >= 0
168 elif TGC_eta ==
"negative":
169 tgc_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.TGC_Sim_stationEta[ntgcHit])) < 0
171 tgc_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.TGC_Sim_stationEta[ntgcHit])) < 9
173 if TGC_sector ==
"None":
174 tgc_sector_sel =
lambda s: MyHistoFiller.Eta(ord(s.TGC_Sim_stationPhi[ntgcHit])) < 51
176 tgc_sector_sel =
lambda s: MyHistoFiller.Eta(ord(s.TGC_Sim_stationPhi[ntgcHit])) == TGC_sector
178 for ntgcHit
in range(0,len(inputTree.TGC_Sim_hitLocalPositionX)):
179 tgchists += [
MyHistoFiller( chamber_name =
"TGC_Sim", eta_sel = tgc_eta_sel, sector_sel = tgc_sector_sel )]
180 tgchists[ntgcHit].
fill(inputTree, ntgcHit)
183 if Options.doMM ==
True:
184 if MM_eta ==
"positive":
185 mm_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.Hits_MM_off_stationEta[nmmHit])) >= 0
186 elif MM_eta ==
"negative":
187 mm_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.Hits_MM_off_stationEta[nmmHit])) < 0
189 mm_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.Hits_MM_off_stationEta[nmmHit])) < 9
191 if MM_sector ==
"None":
192 mm_sector_sel =
lambda s: MyHistoFiller.Eta(ord(s.Hits_MM_off_stationPhi[nmmHit])) < 10
194 mm_sector_sel =
lambda s: MyHistoFiller.Eta(ord(s.Hits_MM_off_stationPhi[nmmHit])) == MM_sector
196 for nmmHit
in range(0,len(inputTree.Hits_MM_GlobalPositionX)):
197 mmhists += [
MyHistoFiller( chamber_name =
"MM_Sim", eta_sel = mm_eta_sel, sector_sel = mm_sector_sel )]
198 mmhists[nmmHit].
fill(inputTree, nmmHit)
201 if Options.doSTGC ==
True:
202 if sTGC_eta ==
"positive":
203 stgc_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.Hits_sTGC_off_stationEta[nstgcHit])) >= 0
204 elif sTGC_eta ==
"negative":
205 stgc_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.Hits_sTGC_off_stationEta[nstgcHit])) < 0
207 stgc_eta_sel =
lambda t: MyHistoFiller.Eta(ord(t.Hits_sTGC_off_stationEta[nstgcHit])) < 9
209 if sTGC_sector ==
"None":
210 stgc_sector_sel =
lambda s: MyHistoFiller.Eta(ord(s.Hits_sTGC_off_stationPhi[nstgcHit])) < 10
212 stgc_sector_sel =
lambda s: MyHistoFiller.Eta(ord(s.Hits_sTGC_off_stationPhi[nstgcHit])) == sTGC_sector
214 for nstgcHit
in range(0,len(inputTree.Hits_sTGC_GlobalPositionX)):
215 stgchists += [
MyHistoFiller( chamber_name =
"sTGC_Sim", eta_sel = stgc_eta_sel, sector_sel = stgc_sector_sel )]
216 stgchists[nstgcHit].
fill(inputTree, nstgcHit)
220 truthhist =
MyHistoFiller( chamber_name =
"TruthInfo", eta_sel =
None, sector_sel =
None )
221 truthhist.write(ODir)
223 rpchist =
MyHistoFiller( chamber_name =
"RPC_Sim", eta_sel =
None, sector_sel =
None )
226 mdthist =
MyHistoFiller( chamber_name =
"MDT_Sim", eta_sel =
None, sector_sel =
None )
229 if Options.doCSC ==
True:
230 cschist =
MyHistoFiller( chamber_name =
"CSC_Sim", eta_sel =
None, sector_sel =
None )
233 tgchist =
MyHistoFiller( chamber_name =
"TGC_Sim", eta_sel =
None, sector_sel =
None )
236 if Options.doMM ==
True:
237 mmhist =
MyHistoFiller( chamber_name =
"MM_Sim", eta_sel =
None, sector_sel =
None )
240 if Options.doSTGC ==
True:
241 stgchist =
MyHistoFiller( chamber_name =
"sTGC_Sim", eta_sel =
None, sector_sel =
None )