ATLAS Offline Software
Loading...
Searching...
No Matches
EventAnalyzer.py
Go to the documentation of this file.
1# Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
2
3import ROOT
4import HistControl
5
6class EventAnalyzer( ROOT.TNamed ):
7
8 def __init__( self, name = "EventAnalyzer", title = "EventAnalyzer", histControlName='histControl', numberOfDarkPhotons=1 ):
9 # Initialise the base class
10 ROOT.TNamed.__init__( self, name, title )
11 self.SetName( name )
12 self.SetTitle( title )
13 self.__nDarkPhotons = numberOfDarkPhotons
14 self.__event = None
18 self.__scalar_index = None
19 histControlName = histControlName.replace('.root', '')
20 self.__histControl = HistControl.HistControl( histControlName, doOutputFile=True )
21
26 self.__doRZPlots = True
27
28 maxPt = 100.e3
29
31 self.__histControl.BookHist1D('Polarization', 'cos(#theta)', 200, -1.1, 1.1)
32 pass
33
35 self.__histControl.BookHist1D('OpeningEtaElectrons', '#Delta #eta', 200, .0, .1)
36 self.__histControl.BookHist1D('OpeningPhiElectrons', '#Delta #phi', 200, .0, .1)
37 self.__histControl.BookHist1D('OpeningRElectrons', '#Delta R', 200, .0, .1)
38
39 self.__histControl.BookHist1D('OpeningEtaMuons', '#Delta #eta', 200, .0, .1)
40 self.__histControl.BookHist1D('OpeningPhiMuons', '#Delta #phi', 200, .0, .1)
41 self.__histControl.BookHist1D('OpeningRMuons', '#Delta R', 200, .0, .1)
42
43 self.__histControl.BookHist1D('OpeningEtaPions', '#Delta #eta', 200, .0, .1)
44 self.__histControl.BookHist1D('OpeningPhiPions', '#Delta #phi', 200, .0, .1)
45 self.__histControl.BookHist1D('OpeningRPions', '#Delta R', 200, .0, .1)
46
47 self.__histControl.BookHist1D('OpeningEtaDarkPhotons', '#Delta #eta', 200, .0, 4.)
48 self.__histControl.BookHist1D('OpeningPhiDarkPhotons', '#Delta #phi', 200, .0, 4.)
49 self.__histControl.BookHist1D('OpeningRDarkPhotons', '#Delta R', 200, .0, 4.)
50 pass
51
52 if self.__doPtBalancePlots:
53 self.__histControl.BookHist1D('leadingPtElectron', 'Leading p_{T};p_{T} [GeV];Entries', 200, .0, maxPt)
54 self.__histControl.BookHist1D('subLeadingPtElectron', 'Sub-Leading p_{T};p_{T} [GeV];Entries', 200, .0, maxPt/2.)
55 self.__histControl.BookHist1D('ptBalanceElectron', 'p_{T} Balance', 200, -.1, 1.1)
56
57 self.__histControl.BookHist1D('leadingPtMuon', 'Leading p_{T};p_{T} [GeV];Entries', 200, .0, maxPt)
58 self.__histControl.BookHist1D('subLeadingPtMuon', 'Sub-Leading p_{T};p_{T} [GeV];Entries', 200, .0, maxPt/2.)
59 self.__histControl.BookHist1D('ptBalanceMuon', 'p_{T} Balance', 200, -.1, 1.1)
60
61 self.__histControl.BookHist1D('leadingPtPion', 'Leading p_{T};p_{T} [GeV];Entries', 200, .0, maxPt)
62 self.__histControl.BookHist1D('subLeadingPtPion', 'Sub-Leading p_{T};p_{T} [GeV];Entries', 200, .0, maxPt/2.)
63 self.__histControl.BookHist1D('ptBalancePion', 'p_{T} Balance', 200, -.1, 1.1)
64
65 self.__histControl.BookHist1D('leadingPtDarkPhoton', 'Leading p_{T};p_{T} [GeV];Entries', 200, .0, maxPt)
66 self.__histControl.BookHist1D('subLeadingPtDarkPhoton', 'Sub-Leading p_{T};p_{T} [GeV];Entries', 200, .0, maxPt/2.)
67 self.__histControl.BookHist1D('ptBalanceDarkPhoton', 'p_{T} Balance', 200, -.1, 1.1)
68 pass
69
70 if self.__doRZPlots:
71 self.__histControl.BookHist1D('decayR1D', 'Dark Photon Decay Radius;R [mm];Entries', 200, 1., -1.)
72 self.__histControl.BookHist1D('decayZ1D', 'Dark Photon Decay Length along Z;Z [mm];Entries', 200, 1., -1.)
73 self.__histControl.BookHist2D('decayZR2D', 'Dark Photon Decay Distance;Z [mm];R [mm]', 150, 1., -1., 150, 1., -1.)
74 self.__histControl.BookHist2D('decayEtaR2D', 'Dark Photon Decay Distance;#eta;R [mm]', 150, 1., -1., 150, 1., -1.)
75 pass
76
77 if self.__doPtEtaPhiPlots:
78 self.__histControl.BookHist1D('darkPhotonEta', 'Dark Photon #eta', 200, 1., -1.)
79 self.__histControl.BookHist1D('darkPhotonPhi', 'Dark Photon #phi', 200, 1., -1.)
80 self.__histControl.BookHist1D('darkPhotonPt', 'Dark Photon p_{T}', 200, 1., -1.)
81 self.__histControl.BookHist1D('scalarEta', 'Scalar #eta', 200, 1., -1.)
82 self.__histControl.BookHist1D('scalarPhi', 'Scalar #phi', 200, 1., -1.)
83 self.__histControl.BookHist1D('scalarPt', 'Scalar p_{T}', 200, 1., -1.)
84 pass
85 return
86
87 def deltaPhi(self, phi1, phi2):
88 dphi = phi1 - phi2
89 while dphi < -ROOT.TMath.Pi():
90 dphi = dphi + 2. * ROOT.TMath.Pi()
91 while dphi > ROOT.TMath.Pi():
92 dphi = dphi - 2. * ROOT.TMath.Pi()
93 return dphi;
94
95 def deltaR(self, eta1, phi1, eta2, phi2):
96 dphi = self.deltaPhi(phi1, phi2)
97 deta = ROOT.TMath.Abs(eta1-eta2)
98 dr = ROOT.TMath.Sqrt(dphi*dphi + deta*deta)
99 return dr
100
101
102 def finalize(self):
103 self.__histControl.SaveAllToFile()
104
106 if self.__nDarkPhotons == 1:
109 self.__dark_photon_index = -999
110 pass
111 elif self.__nDarkPhotons == 2:
114 self.__dark_photon_index = []
115 self.__scalar_index = -999
116 pass
117 else:
118 print 'Cannot have', self.__nDarkPhotons, 'per event. Set to 1 or 2'
119 exit(0)
120 pass
121
122 for index,(mc_child_index,pdgId) in enumerate( zip( self.__event.mc_child_index, self.__event.mc_pdgId ) ):
123 if pdgId==700022:
124 if self.__nDarkPhotons == 1:
125 self.__dark_photon_index = index
126 else:
127 self.__dark_photon_index.append(index)
128 pass
129 if len(mc_child_index) != 2:
130 print 'dark photon has ', len(mc_child_index), 'children instead of 2'
131 return False
132 if self.__nDarkPhotons == 1:
133 self.__dark_photon_child_index_1 = mc_child_index[0]
134 self.__dark_photon_child_index_2 = mc_child_index[1]
135 else:
136 self.__dark_photon_child_index_1.append(mc_child_index[0])
137 self.__dark_photon_child_index_2.append(mc_child_index[1])
138 pass
139 pass
140 if pdgId==700021:
141 if self.__nDarkPhotons == 1:
142 print 'Found a scalar (pdgId==700021). But number of dark photons is set to 1. fix this'
143 exit(0)
144 pass
145 self.__scalar_index = index
146 pass
147 if self.__scalar_index == -999 and self.__nDarkPhotons==2:
148 print 'Didn\'t find a scalar (pdgId==700021). But number of dark photons is set to 2. fix this'
149 exit(0)
150 pass
151 return True
152
153 def processEvent(self, event):
154
155 self.__event = event
156
157 if not self.setParticleIndices():
158 return False
159
160 if self.__doPolarizationPlots:
161 self.cosTheta()
162 pass
163
165 self.openingAngles()
166 pass
167
168 if self.__doPtBalancePlots:
169 self.ptBalance()
170 pass
171
172 if self.__doRZPlots:
173 self.RZ()
174 pass
175
176 if self.__doPtEtaPhiPlots:
177 self.etaPhi()
178 pass
179
180 return True
181
182 def fillPolarization(self, indexDP, indexKid):
183 '''Get cos(theta) between the dark photon vector in the lab frame
184 and the decay product vector in the dark photon\'s rest frame.
185 '''
186 eta = self.__event.mc_eta[indexDP]
187 phi = self.__event.mc_phi[indexDP]
188 pt = self.__event.mc_pt[indexDP]
189
190 theta = 2.*ROOT.TMath.ATan(ROOT.TMath.Exp(-1.*eta))
191 px = pt*ROOT.TMath.Sin(phi)
192 py = pt*ROOT.TMath.Cos(phi)
193 pz = pt/ROOT.TMath.Tan(theta)
194 edp = ROOT.TMath.Sqrt( px**2+py**2+pz**2+self.__event.mc_m[indexDP]**2 )
195
196 eta_1 = self.__event.mc_eta[indexKid]
197 phi_1 = self.__event.mc_phi[indexKid]
198 pt_1 = self.__event.mc_pt[indexKid]
199 theta_1 = 2.*ROOT.TMath.ATan(ROOT.TMath.Exp(-1.*eta_1))
200
201 px_kid_1 = pt_1*ROOT.TMath.Sin(phi_1)
202 py_kid_1 = pt_1*ROOT.TMath.Cos(phi_1)
203 pz_kid_1 = pt_1/ROOT.TMath.Tan(theta_1)
204 e_kid_1 = ROOT.TMath.Sqrt( px_kid_1**2+py_kid_1**2+pz_kid_1**2+self.__event.mc_m[indexKid]**2 )
205
206 tlv_dp = ROOT.TLorentzVector(px, py, pz, edp)
207 bv = tlv_dp.BoostVector()
208
209 tlv_kid_1 = ROOT.TLorentzVector(px_kid_1, py_kid_1, pz_kid_1, e_kid_1)
210
211 tlv_kid_1.Boost(-bv)
212
213 self.__histControl.FillHist1D('Polarization', tlv_dp.Vect().Unit().Dot(tlv_kid_1.Vect().Unit()) )
214 return
215
216 def cosTheta(self):
217 if self.__nDarkPhotons == 1:
219 pass
220 else:
221 for (indexDP,indexKid) in zip (self.__dark_photon_index,self.__dark_photon_child_index_1):
222 self.fillPolarization(indexDP, indexKid)
223 continue
224 pass
225 pass
226 return
227
228 def fillOpeningAngles(self,index1, index2):
229 eta_1 = self.__event.mc_eta[index1]
230 phi_1 = self.__event.mc_phi[index1]
231 pt_1 = self.__event.mc_pt[index1]
232 theta_1 = 2.*ROOT.TMath.ATan(ROOT.TMath.Exp(-1.*eta_1))
233
234 px_kid_1 = pt_1*ROOT.TMath.Sin(phi_1)
235 py_kid_1 = pt_1*ROOT.TMath.Cos(phi_1)
236 pz_kid_1 = pt_1/ROOT.TMath.Tan(theta_1)
237 e_kid_1 = ROOT.TMath.Sqrt( px_kid_1**2+py_kid_1**2+pz_kid_1**2+self.__event.mc_m[index1]**2 )
238
239 eta_2 = self.__event.mc_eta[index2]
240 phi_2 = self.__event.mc_phi[index2]
241 pt_2 = self.__event.mc_pt[index2]
242 theta_2 = 2.*ROOT.TMath.ATan(ROOT.TMath.Exp(-1.*eta_2))
243
244 px_kid_2 = pt_2*ROOT.TMath.Sin(phi_2)
245 py_kid_2 = pt_2*ROOT.TMath.Cos(phi_2)
246 pz_kid_2 = pt_2/ROOT.TMath.Tan(theta_2)
247 e_kid_2 = ROOT.TMath.Sqrt( px_kid_2**2+py_kid_2**2+pz_kid_2**2+self.__event.mc_m[index2]**2 )
248
249 tlv_kid_1 = ROOT.TLorentzVector(px_kid_1, py_kid_1, pz_kid_1, e_kid_1)
250 tlv_kid_2 = ROOT.TLorentzVector(px_kid_2, py_kid_2, pz_kid_2, e_kid_2)
251
252 dPhi = abs( self.deltaPhi( tlv_kid_1.Phi(), tlv_kid_2.Phi() ) )
253 dEta = abs( tlv_kid_1.Eta() - tlv_kid_2.Eta() )
254 dR = self.deltaR( tlv_kid_1.Eta(), tlv_kid_1.Phi(), tlv_kid_2.Eta(), tlv_kid_2.Phi() )
255
256 pdgId = abs(self.__event.mc_pdgId[index1])
257
258 if pdgId == 11:
259 self.__histControl.FillHist1D('OpeningEtaElectrons', dEta )
260 self.__histControl.FillHist1D('OpeningPhiElectrons', dPhi )
261 self.__histControl.FillHist1D('OpeningRElectrons', dR )
262 pass
263 if pdgId == 13:
264 self.__histControl.FillHist1D('OpeningEtaMuons', dEta )
265 self.__histControl.FillHist1D('OpeningPhiMuons', dPhi )
266 self.__histControl.FillHist1D('OpeningRMuons', dR )
267 pass
268 if pdgId == 211:
269 self.__histControl.FillHist1D('OpeningEtaPions', dEta )
270 self.__histControl.FillHist1D('OpeningPhiPions', dPhi )
271 self.__histControl.FillHist1D('OpeningRPions', dR )
272 pass
273 if pdgId == 700022:
274 self.__histControl.FillHist1D('OpeningEtaDarkPhotons', dEta )
275 self.__histControl.FillHist1D('OpeningPhiDarkPhotons', dPhi )
276 self.__histControl.FillHist1D('OpeningRDarkPhotons', dR )
277 pass
278
279 def openingAngles(self):
280 if self.__nDarkPhotons == 1:
282 pass
283 else:
284 for (index1,index2) in zip (self.__dark_photon_child_index_1,self.__dark_photon_child_index_2):
285 self.fillOpeningAngles(index1, index2)
286 pass
288 pass
289
290 return
291
292 def fillPtBalance(self, index1, index2):
293 pt_1 = self.__event.mc_pt[index1]
294 pt_2 = self.__event.mc_pt[index2]
295 pdgId = abs(self.__event.mc_pdgId[index1])
296 if pdgId == 11:
297 self.__histControl.FillHist1D('leadingPtElectron', max(pt_1, pt_2) )
298 self.__histControl.FillHist1D('subLeadingPtElectron', min(pt_1, pt_2) )
299 self.__histControl.FillHist1D('ptBalanceElectron', abs(pt_1 - pt_2)/(pt_1+pt_2) )
300 pass
301 if pdgId == 13:
302 self.__histControl.FillHist1D('leadingPtMuon', max(pt_1, pt_2) )
303 self.__histControl.FillHist1D('subLeadingPtMuon', min(pt_1, pt_2) )
304 self.__histControl.FillHist1D('ptBalanceMuon', abs(pt_1 - pt_2)/(pt_1+pt_2) )
305 pass
306 if pdgId == 211:
307 self.__histControl.FillHist1D('leadingPtPion', max(pt_1, pt_2) )
308 self.__histControl.FillHist1D('subLeadingPtPion', min(pt_1, pt_2) )
309 self.__histControl.FillHist1D('ptBalancePion', abs(pt_1 - pt_2)/(pt_1+pt_2) )
310 pass
311 if pdgId == 700022:
312 self.__histControl.FillHist1D('leadingPtDarkPhoton', max(pt_1, pt_2) )
313 self.__histControl.FillHist1D('subLeadingPtDarkPhoton', min(pt_1, pt_2) )
314 self.__histControl.FillHist1D('ptBalanceDarkPhoton', abs(pt_1 - pt_2)/(pt_1+pt_2) )
315 pass
316 return
317
318 def ptBalance(self):
319 if self.__nDarkPhotons == 1:
321 pass
322 else:
323 for (index1,index2) in zip (self.__dark_photon_child_index_1,self.__dark_photon_child_index_2):
324 self.fillPtBalance(index1, index2)
325 pass
327 pass
328 return
329
330 def fillRZ(self, indexDP, indexKid):
331 decayR = ROOT.TMath.Sqrt( self.__event.mc_vx_x[indexKid]**2 + self.__event.mc_vx_y[indexKid]**2 )
332 decayZ = self.__event.mc_vx_z[indexKid]
333 self.__histControl.FillHist1D('decayR1D', decayR)
334 self.__histControl.FillHist1D('decayZ1D', decayZ)
335 self.__histControl.FillHist2D('decayZR2D', decayZ, decayR)
336 self.__histControl.FillHist2D('decayEtaR2D', self.__event.mc_eta[indexDP], decayR)
337 return
338
339 def RZ(self):
340 if self.__nDarkPhotons == 1:
342 pass
343 else:
344 for (indexDP,indexKid) in zip (self.__dark_photon_index,self.__dark_photon_child_index_1):
345 self.fillRZ(indexDP, indexKid)
346 pass
347 pass
348 return
349
350 def etaPhi(self):
351 if self.__nDarkPhotons == 1:
352 self.__histControl.FillHist1D('darkPhotonEta', self.__event.mc_eta[self.__dark_photon_index])
353 self.__histControl.FillHist1D('darkPhotonPhi', self.__event.mc_phi[self.__dark_photon_index])
354 self.__histControl.FillHist1D('darkPhotonPt', self.__event.mc_pt[self.__dark_photon_index] )
355 else:
356 for index in self.__dark_photon_index:
357 self.__histControl.FillHist1D('darkPhotonEta', self.__event.mc_eta[index])
358 self.__histControl.FillHist1D('darkPhotonPhi', self.__event.mc_phi[index])
359 self.__histControl.FillHist1D('darkPhotonPt', self.__event.mc_pt[index] )
360 pass
361 self.__histControl.FillHist1D('scalarEta', self.__event.mc_eta[self.__scalar_index])
362 self.__histControl.FillHist1D('scalarPhi', self.__event.mc_phi[self.__scalar_index])
363 self.__histControl.FillHist1D('scalarPt', self.__event.mc_pt[self.__scalar_index] )
364 pass
365 return
#define min(a, b)
Definition cfImp.cxx:40
#define max(a, b)
Definition cfImp.cxx:41
fillRZ(self, indexDP, indexKid)
__init__(self, name="EventAnalyzer", title="EventAnalyzer", histControlName='histControl', numberOfDarkPhotons=1)
fillPtBalance(self, index1, index2)
deltaR(self, eta1, phi1, eta2, phi2)
fillOpeningAngles(self, index1, index2)
deltaPhi(self, phi1, phi2)
fillPolarization(self, indexDP, indexKid)