ATLAS Offline Software
EventAnalyzer.py
Go to the documentation of this file.
1 # Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
2 
3 import ROOT
4 import HistControl
5 
6 class 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
17  self.__dark_photon_index = None
18  self.__scalar_index = None
19  histControlName = histControlName.replace('.root', '')
20  self.__histControl = HistControl.HistControl( histControlName, doOutputFile=True )
21 
23  self.__doPtBalancePlots = True
24  self.__doPtEtaPhiPlots = True
26  self.__doRZPlots = True
27 
28  maxPt = 100.e3
29 
30  if self.__doPolarizationPlots:
31  self.__histControl.BookHist1D('Polarization', 'cos(#theta)', 200, -1.1, 1.1)
32  pass
33 
34  if self.__doOpeningAnglesPlots:
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:
107  self.__dark_photon_child_index_1 = -999
108  self.__dark_photon_child_index_2 = -999
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 
164  if self.__doOpeningAnglesPlots:
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
EventAnalyzer.EventAnalyzer.fillOpeningAngles
def fillOpeningAngles(self, index1, index2)
Definition: EventAnalyzer.py:228
EventAnalyzer.EventAnalyzer.__init__
def __init__(self, name="EventAnalyzer", title="EventAnalyzer", histControlName='histControl', numberOfDarkPhotons=1)
Definition: EventAnalyzer.py:8
EventAnalyzer.EventAnalyzer.__nDarkPhotons
__nDarkPhotons
Definition: EventAnalyzer.py:13
max
#define max(a, b)
Definition: cfImp.cxx:41
EventAnalyzer.EventAnalyzer.cosTheta
def cosTheta(self)
Definition: EventAnalyzer.py:216
EventAnalyzer.EventAnalyzer.deltaPhi
def deltaPhi(self, phi1, phi2)
Definition: EventAnalyzer.py:87
EventAnalyzer.EventAnalyzer.__doPtEtaPhiPlots
__doPtEtaPhiPlots
Definition: EventAnalyzer.py:24
EventAnalyzer.EventAnalyzer.__event
__event
Definition: EventAnalyzer.py:14
EventAnalyzer.EventAnalyzer.etaPhi
def etaPhi(self)
Definition: EventAnalyzer.py:350
dumpHVPathFromNtuple.append
bool append
Definition: dumpHVPathFromNtuple.py:91
HistControl.HistControl
Definition: HistControl.py:9
EventAnalyzer.EventAnalyzer.fillPolarization
def fillPolarization(self, indexDP, indexKid)
Definition: EventAnalyzer.py:182
EventAnalyzer.EventAnalyzer.__doPolarizationPlots
__doPolarizationPlots
Definition: EventAnalyzer.py:22
EventAnalyzer.EventAnalyzer.openingAngles
def openingAngles(self)
Definition: EventAnalyzer.py:279
EventAnalyzer.EventAnalyzer.finalize
def finalize(self)
Definition: EventAnalyzer.py:102
EventAnalyzer.EventAnalyzer.__doRZPlots
__doRZPlots
Definition: EventAnalyzer.py:26
EventAnalyzer.EventAnalyzer.__histControl
__histControl
Definition: EventAnalyzer.py:20
calibdata.exit
exit
Definition: calibdata.py:236
EventAnalyzer.EventAnalyzer.setParticleIndices
def setParticleIndices(self)
Definition: EventAnalyzer.py:105
min
#define min(a, b)
Definition: cfImp.cxx:40
EventAnalyzer.EventAnalyzer.__dark_photon_child_index_1
__dark_photon_child_index_1
Definition: EventAnalyzer.py:15
EventAnalyzer.EventAnalyzer.__dark_photon_index
__dark_photon_index
Definition: EventAnalyzer.py:17
EventAnalyzer.EventAnalyzer.ptBalance
def ptBalance(self)
Definition: EventAnalyzer.py:318
EventAnalyzer.EventAnalyzer.__doOpeningAnglesPlots
__doOpeningAnglesPlots
Definition: EventAnalyzer.py:25
EventAnalyzer.EventAnalyzer.deltaR
def deltaR(self, eta1, phi1, eta2, phi2)
Definition: EventAnalyzer.py:95
EventAnalyzer.EventAnalyzer.fillPtBalance
def fillPtBalance(self, index1, index2)
Definition: EventAnalyzer.py:292
EventAnalyzer.EventAnalyzer.fillRZ
def fillRZ(self, indexDP, indexKid)
Definition: EventAnalyzer.py:330
EventAnalyzer.EventAnalyzer.__scalar_index
__scalar_index
Definition: EventAnalyzer.py:18
EventAnalyzer.EventAnalyzer.RZ
def RZ(self)
Definition: EventAnalyzer.py:339
EventAnalyzer.EventAnalyzer
Definition: EventAnalyzer.py:6
EventAnalyzer.EventAnalyzer.__doPtBalancePlots
__doPtBalancePlots
Definition: EventAnalyzer.py:23
EventAnalyzer.EventAnalyzer.__dark_photon_child_index_2
__dark_photon_child_index_2
Definition: EventAnalyzer.py:16
EventAnalyzer.EventAnalyzer.processEvent
def processEvent(self, event)
Definition: EventAnalyzer.py:153