herwig is hosted by Hepforge, IPPP Durham
Herwig  7.2.1
GenericWidthGenerator.h
1 // -*- C++ -*-
2 //
3 // GenericWidthGenerator.h is a part of Herwig - A multi-purpose Monte Carlo event generator
4 // Copyright (C) 2002-2019 The Herwig Collaboration
5 //
6 // Herwig is licenced under version 3 of the GPL, see COPYING for details.
7 // Please respect the MCnet academic guidelines, see GUIDELINES for details.
8 //
9 #ifndef HERWIG_GenericWidthGenerator_H
10 #define HERWIG_GenericWidthGenerator_H
11 //
12 // This is the declaration of the GenericWidthGenerator class.
13 //
14 #include "ThePEG/PDT/WidthGenerator.h"
15 #include "ThePEG/PDT/ParticleData.h"
16 #include "ThePEG/PDT/DecayMode.h"
17 #include "GenericWidthGenerator.fh"
18 #include "Herwig/Decay/DecayIntegrator.h"
19 #include "Herwig/Decay/DecayPhaseSpaceMode.h"
20 #include "Herwig/Utilities/Interpolator.h"
21 #include "GenericMassGenerator.h"
22 #include <iostream>
23 
24 namespace Herwig {
25 using namespace ThePEG;
26 
31 
35  typedef Selector<tDMPtr> DecayMap;
36 
37 
52 
56  friend class ModelGenerator;
57 
58 public:
59 
63  friend class TwoBodyAllOnCalculator;
64 
65 public:
66 
71  : mass_(), prefactor_(1.), initialize_(false),output_(false),
72  BRnorm_(true),npoints_(50),
73  BRminimum_(0.01), intOrder_(1), twoBodyOnly_(false)
74  {}
75 
82  void persistentOutput(PersistentOStream & os) const;
83 
89  void persistentInput(PersistentIStream & is, int version);
91 
95  static void Init();
96 
97 public:
98 
104  virtual bool accept(const ParticleData & part) const {
105  if(!particle_) return false;
106  return part.id() == particle_->id() ||
107  ( part.CC() && part.CC()->id() == particle_->id() );
108  }
109 
118  virtual Energy width(const ParticleData & part, Energy m) const;
119 
125  virtual DecayMap rate(const ParticleData & part) const {
126  return part.decaySelector();
127  }
128 
135  virtual DecayMap rate (const Particle & part);
136 
143  virtual Energy partialWidth(int iloc,Energy m) const;
144 
149  virtual pair<Energy,Energy> width(Energy, const ParticleData &) const;
151 
157  virtual void dataBaseOutput(ofstream & output, bool header=true);
158 
163  virtual Length lifeTime(const ParticleData &, Energy m, Energy w) const;
164 
165 protected:
166 
173  Energy partial2BodyWidth(int iloc,Energy q) const {
174  return partial2BodyWidth(iloc,q,MEmass1_[iloc],MEmass2_[iloc]);
175  }
176 
185  virtual Energy partial2BodyWidth(int iloc,Energy m0,Energy m1,Energy m2) const;
186 
193  virtual void setupMode(tcDMPtr mode, tDecayIntegratorPtr decayer, unsigned int imode);
194 
198  tPDPtr particle() const {return particle_;}
199 
203  void particle(tPDPtr in) {particle_ = in;}
204 
205 protected:
206 
213  virtual IBPtr clone() const {return new_ptr(*this);}
214 
219  virtual IBPtr fullclone() const {return new_ptr(*this);}
221 
222 protected:
223 
231  virtual void doinit();
232 
237  virtual void dofinish();
238 
248  virtual void rebind(const TranslationMap & trans)
249  ;
250 
256  virtual IVector getReferences();
258 
263  int MEcode(int imode) const {return MEcode_[imode];}
264 
269  double MEcoupling(int imode) const {return MEcoupling_[imode];}
270 
271 
275  Energy mass() const {return mass_;}
276 
284  bool initialize() const {return initialize_;}
285 
289  bool output() const {return output_;}
290 
294  vector<tDMPtr> decayModes() const {return decayModes_;}
295 
296 private:
297 
301  void setParticle(string);
302 
306  string getParticle() const;
307 
308 private:
309 
313  GenericWidthGenerator & operator=(const GenericWidthGenerator &) = delete;
314 
315 private:
316 
321 
325  vector<tDMPtr> decayModes_;
326 
330  vector<string> decayTags_;
331 
335  vector<Energy> minMass_;
336 
340  Energy mass_;
341 
345  double prefactor_;
346 
350  vector<int> MEtype_;
351 
355  vector<int> MEcode_;
356 
360  vector<Energy> MEmass1_;
364  vector<Energy> MEmass2_;
365 
369  vector<double> MEcoupling_;
370 
374  vector<bool> modeOn_;
375 
379  vector<Energy> interMasses_;
380 
384  vector<Energy> interWidths_;
385 
389  vector<int> noOfEntries_;
390 
395 
399  bool output_;
400 
404  bool BRnorm_;
405 
409  int npoints_;
410 
414  vector<Interpolator<Energy,Energy>::Ptr> interpolators_;
415 
419  double BRminimum_;
420 
424  unsigned int intOrder_;
425 
431 };
432 
433 }
434 
435 #endif /* HERWIG_GenericWidthGenerator_H */
bool initialize() const
Access to the decay modes.
vector< Interpolator< Energy, Energy >::Ptr > interpolators_
intepolators for the running width
virtual DecayMap rate(const ParticleData &part) const
Initialize the given decay map for the given particle type.
ThePEG::Ptr< DecayMode >::transient_const_pointer tcDMPtr
Energy mass() const
The on-shell mass of the particle.
ThePEG::Ptr< ParticleData >::transient_pointer tPDPtr
ThePEG::Ptr< InterfacedBase >::pointer IBPtr
bool initialize_
initialize the generator using the particle data object
int npoints_
number of points to use for interpolation tables
vector< tDMPtr > decayModes() const
Access to the DecayModes.
tPDPtr particle_
The pointer to the ParticleData object for the particle for this width generator. ...
vector< Energy > interWidths_
storage of the widths to set up the interpolation tables
vector< int > noOfEntries_
the number of entries in the decay table for a particular mode
vector< Energy > MEmass2_
Mass of the second outgoing particle for the simple ME&#39;s.
vector< int > MEtype_
The type of ME, whether it is fixed, calculated by this class or interpolation.
bool output_
Output the parameters.
tPDPtr CC() const
vector< double > MEcoupling_
the coupling for a given me
The GenericWidthGenerator class is designed to automatically calculate the running width for a given ...
double prefactor_
Prefactor to get the on-shell width.
vector< tDMPtr > decayModes_
The decaymodes.
const DecaySelector & decaySelector() const
int MEcode(int imode) const
Matrix element code for a given mode.
double MEcoupling(int imode) const
Coupling for a given mode.
double BRminimum_
minimum branching ratio for the inclusion in the total running width
virtual bool accept(const ParticleData &part) const
Return true if this width generator can handle the given particle type.
bool twoBodyOnly_
Whether or not to only include 2 body modes in the running width calculation, higher modes flat...
Selector< tDMPtr > DecayMap
Typedef to define a DecayMoap.
vector< Energy > MEmass1_
Mass of the first outgoing particle for the simple ME&#39;s.
virtual IBPtr fullclone() const
Make a clone of this object, possibly modifying the cloned object to make it sane.
vector< Energy > minMass_
The minimum mass of the decaying particle for which this decay mode is possible.
GenericWidthGenerator()
Default constructor.
unsigned int intOrder_
Order of the interpolation for the tables.
vector< int > MEcode_
The code for the matrix element.
tPDPtr particle() const
Access to the particle dat for inheriting classes.
vector< Energy > interMasses_
storage of the massesto set up the interpolation tables
-*- C++ -*-
Energy partial2BodyWidth(int iloc, Energy q) const
The width for on-shell particles.
This class is designed to store the particles in some model and then call the appropriate function to...
The TwoBodyAllOnCalculator class is a wrapped around the the simple two body decay matrix elements in...
Energy mass_
The on-shell mass of the particle.
void particle(tPDPtr in)
Set the particle.
vector< IBPtr > IVector
bool BRnorm_
normalise the terms so that the partial widths for an on-shell particle are correct ...
vector< string > decayTags_
The tags for the DecayMode s.
bool output() const
Output option for use by the inheriting classes.
virtual IBPtr clone() const
Make a simple clone of this object.
vector< bool > modeOn_
is this mode used for the running width