Belle II Software development
AafhInputModule.cc
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * See git log for contributors and copyright holders. *
6 * This file is licensed under LGPL-3.0, see LICENSE.md. *
7 **************************************************************************/
8
9#include <generators/modules/aafhinput/AafhInputModule.h>
10#include <generators/utilities/scaleParticleEnergies.h>
11#include <framework/logging/Logger.h>
12
13#include <Math/Vector3D.h>
14
15using namespace Belle2;
16
17//-----------------------------------------------------------------
18// Register the Module
19//-----------------------------------------------------------------
20REG_MODULE(AafhInput);
21
22//-----------------------------------------------------------------
23// Implementation
24//-----------------------------------------------------------------
25
27{
28 // Set module properties
29 setDescription("AAFH Generator to generate non-radiative two-photon events like e+e- -> e+e-e+e-");
30
31 // Parameter definitions
32 addParam("mode", m_mode,
33 "decay mode to generate."
34 " 1: e+e- -> mu+mu-L+L-,"
35 " 2: e+e- -> mu+mu-mu+mu-,"
36 " 3: e+e- -> e+e-mu+mu-,"
37 " 4: e+e- -> e+e-L+L-,"
38 " 5: e+e- -> e+e-e+e-. "
39 "L is a user defined particle (default: tau)",
40 5);
41 addParam("rejection", m_rejection,
42 "Rejection method."
43 " 1: apply once for the final event weight,"
44 " 2: apply for the subgenerator weight and for the final weight",
45 2);
46 addParam("maxTries", m_maxTries,
47 "Maximum tries for event generation (ITOT)",
48 m_generator.getMaxTries());
49 addParam("maxFinalWeight", m_maxFinalWeight,
50 "maximum expected final weight for rejection scheme (ESFT)",
51 m_generator.getMaxFinalWeight());
52 addParam("maxSubgeneratorWeight", m_maxSubgeneratorWeight,
53 "maximum expected subgenerator weight for rejection scheme (ESWE)",
54 m_generator.getMaxSubGeneratorWeight());
55 addParam("subgeneratorWeights", m_subgeneratorWeights,
56 "relative weights of the subgenerators: this must be a list of four "
57 "or eight values (first four are interpreted as WAP, rest as WBP) "
58 "which specify the relative weights for each of the "
59 "four sub generators. The original code states that it the program "
60 "run most efficient when the maximum weight is equal in all sub "
61 "generators and that if one wants to be sure that all peaks in the "
62 "differential cross section are accounted the chance to enter each "
63 "sub generator should be equal. Values which try to fulfill both "
64 "conditions are printed at after generation when the output level "
65 "is set to INFO",
66 m_generator.getGeneratorWeights());
67 addParam("suppressionLimits", m_suppressionLimits,
68 "suppression limits for the matrix element calculations. This must "
69 "be a list of four values, [FACE, FACM, FACL, PROC]. For FACE, FACM "
70 "and FACL these specify the size of the propagator for which the "
71 "corresponding spin configurations will be omitted. A value of 1e3 "
72 "will calculate the dominat part correctly and a very large value "
73 "compared to the mass (i.e. 1e50) will calculate it exactly. PROC "
74 "specifies that feynman diagrams which contribute less than 1/PROC "
75 "of the biggest contribution are omitted. For exact calculation it "
76 "should be very big but 1e9 should be considerably faster without "
77 "affecting the result very much",
78 m_generator.getSuppressionLimits());
79 addParam("minMass", m_minMass,
80 "minimum invariant mass of generated fermion pair. Will "
81 "automatically be set to at least 2 times the generated particle "
82 "mass",
83 0.);
84 addParam("particle", m_particle,
85 "name of the generated particle for mode 1 and 4",
86 m_generator.getParticle());
87}
88
90{
91 m_mcparticles.registerInDataStore();
92
93 //Beam Parameters, initial particle - AAFH cannot handle beam energy spread
94 m_initial.initialize();
95}
96
98{
99
100 // Check if the BeamParameters have changed (if they do, abort the job! otherwise cross section calculation will be a nightmare.)
101 if (m_beamParams.hasChanged()) {
102 if (!m_initialized) {
104 } else {
105 B2FATAL("AafhInputModule::event(): BeamParameters have changed within a job, this is not supported for AAFH!");
106 }
107 }
108
109 // Initial particle from beam parameters (for random vertex)
110 const MCInitialParticles& initial = m_initial.generate();
111
112 // get Lorentz transformation from CMS to LAB
113 ROOT::Math::LorentzRotation boost = initial.getCMSToLab();
114
115 // vertex.
116 ROOT::Math::XYZVector vertex = initial.getVertex();
117
118 MCParticleGraph mpg;
119
120 // generate event
121 m_generator.generateEvent(mpg);
122
123 // scale CMS energy of generated particles to initial.getMass()
124 scaleParticleEnergies(mpg, initial.getMass());
125
126 // transform to lab and set vertex.
127 for (size_t i = 0; i < mpg.size(); ++i) {
128 mpg[i].set4Vector(boost * mpg[i].get4Vector());
129
130 ROOT::Math::XYZVector v3 = mpg[i].getProductionVertex();
131 v3 = v3 + vertex;
132 mpg[i].setProductionVertex(v3);
133 mpg[i].setValidVertex(true);
134 }
135
136 //Fill MCParticle List
138}
139
141{
142 m_generator.finish();
143}
144
146{
147 const BeamParameters& nominal = m_initial.getBeamParameters();
148 const double beamEnergy = nominal.getMass() / 2.;
149
150 //Set Generator options
152 B2ERROR("AafhInputModule::initializeGenerator: 'mode' must be a value between 1 and 5");
153 }
155 B2ERROR("AafhInputModule::initializeGenerator: 'rejection' must be a value between 1 and 2");
156 }
157 m_generator.setMaxTries(m_maxTries);
158 m_generator.setParticle(m_particle);
159 m_generator.setMinimumMass(m_minMass);
160 m_generator.setGeneratorWeights(m_subgeneratorWeights);
161 m_generator.setSupressionLimits(m_suppressionLimits);
163 m_generator.initialize(beamEnergy, (AAFHInterface::EMode) m_mode,
165
166 m_initialized = true;
167
168}
169
170
ERejection
Rejection mode.
EMode
Generation mode.
double m_minMass
minimum invariant mass of generated fermion pair
DBObjPtr< BeamParameters > m_beamParams
BeamParameter.
bool m_initialized
True if generator has been initialized.
int m_rejection
Rejection mode.
std::vector< double > m_suppressionLimits
suppression limits for the matrix element calculations
StoreArray< MCParticle > m_mcparticles
MCParticle collection.
int m_mode
Generator mode.
std::string m_particle
name of the generated particle for modes c_MuonParticle and c_ElectronParticle
virtual void terminate() override
calculate cross section
int m_maxTries
maximum number of tries for generating an event
virtual void generatorInitialize() override
initialize generator
std::vector< double > m_subgeneratorWeights
relative subgenerator weights
AafhInputModule()
Constructor: Sets the description, the properties and the parameters of the module.
double m_maxSubgeneratorWeight
maximum expected subgenerator weight for rejection scheme
InitialParticleGeneration m_initial
initial particle used by BeamParameter class
double m_maxFinalWeight
maximum expected final weight for rejection scheme
AAFHInterface m_generator
interface to the generator
void initializeGenerator()
Method is called to initialize the generator.
virtual void generatorEvent() override
generate event
This class contains the nominal beam parameters and the parameters used for smearing of the primary v...
This class contains the initial state for the given event.
const ROOT::Math::LorentzRotation & getCMSToLab() const
Return the LorentzRotation to convert from CMS to lab frame.
const ROOT::Math::XYZVector & getVertex() const
Get the position of the collision.
double getMass() const
Get the invariant mass of the collision (= energy in CMS)
Class to build, validate and sort a particle decay chain.
@ c_setDecayInfo
Set decay time and vertex.
size_t size() const
Return the number of particles in the graph.
void generateList(const std::string &name="", int options=c_setNothing)
Generates the MCParticle list and stores it in the StoreArray with the given name.
void setDescription(const std::string &description)
Sets the description of the module.
Definition Module.cc:214
void addParam(const std::string &name, T &paramVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
Definition Module.h:559
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition Module.h:649
void scaleParticleEnergies(MCParticleGraph &mpg, double EcmsTarget)
Scale momenta of the particles by a constant factor such that total energy is changed to EcmsTarget.
Abstract base class for different kinds of events.