herwig is hosted by Hepforge, IPPP Durham
Herwig  7.1.5
DecayIntegrator.h
1 // -*- C++ -*-
2 //
3 // DecayIntegrator.h is a part of Herwig - A multi-purpose Monte Carlo event generator
4 // Copyright (C) 2002-2017 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_DecayIntegrator_H
10 #define HERWIG_DecayIntegrator_H
11 //
12 // This is the declaration of the DecayIntegrator class.
13 //
14 #include <ThePEG/PDT/Decayer.h>
15 #include "DecayPhaseSpaceChannel.h"
16 #include <ThePEG/PDT/EnumParticles.h>
17 #include <Herwig/Decay/DecayVertex.h>
18 #include "DecayPhaseSpaceMode.fh"
19 #include "Herwig/PDT/WidthCalculatorBase.fh"
20 #include "Radiation/DecayRadiationGenerator.h"
21 #include "HwDecayerBase.h"
22 #include "DecayIntegrator.fh"
23 
24 namespace Herwig {
25 using namespace ThePEG;
26 
63 
64 public:
65 
69  friend ostream & operator<<(ostream &, const DecayIntegrator &);
70 
74  friend class DecayPhaseSpaceMode;
75 
79  enum MEOption {Initialize,Calculate,Terminate};
80 
81 public:
82 
86  DecayIntegrator() : _niter(10), _npoint(10000), _ntry(500),
87  _generateinter(false), _imode(-1),
88  _realME(false), _virtualME(false) {}
89 
96  virtual bool accept(tcPDPtr parent, const tPDVector & children) const {
97  bool cc;
98  return modeNumber(cc,parent,children)>=0;
99  }
100 
107  virtual ParticleVector decay(const Particle & parent,
108  const tPDVector & children) const;
109 
116  virtual int modeNumber(bool & cc, tcPDPtr parent,
117  const tPDVector & children) const = 0;
118 
122  int imode() const {return _imode;}
123 
130  void addMode(DecayPhaseSpaceModePtr mode,double maxwgt,
131  const vector<double> wgts) const;
132 
143  virtual double me2(const int ichan, const Particle & part,
144  const ParticleVector & decay,MEOption opt) const=0;
145 
149  DecayMEPtr ME() const {return _matrixelement;}
150 
158  virtual bool twoBodyMEcode(const DecayMode &, int & mecode,
159  double & coupling) const {
160  coupling = 1.;
161  mecode = -1;
162  return false;
163  }
164 
170  virtual WidthCalculatorBasePtr threeBodyMEIntegrator(const DecayMode & dm) const;
171 
185  virtual double threeBodyMatrixElement(const int imode, const Energy2 q2,
186  const Energy2 s3, const Energy2 s2,
187  const Energy2 s1, const Energy m1,
188  const Energy m2, const Energy m3) const;
189 
200  virtual InvEnergy threeBodydGammads(const int imode, const Energy2 q2,
201  const Energy2 s,
202  const Energy m1, const Energy m2,
203  const Energy m3) const;
204 
211  void setPartialWidth(const DecayMode & dm, int imode);
212 
217  int findMode(const DecayMode & dm);
218 
224  virtual void dataBaseOutput(ofstream & os,bool header) const;
225 
226 public:
227 
240  return _photongen->generatePhotons(p,children,this);
241  }
242 
246  bool canGeneratePhotons() {return _photongen;}
247 
255  virtual double oneLoopVirtualME(unsigned int imode,
256  const Particle & part,
257  const ParticleVector & products);
258 
262  bool hasOneLoopME() {return _virtualME;}
263 
277  virtual InvEnergy2 realEmissionME(unsigned int imode,
278  const Particle & part,
279  ParticleVector & products,
280  unsigned int iemitter,
281  double ctheta, double stheta,
282  const LorentzRotation & rot1,
283  const LorentzRotation & rot2);
284 
288  bool hasRealEmissionME() {return _realME;}
290 public:
291 
298  void persistentOutput(PersistentOStream & os) const;
299 
305  void persistentInput(PersistentIStream & is, int version);
307 
311  static void Init();
312 
313 protected:
314 
324  ParticleVector generate(bool inter,bool cc, const unsigned int & imode,
325  const Particle & inpart) const;
326 
330  void imode(int in) {_imode=in;}
331 
335  void ME(DecayMEPtr in) const { _matrixelement = in;}
336 
343  void resetIntermediate(tcPDPtr part, Energy mass, Energy width);
344 
348  unsigned int numberModes() const {return _modes.size();}
349 
353  tDecayPhaseSpaceModePtr mode(unsigned int);
354 
358  tcDecayPhaseSpaceModePtr mode(unsigned int) const;
359 
363  bool generateIntermediates() const {return _generateinter;}
364 
368  void generateIntermediates(bool in) {_generateinter=in;}
369 
375  Energy initializePhaseSpaceMode(unsigned int imode,bool init, bool onShell=false) const;
376 
380  void hasOneLoopME(bool in) {_virtualME=in;}
381 
385  void hasRealEmissionME(bool in) {_realME=in;}
386 
387 protected:
388 
394  virtual void doinitrun();
396 
397 private:
398 
403 
407  DecayIntegrator & operator=(const DecayIntegrator &) = delete;
408 
409 private:
410 
414  int _niter;
415 
419  int _npoint;
420 
424  int _ntry;
425 
429  mutable vector<DecayPhaseSpaceModePtr> _modes;
430 
435 
439  DecayRadiationGeneratorPtr _photongen;
440 
444  mutable int _imode;
445 
449  mutable DecayMEPtr _matrixelement;
450 
454  bool _realME;
455 
460 
461 };
465  ostream & operator<<(ostream &, const DecayIntegrator &);
466 
470  class DecayIntegratorError: public Exception {};
471 
472 }
473 
474 
475 namespace ThePEG {
476 
483 template <>
484 struct BaseClassTrait<Herwig::DecayIntegrator,1> {
486  typedef Herwig::HwDecayerBase NthBase;
487 };
488 
493 template <>
494 struct ClassTraits<Herwig::DecayIntegrator>
495  : public ClassTraitsBase<Herwig::DecayIntegrator> {
497  static string className() { return "Herwig::DecayIntegrator"; }
498 };
499 
502 }
503 
504 #endif /* HERWIG_DecayIntegrator_H */
bool _generateinter
Whether to include the intermediates whne outputing the results.
vector< DecayPhaseSpaceModePtr > _modes
List of the decay modes.
bool hasOneLoopME()
Whether or not the one loop matrix element is implemented.
Exception for this class and those inheriting from it.
void generateIntermediates(bool in)
Set whether or not the intermediates are included.
bool generateIntermediates() const
Get whether or not the intermediates are included.
The DecayPhaseSpaceMode class is designed to store a group of phase-space channels for use by the Dec...
DecayMEPtr _matrixelement
The helicity matrix element for the current decay.
void hasOneLoopME(bool in)
Whether or not the one loop matrix element is implemented.
unsigned int numberModes() const
Number of decay modes.
void imode(int in)
Set the mode being use for this decay.
ParticleVector generatePhotons(const Particle &p, ParticleVector children)
Members for the generation of QED radiation in the decays.
MEOption
Enum for the matrix element option.
ostream & operator<<(ostream &, const DecayIntegrator &)
Output information on the DecayIntegrator for debugging purposes.
bool _virtualME
Whether or not the one-loop matrix element exists.
virtual bool accept(tcPDPtr parent, const tPDVector &children) const
Check if this decayer can perfom the decay for a particular mode.
bool canGeneratePhotons()
check if photons can be generated in the decay
vector< tPDPtr > tPDVector
DecayRadiationGeneratorPtr _photongen
Pointer to the object generating the QED radiation in the decay.
void ME(DecayMEPtr in) const
Set the helicity matrix element for the decay.
static AbstractClassDescription< DecayIntegrator > initDecayIntegrator
Describe an abstract base class with persistent data.
Main class for Decayers implementing multi-channel phase space integration.
ThePEG::Ptr< ParticleData >::transient_const_pointer tcPDPtr
bool hasRealEmissionME()
Whether or not the real emission matrix element is implemented.
int _imode
mode currently being generated
void hasRealEmissionME(bool in)
Whether or not the real emission matrix element is implemented.
DecayIntegrator()
Default constructor.
pair< double, double > generate(const Generator< Density > &gen, double r)
Generate a random variable and return its weight.
The HwDecayerBase class is the base class for Decayers in Herwig.
Definition: HwDecayerBase.h:39
-*- C++ -*-
int _npoint
Number of points for initialisation.
virtual bool twoBodyMEcode(const DecayMode &, int &mecode, double &coupling) const
Specify the matrix element to be used in the running width calculation.
vector< PPtr > ParticleVector
bool _realME
Whether or not the real photon emission matrix element exists.
DecayMEPtr ME() const
The helicity amplitude matrix element for spin correlations.
int imode() const
The mode being used for this decay.
int _ntry
number of attempts to generate the decay
int _niter
Number of iterations for th initialization.