PYTHIA  8.312
Classes | Public Member Functions | List of all members
HadronWidths Class Reference

#include <HadronWidths.h>

Inheritance diagram for HadronWidths:
PhysicsBase

Public Member Functions

bool init (string path)
 Load hadron widths data from an xml file. More...
 
bool init (istream &stream)
 Initialize. More...
 
bool check ()
 Check whether input data is valid and matches particle data. More...
 
bool hasResonances (int idA, int idB) const
 Get whether the specified incoming particles can form a resonance. More...
 
set< int > getResonances () const
 Get all implemented resonances.
 
set< int > getResonances (int idA, int idB) const
 Get resonances that can be formed by the specified incoming particles. More...
 
bool hasData (int id) const
 Returns whether the specified particle is handled by HadronWidths.
 
bool canDecay (int id, int prodA, int prodB) const
 Get whether the resonance can decay into the specified products. More...
 
double width (int id, double m) const
 Get the total width of the specified particle at the specified mass. More...
 
double partialWidth (int idR, int prodA, int prodB, double m) const
 Get the partial width for the specified decay channel of the particle. More...
 
double br (int idR, int prodA, int prodB, double m) const
 Get the branching ratio for the specified decay channel of the particle. More...
 
double mDistr (int id, double m) const
 Get the mass distribution density for the particle at the specified mass. More...
 
bool pickMasses (int idA, int idB, double eCM, double &mAOut, double &mBOut, int lType=1)
 Sample masses for the outgoing system with a given eCM. More...
 
bool pickDecay (int idDec, double m, int &idAOut, int &idBOut, double &mAOut, double &mBOut)
 
double widthCalc (int id, double m) const
 Calculate the total width of the particle without using interpolation. More...
 
double widthCalc (int id, int prodA, int prodB, double m) const
 Calculate partial width of the particle without using interpolation. More...
 
bool parameterize (int id, int precision=50)
 Regenerate parameterization for particle and its decay products if needed. More...
 
void parameterizeAll (int precision=50)
 Regenerate parameterization for all particles. More...
 
bool save (ostream &stream) const
 Write all width data to an xml file. More...
 
bool save (string file="HadronWidths.dat") const
 
- Public Member Functions inherited from PhysicsBase
void initInfoPtr (Info &infoPtrIn)
 This function is called from above for physics objects used in a run. More...
 
virtual ~PhysicsBase ()
 Empty virtual destructor.
 
bool flag (string key) const
 Shorthand to read settings values.
 
int mode (string key) const
 
double parm (string key) const
 
string word (string key) const
 
vector< bool > fvec (string key) const
 
vector< int > mvec (string key) const
 
vector< double > pvec (string key) const
 
vector< string > wvec (string key) const
 

Additional Inherited Members

- Public Types inherited from PhysicsBase
enum  Status {
  INCOMPLETE = -1, COMPLETE = 0, CONSTRUCTOR_FAILED, INIT_FAILED,
  LHEF_END, LOWENERGY_FAILED, PROCESSLEVEL_FAILED, PROCESSLEVEL_USERVETO,
  MERGING_FAILED, PARTONLEVEL_FAILED, PARTONLEVEL_USERVETO, HADRONLEVEL_FAILED,
  CHECK_FAILED, OTHER_UNPHYSICAL, HEAVYION_FAILED, HADRONLEVEL_USERVETO
}
 Enumerate the different status codes the event generation can have.
 
- Protected Member Functions inherited from PhysicsBase
 PhysicsBase ()
 Default constructor.
 
virtual void onInitInfoPtr ()
 
virtual void onBeginEvent ()
 This function is called in the very beginning of each Pythia::next call.
 
virtual void onEndEvent (Status)
 
virtual void onStat ()
 This function is called from the Pythia::stat() call.
 
void registerSubObject (PhysicsBase &pb)
 Register a sub object that should have its information in sync with this.
 
- Protected Attributes inherited from PhysicsBase
InfoinfoPtr = {}
 
SettingssettingsPtr = {}
 Pointer to the settings database.
 
ParticleDataparticleDataPtr = {}
 Pointer to the particle data table.
 
LoggerloggerPtr = {}
 Pointer to logger.
 
HadronWidthshadronWidthsPtr = {}
 Pointer to the hadron widths data table.
 
RndmrndmPtr = {}
 Pointer to the random number generator.
 
CoupSMcoupSMPtr = {}
 Pointers to SM and SUSY couplings.
 
CoupSUSYcoupSUSYPtr = {}
 
BeamSetupbeamSetupPtr = {}
 
BeamParticlebeamAPtr = {}
 
BeamParticlebeamBPtr = {}
 
BeamParticlebeamPomAPtr = {}
 
BeamParticlebeamPomBPtr = {}
 
BeamParticlebeamGamAPtr = {}
 
BeamParticlebeamGamBPtr = {}
 
BeamParticlebeamVMDAPtr = {}
 
BeamParticlebeamVMDBPtr = {}
 
PartonSystemspartonSystemsPtr = {}
 Pointer to information on subcollision parton locations.
 
SigmaTotalsigmaTotPtr = {}
 Pointers to the total/elastic/diffractive cross sections.
 
SigmaCombinedsigmaCmbPtr = {}
 
set< PhysicsBase * > subObjects
 
UserHooksPtr userHooksPtr
 

Detailed Description

The HadronWidths class is used to compute mass-dependent widths and branching ratios of hadrons.

Member Function Documentation

double br ( int  idR,
int  prodA,
int  prodB,
double  m 
) const

Get the branching ratio for the specified decay channel of the particle.

Get key and flip idR if necessary.

Find particle data entry.

Check that mass is within range.

If particle is not resonance, get BR from particle database.

If resonance cannot decay into A B, return 0.

For resonances, get branching ratio from parameterization.

bool canDecay ( int  id,
int  prodA,
int  prodB 
) const

Get whether the resonance can decay into the specified products.

Get key and flip idR if necessary.

bool check ( )

Check whether input data is valid and matches particle data.

Check that all resonance entries make sense.

Check that entry id actually corresponds to a particle.

Check that entry id is positive (antiparticles are handled by symmetry).

Check that entry id is hadron.

Check that mass boundaries are same in particle entry and widths entry.

Check that all decay channels make sense.

Check that decay product ids actually correspond to particles.

Check that lType makes sense.

Check that decay conserves charge.

set< int > getResonances ( int  idA,
int  idB 
) const

Get resonances that can be formed by the specified incoming particles.

Get signature for system and look only for resonances that matches it.

For resonances that matches signature, check that decay channel exists.

For pi0pi0 and pi+pi-, add f0(500) explicitly.

Done.

bool hasResonances ( int  idA,
int  idB 
) const

Get whether the specified incoming particles can form a resonance.

Get signature for system and look only for resonances that matches it.

For resonances that matches signature, check that decay channel exists.

No resonances found.

bool init ( string  path)

Load hadron widths data from an xml file.

The HadronWidths class.

Initialize.

bool init ( istream &  stream)

Initialize.

Insert resonance in entries

Insert resonance in signature map

If signature has not been used yet, insert a new vector into the map

If signature has been used already, add id to the existing vector

Generate key to ensure canonical ordering of decay products.

Insert new decay channel.

Search for newly added particles.

Parameterize particle if it is not yet listed and has variable width.

Check that has a resonance decay channel, i.e. into two hadrons.

Done.

double mDistr ( int  id,
double  m 
) const

Get the mass distribution density for the particle at the specified mass.

Get the mass distribution density for the particle at the specified mass, using the Breit-Wigner formula with a mass-dependent width.

bool parameterize ( int  id,
int  precision = 50 
)

Regenerate parameterization for particle and its decay products if needed.

Regenerate parameterization for particle, using the specified number of interpolation points. If needed, its decay products are automatically parameterized as well.

Get particle entry and validate input.

Do recursive parameterization (no validation checks in recursive calls)

void parameterizeAll ( int  precision = 50)

Regenerate parameterization for all particles.

Get all particles with varWidth from particle database.

Clear existing data and parameterize new data.

double partialWidth ( int  idR,
int  prodA,
int  prodB,
double  m 
) const

Get the partial width for the specified decay channel of the particle.

Get key and flip idR if necessary.

Find particle data entry.

Check that mass is within range.

If particle is not resonance, use Breit-Wigner width.

For resonances, get width from parameterization.

bool pickDecay ( int  idDec,
double  m,
int &  idAOut,
int &  idBOut,
double &  mAOut,
double &  mBOut 
)

Pick a decay channel for the specified particle, together with phase space configuration. Returns whether successful.

Pick a decay channel for the specified particle, together with phase space configuration. If successful, the results are written to the output arguments.

Find particle data entry for decaying particle.

If antiparticle, flip the resonance and then later flip decay products.

Find table entry for decaying particle.

Get list of channels that are currently on.

Check that channel is open at this mass.

Check that channel is on in the particle database.

If channel is on, add to list if width is positive.

If no channels are on, the decay fails.

Pick masses of decay products.

Write output values and done.

bool pickMasses ( int  idA,
int  idB,
double  eCM,
double &  mAOut,
double &  mBOut,
int  lType = 1 
)

Sample masses for the outgoing system with a given eCM.

Pick a pair of masses given pair invariant mass and angular momentum.

Minimal masses must be a possible choice.

Done if none of the daughters have a width.

Get width entries for particles with mass-depedent widths.

Parameters for mass selection.

Loop over attempts to pick the two masses simultaneously.

Simplify handling if full procedure does not seem to work.

Initially pick according to simple Breit-Wigner.

Correction given by BW(Gamma_now)/BW(Gamma_fix) for variable width. Note: width not allowed to explode at large masses.

Weight by (p/p_max)^lType.

Give warning message only for more severe cases

Last resort: pick masses within limits known to work, without weight.

Done.

bool save ( ostream &  stream) const

Write all width data to an xml file.

Write all widths data to an xml file.

Counter for number of entries on current line, maximum 8 per line.

Write total width.

Write partial widths.

Done.

double width ( int  id,
double  m 
) const

Get the total width of the specified particle at the specified mass.

Find particle data entry.

Check that mass is within range.

If particle is not resonance, use Breit-Wigner width.

For resonances, get width from parameterization.

double widthCalc ( int  id,
double  m 
) const

Calculate the total width of the particle without using interpolation.

Find particle data entry.

Sum contributions from all channels.

double widthCalc ( int  id,
int  prodA,
int  prodB,
double  m 
) const

Calculate partial width of the particle without using interpolation.

Find particle data entry.

Search for the matching decay channel.

Decay channel not found.


The documentation for this class was generated from the following files: