Belle II Software  release-05-01-25
PhokharaEvtgenAnalyze.py
1 #!/usr/bin/env python3
2 # -*- coding: utf-8 -*-
3 
4 """
5 <header>
6  <input>PhokharaEvtgenData.root</input>
7  <output>PhokharaEvtgenAnalysis.root</output>
8  <contact>Kirill Chilikin (chilikin@lebedev.ru)</contact>
9  <description>Analysis of e+ e- -> J/psi eta_c events.</description>
10 </header>
11 """
12 
13 from basf2 import *
14 import ROOT
15 import numpy
16 from ROOT import Belle2
17 
18 
20  """ Analysis module for PhokharaEvtgen. """
21 
22  def __init__(self):
23  """Initialization."""
24  super(PhokharaEvtgenAnalysisModule, self).__init__()
25 
26  self.output_file = ROOT.TFile('PhokharaEvtgenAnalysis.root', 'recreate')
27 
28  self.tree = ROOT.TTree('tree', '')
29 
30  self.ecms = numpy.zeros(1, dtype=numpy.float32)
31 
32  self.gamma_e = numpy.zeros(1, dtype=numpy.float32)
33 
34  self.gamma_px = numpy.zeros(1, dtype=numpy.float32)
35 
36  self.gamma_py = numpy.zeros(1, dtype=numpy.float32)
37 
38  self.gamma_pz = numpy.zeros(1, dtype=numpy.float32)
39 
40  self.jpsi_e = numpy.zeros(1, dtype=numpy.float32)
41 
42  self.jpsi_px = numpy.zeros(1, dtype=numpy.float32)
43 
44  self.jpsi_py = numpy.zeros(1, dtype=numpy.float32)
45 
46  self.jpsi_pz = numpy.zeros(1, dtype=numpy.float32)
47 
48  self.lepton_e = numpy.zeros(1, dtype=numpy.float32)
49 
50  self.lepton_px = numpy.zeros(1, dtype=numpy.float32)
51 
52  self.lepton_py = numpy.zeros(1, dtype=numpy.float32)
53 
54  self.lepton_pz = numpy.zeros(1, dtype=numpy.float32)
55  self.tree.Branch('ecms', self.ecms, 'ecms/F')
56  self.tree.Branch('gamma_e', self.gamma_e, 'gamma_e/F')
57  self.tree.Branch('gamma_px', self.gamma_px, 'gamma_px/F')
58  self.tree.Branch('gamma_py', self.gamma_py, 'gamma_py/F')
59  self.tree.Branch('gamma_pz', self.gamma_pz, 'gamma_pz/F')
60  self.tree.Branch('jpsi_e', self.jpsi_e, 'jpsi_e/F')
61  self.tree.Branch('jpsi_px', self.jpsi_px, 'jpsi_px/F')
62  self.tree.Branch('jpsi_py', self.jpsi_py, 'jpsi_py/F')
63  self.tree.Branch('jpsi_pz', self.jpsi_pz, 'jpsi_pz/F')
64  self.tree.Branch('lepton_e', self.lepton_e, 'lepton_e/F')
65  self.tree.Branch('lepton_px', self.lepton_px, 'lepton_px/F')
66  self.tree.Branch('lepton_py', self.lepton_py, 'lepton_py/F')
67  self.tree.Branch('lepton_pz', self.lepton_pz, 'lepton_pz/F')
68 
69  def event(self):
70  """ Event function. """
71  mc_initial_particles = Belle2.PyStoreObj('MCInitialParticles')
72  self.ecms[0] = mc_initial_particles.getMass()
73  mc_particles = Belle2.PyStoreArray('MCParticles')
74  for mc_particle in mc_particles:
75  # Select virtual photons.
76  if (mc_particle.getPDG() != 10022):
77  continue
78  p = mc_particle.getMomentum()
79  self.gamma_e[0] = mc_particle.getEnergy()
80  self.gamma_px[0] = p.X()
81  self.gamma_py[0] = p.Y()
82  self.gamma_pz[0] = p.Z()
83  jpsi = mc_particle.getDaughters()[0]
84  p = jpsi.getMomentum()
85  self.jpsi_e[0] = jpsi.getEnergy()
86  self.jpsi_px[0] = p.X()
87  self.jpsi_py[0] = p.Y()
88  self.jpsi_pz[0] = p.Z()
89  lepton = jpsi.getDaughters()[0]
90  p = lepton.getMomentum()
91  self.lepton_e[0] = lepton.getEnergy()
92  self.lepton_px[0] = p.X()
93  self.lepton_py[0] = p.Y()
94  self.lepton_pz[0] = p.Z()
95  self.tree.Fill()
96 
97  def terminate(self):
98  """ Termination function. """
99  self.output_file.cd()
100  self.tree.Write()
101  self.output_file.Close()
102 
103 
104 # Input.
105 root_input = register_module('RootInput')
106 root_input.param('inputFileName', 'PhokharaEvtgenData.root')
107 
108 # Analysis.
109 phokhara_evtgen = PhokharaEvtgenAnalysisModule()
110 
111 # Create main path.
112 main = create_path()
113 
114 # Add modules to main path
115 main.add_module(root_input)
116 main.add_module(phokhara_evtgen)
117 
118 # Run.
119 process(main)
120 
121 print(statistics)
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.terminate
def terminate(self)
Definition: PhokharaEvtgenAnalyze.py:97
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.__init__
def __init__(self)
Definition: PhokharaEvtgenAnalyze.py:22
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.gamma_py
gamma_py
Virtual photon momentum (y component).
Definition: PhokharaEvtgenAnalyze.py:36
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.output_file
output_file
Output file.
Definition: PhokharaEvtgenAnalyze.py:26
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.lepton_py
lepton_py
Lepton momentum (y component).
Definition: PhokharaEvtgenAnalyze.py:52
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule
Definition: PhokharaEvtgenAnalyze.py:19
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.ecms
ecms
Beam energy.
Definition: PhokharaEvtgenAnalyze.py:30
Belle2::PyStoreObj
a (simplified) python wrapper for StoreObjPtr.
Definition: PyStoreObj.h:69
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.lepton_pz
lepton_pz
Lepton momentum (z component).
Definition: PhokharaEvtgenAnalyze.py:54
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.gamma_e
gamma_e
Virtual photon energy.
Definition: PhokharaEvtgenAnalyze.py:32
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.jpsi_px
jpsi_px
J/psi momentum (x component).
Definition: PhokharaEvtgenAnalyze.py:42
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.jpsi_pz
jpsi_pz
J/psi momentum (z component).
Definition: PhokharaEvtgenAnalyze.py:46
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.lepton_e
lepton_e
Lepton energy.
Definition: PhokharaEvtgenAnalyze.py:48
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.tree
tree
Output tree.
Definition: PhokharaEvtgenAnalyze.py:28
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.jpsi_e
jpsi_e
J/psi energy.
Definition: PhokharaEvtgenAnalyze.py:40
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.gamma_pz
gamma_pz
Virtual photon momentum (z component).
Definition: PhokharaEvtgenAnalyze.py:38
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.lepton_px
lepton_px
Lepton momentum (x component).
Definition: PhokharaEvtgenAnalyze.py:50
Belle2::PyStoreArray
a (simplified) python wrapper for StoreArray.
Definition: PyStoreArray.h:58
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.jpsi_py
jpsi_py
J/psi momentum (y component).
Definition: PhokharaEvtgenAnalyze.py:44
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.gamma_px
gamma_px
Virtual photon momentum (x component).
Definition: PhokharaEvtgenAnalyze.py:34
PhokharaEvtgenAnalyze.PhokharaEvtgenAnalysisModule.event
def event(self)
Definition: PhokharaEvtgenAnalyze.py:69