herwig is hosted by Hepforge, IPPP Durham
Herwig 7.3.0
DecayIntegrator.h
1// -*- C++ -*-
2//
3// DecayIntegrator.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_DecayIntegrator_H
10#define Herwig_DecayIntegrator_H
11//
12// This is the declaration of the DecayIntegrator class.
13//
14
15#include "DecayIntegrator.fh"
16#include "HwDecayerBase.h"
17#include "PhaseSpaceMode.fh"
18#include "Herwig/PDT/WidthCalculatorBase.fh"
19#include "Radiation/DecayRadiationGenerator.h"
20#include <Herwig/Decay/DecayVertex.h>
21
22namespace Herwig {
23
24using namespace ThePEG;
25
61
62public:
63
67 friend class PhaseSpaceMode;
68
72 enum MEOption {Initialize,Calculate,Terminate};
73
74public:
75
79 DecayIntegrator() : nIter_(10), nPoint_(10000), nTry_(500),
80 generateInter_(false), iMode_(-1),
81 realME_(false), virtualME_(false), eps_(ZERO), warnings_(false)
82 {}
83
84public:
85
92 virtual bool accept(tcPDPtr parent, const tPDVector & children) const {
93 bool cc;
94 return modeNumber(cc,parent,children)>=0;
95 }
96
103 virtual ParticleVector decay(const Particle & parent,
104 const tPDVector & children) const;
105
112 virtual int modeNumber(bool & cc, tcPDPtr parent,
113 const tPDVector & children) const = 0;
114
118 int imode() const {return iMode_;}
119
124 void addMode(PhaseSpaceModePtr mode) const;
125
135 virtual double me2(const int ichan, const Particle & part,
136 const tPDVector & outgoing,
137 const vector<Lorentz5Momentum> & momenta,
138 MEOption meopt) const = 0;
139
143 virtual void constructSpinInfo(const Particle & part,
144 ParticleVector outgoing) const = 0;
145
151 virtual void dataBaseOutput(ofstream & os,bool header) const;
152
159 void setPartialWidth(const DecayMode & dm, int imode);
167 virtual bool twoBodyMEcode(const DecayMode &, int & mecode,
168 double & coupling) const {
169 coupling = 1.;
170 mecode = -1;
171 return false;
172 }
173
179 virtual WidthCalculatorBasePtr threeBodyMEIntegrator(const DecayMode & dm) const;
180
194 virtual double threeBodyMatrixElement(const int imode, const Energy2 q2,
195 const Energy2 s3, const Energy2 s2,
196 const Energy2 s1, const Energy m1,
197 const Energy m2, const Energy m3) const;
198
209 virtual InvEnergy threeBodydGammads(const int imode, const Energy2 q2,
210 const Energy2 s,
211 const Energy m1, const Energy m2,
212 const Energy m3) const;
213
218 int findMode(const DecayMode & dm);
219
220public:
221
234 return photonGen_->generatePhotons(p,children,this);
235 }
236
241
249 virtual double oneLoopVirtualME(unsigned int imode,
250 const Particle & part,
251 const ParticleVector & products);
252
256 bool hasOneLoopME() {return virtualME_;}
257
271 virtual InvEnergy2 realEmissionME(unsigned int imode,
272 const Particle & part,
273 ParticleVector & products,
274 unsigned int iemitter,
275 double ctheta, double stheta,
276 const LorentzRotation & rot1,
277 const LorentzRotation & rot2);
278
282 bool hasRealEmissionME() {return realME_;}
284
285public:
286
290 friend ostream & operator<<(ostream & os, const DecayIntegrator & decay);
291
292
293public:
294
302
308 void persistentInput(PersistentIStream & is, int version);
310
317 static void Init();
318
319
320protected:
321
328 virtual void doinitrun();
330
331protected:
332
342 ParticleVector generate(bool inter,bool cc, const unsigned int & imode,
343 const Particle & inpart) const;
344
348 void imode(int in) { iMode_ = in;}
349
353 void ME(DecayMEPtr in) const { matrixElement_ = in;}
354
358 DecayMEPtr ME() const {return matrixElement_;}
359
366 void resetIntermediate(tcPDPtr part, Energy mass, Energy width);
367
373 Energy initializePhaseSpaceMode(unsigned int imode,bool init, bool onShell=false) const;
374
375protected:
376
385
390
394 void hasOneLoopME(bool in) {virtualME_=in;}
395
399 void hasRealEmissionME(bool in) {realME_=in;}
400
404 void epsilonPS(Energy in) {eps_=in;}
405
409 void clearModes() {modes_.clear();}
410
411protected:
412
416 unsigned int numberModes() const {return modes_.size();}
417
421 tPhaseSpaceModePtr mode(unsigned int ix) {
422 return modes_[ix];
423 }
427 tcPhaseSpaceModePtr mode(unsigned int ix) const {
428 return modes_[ix];
429 }
430
431public:
432
433bool warnings() const {
434 return warnings_;
435}
436
437private:
438
443
451 unsigned int nIter_;
452
456 unsigned int nPoint_;
457
461 unsigned int nTry_;
462
466 mutable vector<PhaseSpaceModePtr> modes_;
467
469
474
478 DecayRadiationGeneratorPtr photonGen_;
479
483 mutable int iMode_;
484
488 mutable DecayMEPtr matrixElement_;
489
494
499
503 Energy eps_;
504
509
510protected:
511
516
517};
518
522ostream & operator<<(ostream &, const DecayIntegrator &);
523
524}
525
526#endif /* Herwig_DecayIntegrator_H */
Exception for this class and those inheriting from it.
Main class for Decayers implementing multi-channel phase space integration.
bool warnings_
option for turinh on/off log warnings in Phase class
int imode() const
The mode being used for this decay.
unsigned int nTry_
number of attempts to generate the decay
void clearModes()
Clear the models.
void epsilonPS(Energy in)
Set the epsilon parameter.
tcPhaseSpaceModePtr mode(unsigned int ix) const
Pointer to a mode.
virtual double me2(const int ichan, const Particle &part, const tPDVector &outgoing, const vector< Lorentz5Momentum > &momenta, MEOption meopt) const =0
Return the matrix element squared for a given mode and phase-space channel.
friend ostream & operator<<(ostream &os, const DecayIntegrator &decay)
The output operator is a friend, this is mainly for debugging.
DecayIntegrator()
The default constructor.
void resetIntermediate(tcPDPtr part, Energy mass, Energy width)
Reset the properities of all intermediates.
void hasRealEmissionME(bool in)
Whether or not the real emission matrix element is implemented.
virtual InvEnergy2 realEmissionME(unsigned int imode, const Particle &part, ParticleVector &products, unsigned int iemitter, double ctheta, double stheta, const LorentzRotation &rot1, const LorentzRotation &rot2)
The real emission matrix element.
int iMode_
mode currently being generated
virtual void doinitrun()
Initialize this object.
void addMode(PhaseSpaceModePtr mode) const
Add a phase-space mode to the list.
virtual InvEnergy threeBodydGammads(const int imode, const Energy2 q2, const Energy2 s, const Energy m1, const Energy m2, const Energy m3) const
The differential three body decay rate with one integral performed.
virtual double oneLoopVirtualME(unsigned int imode, const Particle &part, const ParticleVector &products)
The one-loop virtual correction.
vector< PhaseSpaceModePtr > modes_
List of the decay modes.
void persistentOutput(PersistentOStream &os) const
Function used to write out object persistently.
void setPartialWidth(const DecayMode &dm, int imode)
Set the code for the partial width.
bool hasRealEmissionME()
Whether or not the real emission matrix element is implemented.
Energy initializePhaseSpaceMode(unsigned int imode, bool init, bool onShell=false) const
Initialize the phase-space mode.
DecayRadiationGeneratorPtr photonGen_
Pointer to the object generating the QED radiation in the decay.
ParticleVector generate(bool inter, bool cc, const unsigned int &imode, const Particle &inpart) const
Generate the momenta for the decay.
DecayMEPtr matrixElement_
The helicity matrix element for the current decay.
bool canGeneratePhotons()
check if photons can be generated in the decay
bool hasOneLoopME()
Whether or not the one loop matrix element is implemented.
void imode(int in)
Set the mode being use for this decay.
virtual double threeBodyMatrixElement(const int imode, const Energy2 q2, const Energy2 s3, const Energy2 s2, const Energy2 s1, const Energy m1, const Energy m2, const Energy m3) const
The matrix element to be integrated for the three-body decays as a function of the invariant masses o...
virtual int modeNumber(bool &cc, tcPDPtr parent, const tPDVector &children) const =0
Which of the possible decays is required.
int findMode(const DecayMode &dm)
Finds the phase-space mode corresponding to a given decay mode.
virtual WidthCalculatorBasePtr threeBodyMEIntegrator(const DecayMode &dm) const
Method to return an object to calculate the 3 (or higher body) partial width.
void hasOneLoopME(bool in)
Whether or not the one loop matrix element is implemented.
unsigned 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.
DecayIntegrator & operator=(const DecayIntegrator &)=delete
Private and non-existent assignment operator.
bool virtualME_
Whether or not the one-loop matrix element exists.
bool generateInter_
Whether to include the intermediates whne outputing the results.
void ME(DecayMEPtr in) const
Set the helicity matrix element for the decay.
Energy eps_
Epsilon parameter for phase-space integration.
bool generateIntermediates() const
Set whether or not the intermediates are included.
unsigned int nIter_
Parameters for the integration.
static void Init()
The standard Init function used to initialize the interfaces.
virtual void dataBaseOutput(ofstream &os, bool header) const
Output the setup information for the particle database.
unsigned int numberModes() const
Number of decay modes.
DecayMEPtr ME() const
The helicity amplitude matrix element for spin correlations.
void persistentInput(PersistentIStream &is, int version)
Function used to read in object persistently.
virtual bool accept(tcPDPtr parent, const tPDVector &children) const
Check if this decayer can perfom the decay for a particular mode.
virtual ParticleVector decay(const Particle &parent, const tPDVector &children) const
For a given decay mode and a given particle instance, perform the decay and return the decay products...
virtual void constructSpinInfo(const Particle &part, ParticleVector outgoing) const =0
Construct the SpinInfos for the particles produced in the decay.
MEOption
Enum for the matrix element option.
void generateIntermediates(bool in)
Methods to set variables in inheriting classes.
tPhaseSpaceModePtr mode(unsigned int ix)
Pointer to a mode.
bool realME_
Whether or not the real photon emission matrix element exists.
ParticleVector generatePhotons(const Particle &p, ParticleVector children)
Members for the generation of QED radiation in the decays.
The HwDecayerBase class is the base class for Decayers in Herwig.
Definition: HwDecayerBase.h:37
The PhaseSpaceMode class is designed to store a group of phase-space channels for use by the DecayInt...
-*- C++ -*-
vector< tPDPtr > tPDVector
vector< PPtr > ParticleVector
constexpr ZeroUnit ZERO
ThePEG::Ptr< ParticleData >::transient_const_pointer tcPDPtr
ostream & operator<<(ostream &, const Collision &)