herwig
is hosted by
Hepforge
,
IPPP Durham
Herwig
7.2.1
Decay
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
#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
62
class
DecayIntegrator
:
public
HwDecayerBase
{
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), _eps(
ZERO
) {}
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
229
Energy
epsilonPS
()
const
{
return
_eps;}
230
231
public
:
232
244
ParticleVector
generatePhotons
(
const
Particle
& p,
ParticleVector
children) {
245
return
_photongen->generatePhotons(p,children,
this
);
246
}
247
251
bool
canGeneratePhotons
() {
return
_photongen;}
252
260
virtual
double
oneLoopVirtualME(
unsigned
int
imode,
261
const
Particle
& part,
262
const
ParticleVector
& products);
263
267
bool
hasOneLoopME
() {
return
_virtualME;}
268
282
virtual
InvEnergy2 realEmissionME(
unsigned
int
imode,
283
const
Particle
& part,
284
ParticleVector
& products,
285
unsigned
int
iemitter,
286
double
ctheta,
double
stheta,
287
const
LorentzRotation
& rot1,
288
const
LorentzRotation
& rot2);
289
293
bool
hasRealEmissionME
() {
return
_realME;}
295
public
:
296
303
void
persistentOutput(
PersistentOStream
& os)
const
;
304
310
void
persistentInput(
PersistentIStream
& is,
int
version);
312
316
static
void
Init();
317
318
protected
:
319
329
ParticleVector
generate
(
bool
inter,
bool
cc,
const
unsigned
int
& imode,
330
const
Particle
& inpart)
const
;
331
335
void
imode
(
int
in) {_imode=in;}
336
340
void
ME
(DecayMEPtr in)
const
{ _matrixelement = in;}
341
348
void
resetIntermediate(
tcPDPtr
part, Energy mass, Energy width);
349
353
unsigned
int
numberModes
()
const
{
return
_modes.size();}
354
358
tDecayPhaseSpaceModePtr mode(
unsigned
int
);
359
363
tcDecayPhaseSpaceModePtr mode(
unsigned
int
)
const
;
364
368
bool
generateIntermediates
()
const
{
return
_generateinter;}
369
373
void
generateIntermediates
(
bool
in) {_generateinter=in;}
374
380
Energy initializePhaseSpaceMode(
unsigned
int
imode,
bool
init,
bool
onShell=
false
)
const
;
381
385
void
hasOneLoopME
(
bool
in) {_virtualME=in;}
386
390
void
hasRealEmissionME
(
bool
in) {_realME=in;}
391
395
void
epsilonPS
(Energy in) {_eps=in;}
396
400
void
clearModes
() {_modes.clear();}
401
402
protected
:
403
409
virtual
void
doinitrun();
411
412
private
:
413
417
DecayIntegrator
& operator=(
const
DecayIntegrator
&) =
delete
;
418
419
private
:
420
424
int
_niter
;
425
429
int
_npoint
;
430
434
int
_ntry
;
435
439
mutable
vector<DecayPhaseSpaceModePtr>
_modes
;
440
444
bool
_generateinter
;
445
449
DecayRadiationGeneratorPtr
_photongen
;
450
454
mutable
int
_imode
;
455
459
mutable
DecayMEPtr
_matrixelement
;
460
464
bool
_realME
;
465
469
bool
_virtualME
;
470
474
Energy
_eps
;
475
476
};
480
ostream &
operator<<
(ostream &,
const
DecayIntegrator
&);
481
485
class
DecayIntegratorError
:
public
Exception
{};
486
487
}
488
489
490
#endif
/* HERWIG_DecayIntegrator_H */
ThePEG::PersistentIStream
Herwig::DecayIntegrator::_generateinter
bool _generateinter
Whether to include the intermediates whne outputing the results.
Definition:
DecayIntegrator.h:444
Herwig::DecayIntegrator::epsilonPS
void epsilonPS(Energy in)
Set the epsilon parameter.
Definition:
DecayIntegrator.h:395
Herwig::DecayIntegrator::_modes
vector< DecayPhaseSpaceModePtr > _modes
List of the decay modes.
Definition:
DecayIntegrator.h:439
ThePEG::Particle
Herwig::DecayIntegrator::hasOneLoopME
bool hasOneLoopME()
Whether or not the one loop matrix element is implemented.
Definition:
DecayIntegrator.h:267
Herwig::DecayIntegratorError
Exception for this class and those inheriting from it.
Definition:
DecayIntegrator.h:485
Herwig::DecayIntegrator::generateIntermediates
void generateIntermediates(bool in)
Set whether or not the intermediates are included.
Definition:
DecayIntegrator.h:373
Herwig::DecayIntegrator::generateIntermediates
bool generateIntermediates() const
Get whether or not the intermediates are included.
Definition:
DecayIntegrator.h:368
Herwig::DecayPhaseSpaceMode
The DecayPhaseSpaceMode class is designed to store a group of phase-space channels for use by the Dec...
Definition:
DecayPhaseSpaceMode.h:43
Herwig::DecayIntegrator::_matrixelement
DecayMEPtr _matrixelement
The helicity matrix element for the current decay.
Definition:
DecayIntegrator.h:459
ThePEG::PersistentOStream
Herwig::DecayIntegrator::hasOneLoopME
void hasOneLoopME(bool in)
Whether or not the one loop matrix element is implemented.
Definition:
DecayIntegrator.h:385
ThePEG::LorentzRotation
Herwig::DecayIntegrator::numberModes
unsigned int numberModes() const
Number of decay modes.
Definition:
DecayIntegrator.h:353
ThePEG
Herwig::DecayIntegrator::imode
void imode(int in)
Set the mode being use for this decay.
Definition:
DecayIntegrator.h:335
Herwig::DecayIntegrator::_eps
Energy _eps
Epsilon parameter for phase-space integration.
Definition:
DecayIntegrator.h:474
Herwig::DecayIntegrator::generatePhotons
ParticleVector generatePhotons(const Particle &p, ParticleVector children)
Members for the generation of QED radiation in the decays.
Definition:
DecayIntegrator.h:244
Herwig::DecayIntegrator::MEOption
MEOption
Enum for the matrix element option.
Definition:
DecayIntegrator.h:79
Herwig::operator<<
ostream & operator<<(ostream &, const DecayIntegrator &)
Output information on the DecayIntegrator for debugging purposes.
Herwig::DecayIntegrator::clearModes
void clearModes()
Clear the models.
Definition:
DecayIntegrator.h:400
Herwig::DecayIntegrator::_virtualME
bool _virtualME
Whether or not the one-loop matrix element exists.
Definition:
DecayIntegrator.h:469
Herwig::DecayIntegrator::accept
virtual bool accept(tcPDPtr parent, const tPDVector &children) const
Check if this decayer can perfom the decay for a particular mode.
Definition:
DecayIntegrator.h:96
Herwig::DecayIntegrator::canGeneratePhotons
bool canGeneratePhotons()
check if photons can be generated in the decay
Definition:
DecayIntegrator.h:251
ThePEG::tPDVector
vector< tPDPtr > tPDVector
Herwig::DecayIntegrator::_photongen
DecayRadiationGeneratorPtr _photongen
Pointer to the object generating the QED radiation in the decay.
Definition:
DecayIntegrator.h:449
Herwig::DecayIntegrator::ME
void ME(DecayMEPtr in) const
Set the helicity matrix element for the decay.
Definition:
DecayIntegrator.h:340
Herwig::DecayIntegrator
Main class for Decayers implementing multi-channel phase space integration.
Definition:
DecayIntegrator.h:62
ThePEG::tcPDPtr
ThePEG::Ptr< ParticleData >::transient_const_pointer tcPDPtr
Herwig::DecayIntegrator::hasRealEmissionME
bool hasRealEmissionME()
Whether or not the real emission matrix element is implemented.
Definition:
DecayIntegrator.h:293
Herwig::DecayIntegrator::_imode
int _imode
mode currently being generated
Definition:
DecayIntegrator.h:454
Herwig::DecayIntegrator::epsilonPS
Energy epsilonPS() const
Access to the epsilon parameter.
Definition:
DecayIntegrator.h:229
Herwig::DecayIntegrator::hasRealEmissionME
void hasRealEmissionME(bool in)
Whether or not the real emission matrix element is implemented.
Definition:
DecayIntegrator.h:390
Herwig::DecayIntegrator::DecayIntegrator
DecayIntegrator()
Default constructor.
Definition:
DecayIntegrator.h:86
Herwig::RandomHelpers::generate
pair< double, double > generate(const Generator< Density > &gen, double r)
Generate a random variable and return its weight.
Definition:
RandomHelpers.h:182
ThePEG::Exception
Herwig::HwDecayerBase
The HwDecayerBase class is the base class for Decayers in Herwig.
Definition:
HwDecayerBase.h:37
Herwig
-*- C++ -*-
Definition:
BasicConsistency.h:17
ThePEG::DecayMode
Herwig::DecayIntegrator::_npoint
int _npoint
Number of points for initialisation.
Definition:
DecayIntegrator.h:429
Herwig::DecayIntegrator::twoBodyMEcode
virtual bool twoBodyMEcode(const DecayMode &, int &mecode, double &coupling) const
Specify the matrix element to be used in the running width calculation.
Definition:
DecayIntegrator.h:158
ThePEG::ParticleVector
vector< PPtr > ParticleVector
Herwig::DecayIntegrator::_realME
bool _realME
Whether or not the real photon emission matrix element exists.
Definition:
DecayIntegrator.h:464
Herwig::DecayIntegrator::ME
DecayMEPtr ME() const
The helicity amplitude matrix element for spin correlations.
Definition:
DecayIntegrator.h:149
Herwig::DecayIntegrator::imode
int imode() const
The mode being used for this decay.
Definition:
DecayIntegrator.h:122
Herwig::DecayIntegrator::_ntry
int _ntry
number of attempts to generate the decay
Definition:
DecayIntegrator.h:434
Herwig::DecayIntegrator::_niter
int _niter
Number of iterations for th initialization.
Definition:
DecayIntegrator.h:424
ZERO
constexpr ZeroUnit ZERO
Generated on Sat Apr 11 2020 14:50:29 for Herwig by
1.8.13