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