Belle II Software  release-05-02-19
AAFHInterface.h
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010 - 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 #pragma once
12 #ifndef GENERATORS_AAFH_INTERFACE_H
13 #define GENERATORS_AAFH_INTERFACE_H
14 
15 #include <vector>
16 #include <mdst/dataobjects/MCParticleGraph.h>
17 
18 
19 namespace Belle2 {
31  class AAFHInterface {
32  public:
34  enum EMode {
45  };
46 
48  enum ERejection {
50  c_Once = 1,
53  c_Twice = 2,
54  };
55 
57  AAFHInterface() = default;
59  void setMaxTries(int tries);
61  void setGeneratorWeights(const std::vector<double>& weights);
65  void setMaxWeights(double subgeneratorWeight, double finalWeight);
67  void setSupressionLimits(std::vector<double> limits);
69  void setMinimumMass(double minMass) { m_minMass = minMass; }
71  void setParticle(const std::string& particle)
72  {
73  m_particle = particle;
74  }
75 
77  int getMaxTries() const;
79  std::vector<double> getGeneratorWeights() const;
81  double getMaxSubGeneratorWeight() const;
83  double getMaxFinalWeight() const;
85  std::vector<double> getSuppressionLimits() const;
87  double getMinimumMass() const { return m_minMass; }
89  std::string getParticle() const { return m_particle; }
91  double getTotalCrossSection() const;
93  double getTotalCrossSectionError() const;
94 
96  void initialize(double beamEnergy, EMode mode, ERejection rejection);
98  void generateEvent(MCParticleGraph& mpg);
100  void finish();
101  private:
108  void updateParticle(MCParticleGraph::GraphParticle& p, double* q, int pdg, bool isInitial = false);
110  double m_minMass {0};
112  std::string m_particle {"tau-"};
115  int m_particlePDG {0};
116  };
117 
119 } //Belle2 namespace
120 #endif // GENERATORS_AAFH_INTERFACE_H
Belle2::AAFHInterface::getTotalCrossSectionError
double getTotalCrossSectionError() const
Return error on the total cross section.
Definition: AAFHInterface.cc:222
Belle2::AAFHInterface::m_minMass
double m_minMass
minimum invariant mass
Definition: AAFHInterface.h:110
Belle2::AAFHInterface::c_ElectonMuon
@ c_ElectonMuon
E+E- -> E+E-Mu+Mu-.
Definition: AAFHInterface.h:40
Belle2::MCParticleGraph
Class to build, validate and sort a particle decay chain.
Definition: MCParticleGraph.h:48
Belle2::AAFHInterface::getParticle
std::string getParticle() const
Get the particle type for modes c_MuonParticle and c_ElectronParticle.
Definition: AAFHInterface.h:89
Belle2::AAFHInterface::getSuppressionLimits
std::vector< double > getSuppressionLimits() const
Get suppression limits.
Definition: AAFHInterface.cc:212
Belle2::AAFHInterface::updateParticle
void updateParticle(MCParticleGraph::GraphParticle &p, double *q, int pdg, bool isInitial=false)
update particle with generated values
Definition: AAFHInterface.cc:252
Belle2::AAFHInterface::setSupressionLimits
void setSupressionLimits(std::vector< double > limits)
Set the suppression limits when calculation the matrix elements.
Definition: AAFHInterface.cc:175
Belle2::AAFHInterface::m_particle
std::string m_particle
name of the particle for modes c_MuonParticle and c_ElectronParticle
Definition: AAFHInterface.h:112
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::AAFHInterface
Class to inferface AAFH/DIAG36 Generator written in Fortran.
Definition: AAFHInterface.h:31
Belle2::AAFHInterface::getTotalCrossSection
double getTotalCrossSection() const
Return total cross section.
Definition: AAFHInterface.cc:217
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::AAFHInterface::c_ElectonElectron
@ c_ElectonElectron
E+E- -> E+E-E+E-.
Definition: AAFHInterface.h:44
Belle2::AAFHInterface::getGeneratorWeights
std::vector< double > getGeneratorWeights() const
Get the relative weights for the sub generators.
Definition: AAFHInterface.cc:194
Belle2::AAFHInterface::c_MuonMuon
@ c_MuonMuon
E+E- -> Mu+Mu-Mu+Mu-.
Definition: AAFHInterface.h:38
Belle2::AAFHInterface::c_MuonParticle
@ c_MuonParticle
E+E- -> Mu+Mu-L+L- where L can be anything, defaults to tau.
Definition: AAFHInterface.h:36
Belle2::AAFHInterface::getMinimumMass
double getMinimumMass() const
Get the minimum invariant mass for the generated event.
Definition: AAFHInterface.h:87
Belle2::AAFHInterface::setGeneratorWeights
void setGeneratorWeights(const std::vector< double > &weights)
Set the relative weights for the sub generators.
Definition: AAFHInterface.cc:152
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::AAFHInterface::c_Once
@ c_Once
Apply rejection only once for the final event.
Definition: AAFHInterface.h:50
Belle2::AAFHInterface::initialize
void initialize(double beamEnergy, EMode mode, ERejection rejection)
initialize the generator
Definition: AAFHInterface.cc:227
Belle2::AAFHInterface::m_particlePDG
int m_particlePDG
pdg of the particle for modes c_MuonParticle and c_ElectronParticle, gets set in initialize()
Definition: AAFHInterface.h:115
Belle2::AAFHInterface::ERejection
ERejection
Rejection mode.
Definition: AAFHInterface.h:48
Belle2::AAFHInterface::EMode
EMode
Generation mode.
Definition: AAFHInterface.h:34
Belle2::AAFHInterface::getMaxFinalWeight
double getMaxFinalWeight() const
Get the maximum expected subgenerator weight for the rejection method.
Definition: AAFHInterface.cc:207
Belle2::AAFHInterface::getMaxTries
int getMaxTries() const
Get the maximum number of tries to generate an event.
Definition: AAFHInterface.cc:189
Belle2::AAFHInterface::c_ElectronParticle
@ c_ElectronParticle
E+E- -> E+E-L+L- where L can be anything, defaults to tau.
Definition: AAFHInterface.h:42
Belle2::AAFHInterface::AAFHInterface
AAFHInterface()=default
Default Constructor.
Belle2::AAFHInterface::generateEvent
void generateEvent(MCParticleGraph &mpg)
generate one event and add it to the graph in CMS
Definition: AAFHInterface.cc:269
Belle2::AAFHInterface::setMinimumMass
void setMinimumMass(double minMass)
Set the minimum invariant mass for the generated event.
Definition: AAFHInterface.h:69
Belle2::AAFHInterface::getMaxSubGeneratorWeight
double getMaxSubGeneratorWeight() const
Get the maximum expected final weight for the rejection method.
Definition: AAFHInterface.cc:202
Belle2::AAFHInterface::setMaxWeights
void setMaxWeights(double subgeneratorWeight, double finalWeight)
Set the maximum expected weights for the rejection method.
Definition: AAFHInterface.cc:169
Belle2::AAFHInterface::c_Twice
@ c_Twice
Apply rejection first for each sub generator and then for the final event.
Definition: AAFHInterface.h:53
Belle2::MCParticleGraph::GraphParticle
Class to represent Particle data in graph.
Definition: MCParticleGraph.h:86