SND@LHC Software
Loading...
Searching...
No Matches
MufluxReco.py
Go to the documentation of this file.
1#!/usr/bin/env python
2#inputFile = '/eos/experiment/ship/data/muflux/run_fixedtarget/19april2018/pythia.root'
3#geoFile = '/eos/experiment/ship/data/muflux/run_fixedtarget/19april2018/geofile_full.root'
4from __future__ import print_function
5from __future__ import division
6import global_variables
7debug = False#False
8
9withNoStrawSmearing = None # True for debugging purposes
10withDist2Wire = False
11withNTaggerHits = 0
12nEvents = 10000
13firstEvent = 0
14withHists = True
15dy = None
16saveDisk = False # remove input file
17realPR = ''
18withT0 = False
19
20import resource
21def mem_monitor():
22 # Getting virtual memory size
23 pid = os.getpid()
24 with open(os.path.join("/proc", str(pid), "status")) as f:
25 lines = f.readlines()
26 _vmsize = [l for l in lines if l.startswith("VmSize")][0]
27 vmsize = int(_vmsize.split()[1])
28 #Getting physical memory size
29 pmsize = resource.getrusage(resource.RUSAGE_SELF).ru_maxrss
30 print("memory: virtuell = %5.2F MB physical = %5.2F MB"%(vmsize/1.0E3,pmsize/1.0E3))
31
32import ROOT,os,sys,getopt
33import rootUtils as ut
34import shipunit as u
35import shipRoot_conf
36
37# geoMat = ROOT.genfit.TGeoMaterialInterface()
38shipRoot_conf.configure()
39
40try:
41 opts, args = getopt.getopt(sys.argv[1:], "o:D:FHPu:n:f:g:c:hqv:sl:A:Y:i:t:",\
42 ["ecalDebugDraw","inputFile=","geoFile=","nEvents=","noStrawSmearing","noVertexing","saveDisk","realPR","withT0", "withNTaggerHits=", "withDist2Wire"])
43except getopt.GetoptError:
44 # print help information and exit:
45 print(' enter --inputFile= --geoFile= --nEvents= --firstEvent=,')
46 print(' noStrawSmearing: no smearing of distance to wire, default on')
47 print(' outputfile will have same name with _rec added')
48 sys.exit()
49for o, a in opts:
50 if o in ("--noVertexing",):
51 print("WARNING: --noVertexing option not currently used by script.")
52 if o in ("--noStrawSmearing",):
53 withNoStrawSmearing = True
54 if o in ("--withT0",):
55 withT0 = True
56 if o in ("--withDist2Wire",):
57 withDist2Wire = True
58 if o in ("-t", "--withNTaggerHits"):
59 withNTaggerHits = int(a)
60 if o in ("-f", "--inputFile"):
61 inputFile = a
62 if o in ("-g", "--geoFile"):
63 geoFile = a
64 if o in ("-n", "--nEvents="):
65 nEvents = int(a)
66 if o in ("-Y",):
67 dy = float(a)
68 if o in ("--ecalDebugDraw",):
69 print("WARNING: --ecalDebugDraw option not currently used by script.")
70 if o in ("--saveDisk",):
71 saveDisk = True
72 if o in ("--realPR",):
73 realPR = "_PR"
74
75
76# need to figure out which geometry was used
77if not dy:
78 # try to extract from input file name
79 tmp = inputFile.split('.')
80 try:
81 dy = float( tmp[1]+'.'+tmp[2] )
82 except:
83 dy = None
84print('configured to process ', nEvents, ' events from ', inputFile,
85 ' starting with event ', firstEvent, ' with option Yheight = ', dy,
86 ' with vertexing', False,' and real pattern reco',realPR=="_PR")
87if not inputFile.find('_rec.root') < 0:
88 outFile = inputFile
89 inputFile = outFile.replace('_rec.root','.root')
90else:
91 outFile = inputFile.replace('.root','_rec.root')
92# outfile should be in local directory
93 tmp = outFile.split('/')
94 outFile = tmp[len(tmp)-1]
95 if saveDisk: os.system('mv '+inputFile+' '+outFile)
96 else : os.system('cp '+inputFile+' '+outFile)
97
98# check if simulation or raw data
99f=ROOT.TFile.Open(outFile)
100if f.cbmsim.FindBranch('MCTrack'): simulation = True
101else: simulation = False
102f.Close()
103
104if simulation and geoFile:
105 fgeo = ROOT.TFile.Open(geoFile)
106 from ShipGeoConfig import ConfigRegistry
107 from rootpyPickler import Unpickler
108#load Shipgeo dictionary
109 upkl = Unpickler(fgeo)
110 ShipGeo = upkl.load('ShipGeo')
111else:
112 ShipGeo = ConfigRegistry.loadpy("$FAIRSHIP/geometry/charm-geometry_config.py", Yheight = dy)
113
114h={}
115log={}
116if withHists:
117 ut.bookHist(h,'distu','distance to wire',100,0.,5.)
118 ut.bookHist(h,'distv','distance to wire',100,0.,5.)
119 ut.bookHist(h,'disty','distance to wire',100,0.,5.)
120 ut.bookHist(h,'nmeas','nr measurements',100,0.,50.)
121 ut.bookHist(h,'chi2','Chi2/DOF',100,0.,20.)
122 ut.bookHist(h,'p-fittedtracks','p of fitted tracks',40,0.,400.)
123 ut.bookHist(h,'1/p-fittedtracks','1/p of fitted tracks',120,-0.2,1.)
124 ut.bookHist(h,'pt-fittedtracks','pt of fitted tracks',100,0.,10.)
125 ut.bookHist(h,'1/pt-fittedtracks','1/pt of fitted tracks',120,-0.2,1.)
126 ut.bookHist(h,'ptruth','ptruth',40,0.,400.)
127 ut.bookHist(h,'delPOverP','Pfitted/Ptrue-1 vs Ptrue',40,0.,400.,50,-2.0,2.0)
128 ut.bookHist(h,'invdelPOverP','1/Pfitted-1/Ptrue)/(1/Ptrue) vs Ptrue',40,0.,400.,50,-2.0,2.0)
129 ut.bookProf(h,'deltaPOverP','Pfitted/Ptrue-1 vs Ptrue',40,0.,400.,-10.,10.0)
130 ut.bookHist(h,'Pfitted-Pgun','P-fitted vs P-gun',40,0.,400.,50,0.,500.0)
131 ut.bookHist(h,'Px/Pzfitted','Px/Pz-fitted',100,-0.04,0.04)
132 ut.bookHist(h,'Py/Pzfitted','Py/Pz-fitted',100,-0.04,0.04)
133 ut.bookHist(h,'Px/Pztrue','Px/Pz-true',100,-0.04,0.04)
134 ut.bookHist(h,'Py/Pztrue','Py/Pz-true',100,-0.04,0.04)
135 ut.bookHist(h,'Px/Pzfitted-noT4','Px/Pz-fitted only T1,T2,T3 ',100,-0.04,0.04)
136 ut.bookHist(h,'Py/Pzfitted-noT4','Py/Pz-fitted only T1,T2,T3',100,-0.04,0.04)
137 ut.bookHist(h,'Px/Pztrue-noT4','Px/Pz-true only T1,T2,T3',100,-0.04,0.04)
138 ut.bookHist(h,'Py/Pztrue-noT4','Py/Pz-true only T1,T2,T3',100,-0.04,0.04)
139 ut.bookHist(h,'Px/Pzfitted-all','Px/Pz-fitted',100,-0.04,0.04)
140 ut.bookHist(h,'Py/Pzfitted-all','Py/Pz-fitted',100,-0.04,0.04)
141 ut.bookHist(h,'Px/Pztrue-all','Px/Pz-true',100,-0.04,0.04)
142 ut.bookHist(h,'Py/Pztrue-all','Py/Pz-true',100,-0.04,0.04)
143 ut.bookHist(h,'Px/Pzfitted-Px/Pztruth','Px/Pz-fitted - Px/Pz-true vs P-true',40,0.,400.,100,-0.002,0.002)
144 ut.bookHist(h,'Py/Pzfitted-Py/Pztruth','Py/Pz-fitted - Py/Pz-true vs P-true',40,0.,400.,50,-0.02,0.02)
145 ut.bookHist(h,'Px/Pzfitted-Px/Pztruth-noT4','Px/Pz-fitted - Px/Pz-true vs P-true (no stereo layers)',40,0.,400.,100,-0.002,0.002)
146 ut.bookHist(h,'Py/Pzfitted-Py/Pztruth-noT4','Py/Pz-fitted - Py/Pz-true vs P-true (no stereo layers)',40,0.,400.,50,-0.02,0.02)
147 ut.bookHist(h,'Px/Pzfitted-Px/Pztruth-all','Px/Pz-fitted - Px/Pz-true vs P-true (no stereo layers)',40,0.,400.,100,-0.002,0.002)
148 ut.bookHist(h,'Py/Pzfitted-Py/Pztruth-all','Py/Pz-fitted - Py/Pz-true vs P-true (no stereo layers)',40,0.,400.,50,-0.02,0.02)
149
150 ut.bookHist(h,'p-value','p-value of fit',100,0.,1.)
151 ut.bookHist(h,'pt-kick','pt-kick',100,-2.,2.)
152 ut.bookHist(h,'nmeas-noT4','nr measurements only T1,T2,T3',100,0.,50.)
153 ut.bookHist(h,'chi2-noT4','Chi2/DOF only T1,T2,T3',100,0.,20.)
154 ut.bookHist(h,'nmeas-all','nr measurements',100,0.,50.)
155 ut.bookHist(h,'chi2-all','Chi2/DOF',100,0.,20.)
156 ut.bookHist(h,'p-fittedtracks-noT4','p of fitted track only T1,T2,T3',40,0.,400.)
157 ut.bookHist(h,'1/p-fittedtracks-noT4','1/p of fitted tracks only T1,T2,T3',120,-0.2,1.)
158 ut.bookHist(h,'pt-fittedtracks-noT4','pt of fitted tracks only T1,T2,T3',100,0.,10.)
159 ut.bookHist(h,'1/pt-fittedtracks-noT4','1/pt of fitted tracks only T1,T2,T3',120,-0.2,1.)
160 ut.bookHist(h,'ptruth-noT4','ptruth only T1,T2,T3',40,0.,400.)
161 ut.bookHist(h,'delPOverP-noT4','Pfitted/Ptrue-1 vs Ptrue only T1,T2,T3',40,0.,400.,50,-2.0,2.0)
162 ut.bookHist(h,'invdelPOverP-noT4','1/Pfitted-1/Ptrue)/(1/Ptrue) vs Ptrue only T1,T2,T3',40,0.,400.,50,-2.0,2.0)
163 ut.bookProf(h,'deltaPOverP-noT4','Pfitted/Ptrue-1 vs Ptrue only T1,T2,T3',40,0.,400.,-10.,10.0)
164 ut.bookHist(h,'Pfitted-Pgun-noT4','P-fitted vs P-gun only T1,T2,T3',40,0.,400.,50,0.,500.0)
165 ut.bookHist(h,'p-value-noT4','p-value of fit only T1,T2,T3',100,0.,1.)
166 ut.bookHist(h,'p-fittedtracks-all','p of fitted track',40,0.,400.)
167 ut.bookHist(h,'1/p-fittedtracks-all','1/p of fitted tracks',120,-0.2,1.)
168 ut.bookHist(h,'pt-fittedtracks-all','pt of fitted tracks',100,0.,10.)
169 ut.bookHist(h,'1/pt-fittedtracks-all','1/pt of fitted tracks',120,-0.2,1.)
170 ut.bookHist(h,'ptruth-all','ptruth',40,0.,400.)
171 ut.bookHist(h,'delPOverP-all','Pfitted/Ptrue-1 vs Ptrue',40,0.,400.,50,-2.0,2.0)
172 ut.bookHist(h,'invdelPOverP-all','1/Pfitted-1/Ptrue)/(1/Ptrue) vs Ptrue',40,0.,400.,50,-2.0,2.0)
173 ut.bookProf(h,'deltaPOverP-all','Pfitted/Ptrue-1 vs Ptrue',40,0.,400.,-10.,10.0)
174 ut.bookHist(h,'Pfitted-Pgun-all','P-fitted vs P-gun only',40,0.,400.,50,0.,500.0)
175 ut.bookHist(h,'p-value-all','p-value of fit',100,0.,1.)
176 ut.bookHist(h,'hits-T1','x vs y hits in T1',50,-25.,25.,100,-50.,50)
177 ut.bookHist(h,'hits-T2','x vs y hits in T2',50,-25.,25.,100,-50.,50)
178 ut.bookHist(h,'hits-T1x','x vs y hits in T1 x plane',50,-25.,25.,100,-50.,50)
179 ut.bookHist(h,'hits-T1u','x vs y hits in T1 u plane',50,-25.,25.,100,-50.,50)
180 ut.bookHist(h,'hits-T2x','x vs y hits in T2 x plane',50,-25.,25.,100,-50.,50)
181 ut.bookHist(h,'hits-T2v','x vs y hits in T2 v plane',50,-25.,25.,100,-50.,50)
182 ut.bookHist(h,'hits-T3','x vs y hits in T3',200,-100.,100.,160,-80.,80)
183 ut.bookHist(h,'hits-T4','x vs y hits in T4',200,-100.,100.,160,-80.,80)
184
185 ut.bookHist(h,'muontaggerhits', 'Muon Tagger Points', 300, -150, 150, 200, -100, 100)
186 h['muontaggerhits'].SetMarkerSize(15)
187
188 ut.bookHist(h, 'muontagger_z', 'Z Hits', 600, 850, 2500)
189 ut.bookHist(h, 'muontaggerdist', 'Muontagger Hits', 300, -150, 150, 200, -100, 100, 600, 850, 2000)
190
191
192 ut.bookHist(h, 'muontagger_clusters', 'Clusters', 50, 0, 50)
193
194 ut.bookHist(h,'NTrueTracks','Number of tracks.', 3, -0.5, 2.5)
195 h['NTrueTracks'].GetXaxis().SetBinLabel(1,"Stations 1&2, Y views")
196 h['NTrueTracks'].GetXaxis().SetBinLabel(2,"Stations 1&2, Stereo views")
197 h['NTrueTracks'].GetXaxis().SetBinLabel(3,"Stations 3&4")
198
199 ut.bookHist(h,'Reco_y12','Number of recognized tracks, clones and ghosts in stations 1&2, Y views', 5, -0.5, 4.5)
200 h['Reco_y12'].GetXaxis().SetBinLabel(1,"N total")
201 h['Reco_y12'].GetXaxis().SetBinLabel(2,"N recognized tracks")
202 h['Reco_y12'].GetXaxis().SetBinLabel(3,"N clones")
203 h['Reco_y12'].GetXaxis().SetBinLabel(4,"N ghosts")
204 h['Reco_y12'].GetXaxis().SetBinLabel(5,"N others")
205
206 ut.bookHist(h,'Reco_stereo12','Number of recognized tracks, clones and ghosts in stations 1&2, Stereo views', 5, -0.5, 4.5)
207 h['Reco_stereo12'].GetXaxis().SetBinLabel(1,"N total")
208 h['Reco_stereo12'].GetXaxis().SetBinLabel(2,"N recognized tracks")
209 h['Reco_stereo12'].GetXaxis().SetBinLabel(3,"N clones")
210 h['Reco_stereo12'].GetXaxis().SetBinLabel(4,"N ghosts")
211 h['Reco_stereo12'].GetXaxis().SetBinLabel(5,"N others")
212
213 ut.bookHist(h,'Reco_34','Number of recognized tracks, clones and ghosts in stations 3&4', 5, -0.5, 4.5)
214 h['Reco_34'].GetXaxis().SetBinLabel(1,"N total")
215 h['Reco_34'].GetXaxis().SetBinLabel(2,"N recognized tracks")
216 h['Reco_34'].GetXaxis().SetBinLabel(3,"N clones")
217 h['Reco_34'].GetXaxis().SetBinLabel(4,"N ghosts")
218 h['Reco_34'].GetXaxis().SetBinLabel(5,"N others")
219
220
221
222 ut.bookHist(h,'NTrueTracks_3hits','Number of tracks with more than 3 hits.', 3, -0.5, 2.5)
223 h['NTrueTracks_3hits'].GetXaxis().SetBinLabel(1,"Stations 1&2, Y views")
224 h['NTrueTracks_3hits'].GetXaxis().SetBinLabel(2,"Stations 1&2, Stereo views")
225 h['NTrueTracks_3hits'].GetXaxis().SetBinLabel(3,"Stations 3&4")
226
227 ut.bookHist(h,'Reco_y12_3hits','Number of recognized tracks, clones and ghosts with more than 3 hits in stations 1&2, Y views', 5, -0.5, 4.5)
228 h['Reco_y12_3hits'].GetXaxis().SetBinLabel(1,"N total")
229 h['Reco_y12_3hits'].GetXaxis().SetBinLabel(2,"N recognized tracks")
230 h['Reco_y12_3hits'].GetXaxis().SetBinLabel(3,"N clones")
231 h['Reco_y12_3hits'].GetXaxis().SetBinLabel(4,"N ghosts")
232 h['Reco_y12_3hits'].GetXaxis().SetBinLabel(5,"N others")
233
234 ut.bookHist(h,'Reco_stereo12_3hits','Number of recognized tracks, clones and ghosts with more than 3 hits in stations 1&2, Stereo views', 5, -0.5, 4.5)
235 h['Reco_stereo12_3hits'].GetXaxis().SetBinLabel(1,"N total")
236 h['Reco_stereo12_3hits'].GetXaxis().SetBinLabel(2,"N recognized tracks")
237 h['Reco_stereo12_3hits'].GetXaxis().SetBinLabel(3,"N clones")
238 h['Reco_stereo12_3hits'].GetXaxis().SetBinLabel(4,"N ghosts")
239 h['Reco_stereo12_3hits'].GetXaxis().SetBinLabel(5,"N others")
240
241 ut.bookHist(h,'Reco_34_3hits','Number of recognized tracks, clones and ghosts with more than 3 hits in stations 3&4', 5, -0.5, 4.5)
242 h['Reco_34_3hits'].GetXaxis().SetBinLabel(1,"N total")
243 h['Reco_34_3hits'].GetXaxis().SetBinLabel(2,"N recognized tracks")
244 h['Reco_34_3hits'].GetXaxis().SetBinLabel(3,"N clones")
245 h['Reco_34_3hits'].GetXaxis().SetBinLabel(4,"N ghosts")
246 h['Reco_34_3hits'].GetXaxis().SetBinLabel(5,"N others")
247
248
249
250 ut.bookHist(h,'NTrueTracks_Tr4','Number of tracks. At least one hit in stations 1-4.', 3, -0.5, 2.5)
251 h['NTrueTracks_Tr4'].GetXaxis().SetBinLabel(1,"Stations 1&2, Y views")
252 h['NTrueTracks_Tr4'].GetXaxis().SetBinLabel(2,"Stations 1&2, Stereo views")
253 h['NTrueTracks_Tr4'].GetXaxis().SetBinLabel(3,"Stations 3&4")
254
255 ut.bookHist(h,'Reco_y12_Tr4','Number of recognized tracks, clones and ghosts in stations 1&2, Y views. At least one hit in stations 1-4.', 5, -0.5, 4.5)
256 h['Reco_y12_Tr4'].GetXaxis().SetBinLabel(1,"N total")
257 h['Reco_y12_Tr4'].GetXaxis().SetBinLabel(2,"N recognized tracks")
258 h['Reco_y12_Tr4'].GetXaxis().SetBinLabel(3,"N clones")
259 h['Reco_y12_Tr4'].GetXaxis().SetBinLabel(4,"N ghosts")
260 h['Reco_y12_Tr4'].GetXaxis().SetBinLabel(5,"N others")
261
262 ut.bookHist(h,'Reco_stereo12_Tr4','Number of recognized tracks, clones and ghosts in stations 1&2, Stereo views. At least one hit in stations 1-4.', 5, -0.5, 4.5)
263 h['Reco_stereo12_Tr4'].GetXaxis().SetBinLabel(1,"N total")
264 h['Reco_stereo12_Tr4'].GetXaxis().SetBinLabel(2,"N recognized tracks")
265 h['Reco_stereo12_Tr4'].GetXaxis().SetBinLabel(3,"N clones")
266 h['Reco_stereo12_Tr4'].GetXaxis().SetBinLabel(4,"N ghosts")
267 h['Reco_stereo12_Tr4'].GetXaxis().SetBinLabel(5,"N others")
268
269 ut.bookHist(h,'Reco_34_Tr4','Number of recognized tracks, clones and ghosts in stations 3&4. At least one hit in stations 1-4.', 5, -0.5, 4.5)
270 h['Reco_34_Tr4'].GetXaxis().SetBinLabel(1,"N total")
271 h['Reco_34_Tr4'].GetXaxis().SetBinLabel(2,"N recognized tracks")
272 h['Reco_34_Tr4'].GetXaxis().SetBinLabel(3,"N clones")
273 h['Reco_34_Tr4'].GetXaxis().SetBinLabel(4,"N ghosts")
274 h['Reco_34_Tr4'].GetXaxis().SetBinLabel(5,"N others")
275
276 ut.bookHist(h,'Reco_target','Number of recognized target tracks, clones and ghosts.', 5, -0.5, 4.5)
277 h['Reco_target'].GetXaxis().SetBinLabel(1,"N total")
278 h['Reco_target'].GetXaxis().SetBinLabel(2,"N recognized tracks")
279 h['Reco_target'].GetXaxis().SetBinLabel(3,"N clones")
280 h['Reco_target'].GetXaxis().SetBinLabel(4,"N ghosts")
281 h['Reco_target'].GetXaxis().SetBinLabel(5,"N others")
282
283 ut.bookHist(h,'Reco_muon','Number of recognized muon tracks, clones and ghosts.', 5, -0.5, 4.5)
284 h['Reco_muon'].GetXaxis().SetBinLabel(1,"N total")
285 h['Reco_muon'].GetXaxis().SetBinLabel(2,"N recognized tracks")
286 h['Reco_muon'].GetXaxis().SetBinLabel(3,"N clones")
287 h['Reco_muon'].GetXaxis().SetBinLabel(4,"N ghosts")
288 h['Reco_muon'].GetXaxis().SetBinLabel(5,"N others")
289
290 ut.bookHist(h,'Reco_muon_with_tagger','Number of recognized muon tracks, clones and ghosts with tagger hits.', 5, -0.5, 4.5)
291 h['Reco_muon_with_tagger'].GetXaxis().SetBinLabel(1,"N total")
292 h['Reco_muon_with_tagger'].GetXaxis().SetBinLabel(2,"N recognized tracks")
293 h['Reco_muon_with_tagger'].GetXaxis().SetBinLabel(3,"N clones")
294 h['Reco_muon_with_tagger'].GetXaxis().SetBinLabel(4,"N ghosts")
295 h['Reco_muon_with_tagger'].GetXaxis().SetBinLabel(5,"N others")
296
297 ut.bookHist(h,'Reco_all_tracks','Number of recognized all tracks, clones and ghosts.', 5, -0.5, 4.5)
298 h['Reco_all_tracks'].GetXaxis().SetBinLabel(1,"N total")
299 h['Reco_all_tracks'].GetXaxis().SetBinLabel(2,"N recognized tracks")
300 h['Reco_all_tracks'].GetXaxis().SetBinLabel(3,"N clones")
301 h['Reco_all_tracks'].GetXaxis().SetBinLabel(4,"N ghosts")
302 h['Reco_all_tracks'].GetXaxis().SetBinLabel(5,"N others")
303
304
305 ut.bookHist(h,'NHits_true_y12','Number of hits per MC track. Stations 1&2, Y views', 9, -0.5, 8.5)
306 ut.bookHist(h,'NHits_true_stereo12','Number of hits per MC track. Stations 1&2, Stereo views', 9, -0.5, 8.5)
307 ut.bookHist(h,'NHits_true_34','Number of hits per MC track. Stations 3&4, Y views', 9, -0.5, 8.5)
308
309 ut.bookHist(h,'NHits_reco_y12','Number of hits per Reco track. Stations 1&2, Y views', 9, -0.5, 8.5)
310 ut.bookHist(h,'NHits_reco_stereo12','Number of hits per Reco track. Stations 1&2, Stereo views', 9, -0.5, 8.5)
311 ut.bookHist(h,'NHits_reco_34','Number of hits per Reco track. Stations 3&4, Y views', 9, -0.5, 8.5)
312
313
314 ut.bookHist(h,'p/pt','P vs Pt (GeV), Reco',100,0.,400.,100,0.,10.)
315 ut.bookHist(h,'p/pt_truth','P vs Pt (GeV), MC Truth',100,0.,400.,100,0.,10.)
316 ut.bookHist(h,'p/pt_noT4','P vs Pt (GeV), Reco',100,0.,400.,100,0.,10.)
317 ut.bookHist(h,'p/pt_truth_noT4','P vs Pt (GeV), MC Truth',100,0.,400.,100,0.,10.)
318 ut.bookHist(h,'p/pt_all','P vs Pt (GeV), Reco',100,0.,400.,100,0.,10.)
319 ut.bookHist(h,'p/pt_truth_all','P vs Pt (GeV), MC Truth',100,0.,400.,100,0.,10.)
320
321 ut.bookHist(h,'p_rel_error','(P_reco - P_true) / P_true',200,-2.,2.)
322 ut.bookHist(h,'pt_rel_error','(Pt_reco - Pt_true) / Pt_true',200,-2.,2.)
323 ut.bookHist(h,'p_rel_error_noT4','(P_reco - P_true) / P_true',200,-2.,2.)
324 ut.bookHist(h,'pt_rel_error_noT4','(Pt_reco - Pt_true) / Pt_true',200,-2.,2.)
325 ut.bookHist(h,'p_rel_error_all','(P_reco - P_true) / P_true',200,-2.,2.)
326 ut.bookHist(h,'pt_rel_error_all','(Pt_reco - Pt_true) / Pt_true',200,-2.,2.)
327
328
329 ut.bookHist(h,'Reco_muons_vs_p_true','Number of recognized muons vs P MC truth',40,0.,400.)
330 ut.bookHist(h,'Ghosts_muons_vs_p_true','Number of ghosts vs P MC truth',40,0.,400.)
331 ut.bookHist(h,'True_muons_vs_p_true','Number of muons vs P MC truth',40,0.,400.)
332 ut.bookHist(h,'True_all_tracks_vs_p_true','Number of muons vs P MC truth',40,0.,400.)
333
334# -----Create geometry----------------------------------------------
335import charmDet_conf
336run = ROOT.FairRunSim()
337run.SetName("TGeant4") # Transport engine
338run.SetOutputFile(ROOT.TMemFile('output', 'recreate')) # Output file
339run.SetUserConfig("g4Config_basic.C") # geant4 transport not used, only needed for creating VMC field
340rtdb = run.GetRuntimeDb()
341modules = charmDet_conf.configure(run,ShipGeo)
342# -----Create geometry----------------------------------------------
343fgeo.FAIRGeom
344
345# make global variables
346global_variables.debug = debug
347global_variables.withT0 = withT0
348global_variables.realPR = realPR
349
350global_variables.ShipGeo = ShipGeo
351global_variables.modules = modules
352
353global_variables.withNoStrawSmearing = withNoStrawSmearing
354global_variables.withNTaggerHits = withNTaggerHits
355global_variables.withDist2Wire = withDist2Wire
356global_variables.h = h
357global_variables.log = log
358global_variables.simulation = simulation
359global_variables.iEvent = 0
360
361# import reco tasks
362import MufluxDigiReco
363SHiP = MufluxDigiReco.MufluxDigiReco(outFile)
364
365nEvents = min(SHiP.sTree.GetEntries(),nEvents)
366# main loop
367for global_variables.iEvent in range(firstEvent, nEvents):
368 if global_variables.iEvent % 1000 == 0 or global_variables.debug:
369 print('event ', global_variables.iEvent)
370 SHiP.iEvent = global_variables.iEvent
371 rc = SHiP.sTree.GetEvent(global_variables.iEvent)
372 if global_variables.simulation:
373 SHiP.digitize()
374 # IS BROKEN SHiP.reconstruct()
375 # memory monitoring
376 # mem_monitor()
377
378# end loop over events
379SHiP.finish()
charmDet_conf.configure
configure(run, ship_geo, Gfield='')
Definition charmDet_conf.py:22
shipRoot_conf.configure
configure(darkphoton=None)
Definition shipRoot_conf.py:67