1 // ParticleDecays.h is a part of the PYTHIA event generator.
2 // Copyright (C) 2008 Torbjorn Sjostrand.
3 // PYTHIA is licenced under the GNU GPL version 2, see COPYING for details.
4 // Please respect the MCnet Guidelines, see GUIDELINES for details.
6 // This file contains the classes to perform a particle decay.
7 // DecayHandler: base class for external handling of decays.
8 // ParticleDecays: decay a particle.
10 #ifndef Pythia8_ParticleDecays_H
11 #define Pythia8_ParticleDecays_H
15 #include "FragmentationFlavZpT.h"
17 #include "ParticleData.h"
18 #include "PythiaStdlib.h"
20 #include "TimeShower.h"
24 //**************************************************************************
26 // DecayHandler is base class for the external handling of decays.
27 // There is only one pure virtual method, that should do the decay.
33 // A pure virtual method, wherein the derived class method does a decay.
34 virtual bool decay(vector<int>& idProd, vector<double>& mProd,
35 vector<Vec4>& pProd, int iDec, const Event& event) = 0;
40 virtual ~DecayHandler() {}
44 //**************************************************************************
46 // The ParticleDecays class contains the routines to decay a particle.
48 class ParticleDecays {
55 // Initialize: store pointers and find settings
56 void init(Info* infoPtrIn, TimeShower* timesDecPtrIn,
57 StringFlav* flavSelPtrIn, DecayHandler* decayHandlePtrIn,
58 vector<int> handledParticles);
60 // Perform a decay of a single particle.
61 bool decay(int iDec, Event& event);
63 // Did decay result in new partons to hadronize?
64 bool moreToDo() const {return hasPartons && keepPartons;}
68 // Constants: could only be changed in the code itself.
69 static const int NTRYDECAY, NTRYPICK;
70 static const double MSAFEDALITZ, WTCORRECTION[11];
72 // Pointer to various information on the generation.
75 // Pointers to timelike showers, for decays to partons (e.g. Upsilon).
76 TimeShower* timesDecPtr;
78 // Pointer to class for flavour generation; needed when to pick hadrons.
79 StringFlav* flavSelPtr;
81 // Pointer to a handler of external decays.
82 DecayHandler* decayHandlePtr;
84 // Initialization data, read from Settings..
85 bool limitTau0, limitTau, limitRadius, limitCylinder, limitDecay,
87 double mSafety, tau0Max, tauMax, rMax, xyMax, zMax, xBdMix, xBsMix,
88 sigmaSoft, multIncrease, multRefMass, multGoffset, colRearrange,
89 stopMass, sRhoDal, wRhoDal;
91 // Multiplicity. Decay products positions and masses.
92 bool hasPartons, keepPartons;
93 int idDec, meMode, mult;
94 vector<int> iProd, idProd, cols, acols, idPartons;
95 vector<double> mProd, mInv, rndmOrd;
96 vector<Vec4> pInv, pProd;
97 vector<FlavContainer> flavEnds;
99 // Pointer to particle data for currently decaying particle
100 ParticleDataEntry* decDataPtr;
102 // Check whether a decay is allowed, given the upcoming decay vertex.
103 bool checkVertex(Particle& decayer);
105 // Check for oscillations B0 <-> B0bar or B_s0 <-> B_s0bar.
106 bool oscillateB(Particle& decayer);
108 // Do a one-body decay.
109 bool oneBody(Event& event);
111 // Do a two-body decay;
112 bool twoBody(Event& event);
114 // Do a three-body decay;
115 bool threeBody(Event& event);
117 // Do a multibody decay using the M-generator algorithm.
118 bool mGenerator(Event& event);
120 // Select mass of lepton pair in a Dalitz decay.
123 // Do kinematics of gamma* -> l- l+ in Dalitz decay.
124 bool dalitzKinematics(Event& event);
126 // Translate a partonic content into a set of actual hadrons.
129 // Set colour flow and scale in a decay explicitly to partons.
130 bool setColours(Event& event);
134 //**************************************************************************
136 } // end namespace Pythia8
138 #endif // Pythia8_ParticleDecays_H