PYTHIA  8.312
Public Member Functions | Protected Member Functions | Protected Attributes | Static Protected Attributes | List of all members
PhaseSpace Class Referenceabstract

#include <PhaseSpace.h>

Inheritance diagram for PhaseSpace:
PhysicsBase PhaseSpace2to1tauy PhaseSpace2to2diffractive PhaseSpace2to2elastic PhaseSpace2to2nondiffractive PhaseSpace2to2tauyz PhaseSpace2to3diffractive PhaseSpace2to3tauycyl PhaseSpace2to3yyycyl PhaseSpaceLHA

Public Member Functions

virtual ~PhaseSpace ()
 Destructor.
 
void init (bool isFirst, SigmaProcessPtr sigmaProcessPtrIn)
 Perform simple initialization and store pointers. More...
 
void updateBeamIDs ()
 Switch to new beam particle identities; for similar hadrons only.
 
void newECM (double eCMin)
 Update the CM energy of the event.
 
void setLHAPtr (LHAupPtr lhaUpPtrIn)
 Store or replace Les Houches pointer.
 
virtual bool setupSampling ()=0
 
virtual bool trialKin (bool inEvent=true, bool repeatSame=false)=0
 
virtual bool finalKin ()=0
 
void decayKinematics (Event &process)
 Allow for nonisotropic decays when ME's available. More...
 
double sigmaNow () const
 Give back current or maximum cross section, or set latter.
 
double sigmaMax () const
 
double biasSelectionWeight () const
 
bool newSigmaMax () const
 
void setSigmaMax (double sigmaMaxIn)
 
double sigmaMaxSwitch ()
 New value for switched beam identity or energy (for SoftQCD processes).
 
virtual double sigmaSumSigned () const
 For Les Houches with negative event weight needs.
 
Vec4 p (int i) const
 Give back constructed four-vectors and known masses.
 
double m (int i) const
 
void setP (int i, Vec4 pNew)
 Reset the four-momentum.
 
double ecm () const
 Give back other event properties.
 
double x1 () const
 
double x2 () const
 
double sHat () const
 
double tHat () const
 
double uHat () const
 
double pTHat () const
 
double thetaHat () const
 
double phiHat () const
 
double runBW3 () const
 
double runBW4 () const
 
double runBW5 () const
 
virtual bool isResolved () const
 Inform whether beam particles are resolved in partons or scatter directly.
 
virtual void rescaleSigma (double)
 
virtual void rescaleMomenta (double)
 
virtual double weightGammaPDFApprox ()
 Calculate the weight for over-estimated cross section.
 
virtual void setGammaKinPtr (GammaKinematics *gammaKinPtrIn)
 Set the GammaKinematics pointer needed for soft photoproduction.
 
- 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
 

Protected Member Functions

 PhaseSpace ()
 Constructor.
 
void decayKinematicsStep (Event &process, int iRes)
 Reselect decay products momenta isotropically in phase space. More...
 
void setup3Body ()
 Much common code for normal 2 -> 1, 2 -> 2 and 2 -> 3 cases: More...
 
bool setupSampling123 (bool is2, bool is3)
 Determine how phase space should be sampled. More...
 
bool trialKin123 (bool is2, bool is3, bool inEvent=true)
 Select a trial kinematics phase space point. More...
 
bool limitTau (bool is2, bool is3)
 Calculate kinematical limits for 2 -> 1/2/3. More...
 
bool limitY ()
 Find range of allowed y values. More...
 
bool limitZ ()
 Find range of allowed z = cos(theta) values. More...
 
void selectTau (int iTau, double tauVal, bool is2)
 Select kinematical variable between defined limits for 2 -> 1/2/3. More...
 
void selectY (int iY, double yVal)
 Select y according to a choice of shapes. More...
 
void selectZ (int iZ, double zVal)
 
bool select3Body ()
 
void solveSys (int n, int bin[8], double vec[8], double mat[8][8], double coef[8])
 Solve equation system for better phase space coefficients in 2 -> 1/2/3. More...
 
void setupMass1 (int iM)
 Setup mass selection for one resonance at a time. Split in two parts. More...
 
void setupMass2 (int iM, double distToThresh)
 Setup mass selection for one resonance at a time - part 2. More...
 
void trialMass (int iM)
 Do mass selection and find the associated weight. More...
 
double weightMass (int iM)
 
pair< double, double > tRange (double sIn, double s1In, double s2In, double s3In, double s4In)
 
bool tInRange (double tIn, double sIn, double s1In, double s2In, double s3In, double s4In)
 
- 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

SigmaProcessPtr sigmaProcessPtr
 Pointer to cross section.
 
LHAupPtr lhaUpPtr
 Pointer to LHAup for generating external events.
 
GammaKinematicsgammaKinPtr
 Pointer to object that samples photon kinematics from leptons.
 
bool useBreitWigners
 Initialization data, normally only set once.
 
bool doEnergySpread
 
bool showSearch
 
bool showViolation
 
bool increaseMaximum
 
bool hasQ2Min
 
int gmZmodeGlobal
 
double mHatGlobalMin
 
double mHatGlobalMax
 
double pTHatGlobalMin
 
double pTHatGlobalMax
 
double Q2GlobalMin
 
double pTHatMinDiverge
 
double minWidthBreitWigners
 
double minWidthNarrowBW
 
int idA
 Information on incoming beams.
 
int idB
 
int idAold
 
int idBold
 
int idAgm
 
int idBgm
 
double mA
 
double mB
 
double eCM
 
double s
 
double sigmaMxGm
 
bool hasLeptonBeamA
 
bool hasLeptonBeamB
 
bool hasOneLeptonBeam
 
bool hasTwoLeptonBeams
 
bool hasPointGammaA
 
bool hasPointGammaB
 
bool hasOnePointParticle
 
bool hasTwoPointParticles
 
bool hasGamma
 
bool hasVMD
 
bool newSigmaMx
 Cross section information.
 
bool canModifySigma
 
bool canBiasSelection
 
bool canBias2Sel
 
int gmZmode
 
double bias2SelPow
 
double bias2SelRef
 
double wtBW
 
double sigmaNw
 
double sigmaMx
 
double sigmaPos
 
double sigmaNeg
 
double biasWt
 
double mHatMin
 Process-specific kinematics properties, almost always available.
 
double mHatMax
 
double sHatMin
 
double sHatMax
 
double pTHatMin
 
double pTHatMax
 
double pT2HatMin
 
double pT2HatMax
 
double x1H
 Event-specific kinematics properties, almost always available.
 
double x2H
 
double m3
 
double m4
 
double m5
 
double s3
 
double s4
 
double s5
 
double mHat
 
double sH
 
double tH
 
double uH
 
double pAbs
 
double p2Abs
 
double pTH
 
double theta
 
double phi
 
double betaZ
 
Vec4 pH [12]
 
double mH [12]
 
int idResA
 Presence and properties of any s-channel resonances.
 
int idResB
 
double mResA
 
double mResB
 
double GammaResA
 
double GammaResB
 
double tauResA
 
double tauResB
 
double widResA
 
double widResB
 
bool sameResMass
 
bool useMirrorWeight
 Kinematics properties specific to 2 -> 1/2/3.
 
bool hasNegZ
 
bool hasPosZ
 
double tau
 
double y
 
double z
 
double tauMin
 
double tauMax
 
double yMax
 
double zMin
 
double zMax
 
double ratio34
 
double unity34
 
double zNeg
 
double zPos
 
double wtTau
 
double wtY
 
double wtZ
 
double wt3Body
 
double runBW3H
 
double runBW4H
 
double runBW5H
 
double intTau0
 
double intTau1
 
double intTau2
 
double intTau3
 
double intTau4
 
double intTau5
 
double intTau6
 
double intY0
 
double intY12
 
double intY34
 
double intY56
 
double mTchan1
 
double sTchan1
 
double mTchan2
 
double sTchan2
 
double frac3Flat
 
double frac3Pow1
 
double frac3Pow2
 
double zNegMin
 
double zNegMax
 
double zPosMin
 
double zPosMax
 
Vec4 p3cm
 
Vec4 p4cm
 
Vec4 p5cm
 
int nTau
 Coefficients for optimized selection in 2 -> 1/2/3.
 
int nY
 
int nZ
 
double tauCoef [8]
 
double yCoef [8]
 
double zCoef [8]
 
double tauCoefSum [8]
 
double yCoefSum [8]
 
double zCoefSum [8]
 
bool useBW [6]
 Properties specific to resonance mass selection in 2 -> 2 and 2 -> 3.
 
bool useNarrowBW [6]
 
int idMass [6]
 
double mPeak [6]
 
double sPeak [6]
 
double mWidth [6]
 
double mMin [6]
 
double mMax [6]
 
double mw [6]
 
double wmRat [6]
 
double mLower [6]
 
double mUpper [6]
 
double sLower [6]
 
double sUpper [6]
 
double fracFlatS [6]
 
double fracFlatM [6]
 
double fracInv [6]
 
double fracInv2 [6]
 
double atanLower [6]
 
double atanUpper [6]
 
double intBW [6]
 
double intFlatS [6]
 
double intFlatM [6]
 
double intInv [6]
 
double intInv2 [6]
 
- 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
 

Static Protected Attributes

static const int NMAXTRY = 2
 Constants: could only be changed in the code itself. More...
 
static const int NTRY3BODY = 20
 Number of three-body trials in phase space optimization.
 
static const double SAFETYMARGIN = 1.05
 Maximum cross section increase, just in case true maximum not found.
 
static const double TINY = 1e-20
 Small number to avoid division by zero.
 
static const double EVENFRAC = 0.4
 Fraction of total weight that is shared evenly between all shapes.
 
static const double SAMESIGMA = 1e-6
 Two cross sections with a small relative error are assumed same.
 
static const double MRESMINABS = 0.001
 Do not allow resonance to have mass below 2 m_e.
 
static const double WIDTHMARGIN = 20.
 Do not include resonances peaked too far outside allowed mass region.
 
static const double SAMEMASS = 0.01
 Special optimization treatment when two resonances at almost same mass.
 
static const double MASSMARGIN = 0.01
 Minimum phase space left when kinematics constraints are combined.
 
static const double EXTRABWWTMAX = 1.25
 When using Breit-Wigners in 2 -> 2 raise maximum weight estimate.
 
static const double THRESHOLDSIZE = 3.
 Size of Breit-Wigner threshold region, for mass selection biasing.
 
static const double THRESHOLDSTEP = 0.2
 Step size in optimal-mass search, for mass selection biasing.
 
static const double YRANGEMARGIN = 1e-6
 Minimal rapidity range for allowed open range (in 2 -> 3).
 
static const double LEPTONXMIN = 1e-10
 
static const double LEPTONXMAX = 1. - 1e-10
 
static const double LEPTONXLOGMIN = log(1e-10)
 
static const double LEPTONXLOGMAX = log(1. - 1e-10)
 
static const double LEPTONTAUMIN = 2e-10
 
static const double SHATMINZ = 1.
 Safety to avoid division with unreasonably small value for z selection.
 
static const double PT2RATMINZ = 0.0001
 Regularization for small pT2min in z = cos(theta) selection.
 
static const double WTCORRECTION [11]
 

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.
 

Detailed Description

PhaseSpace is a base class for phase space generators used in the selection of hard-process kinematics.

Member Function Documentation

void decayKinematics ( Event process)

Allow for nonisotropic decays when ME's available.

Identify sets of sister partons.

Check that at least one of them is a resonance.

Evaluate matrix element and decide whether to keep kinematics.

Find resonances for which to redo decay angles.

Do decay of this mother isotropically in phase space.

End loop over resonance decay chains.

Ready to allow new test of matrix element.

End loop over sets of sister resonances/partons.

void decayKinematicsStep ( Event process,
int  iRes 
)
protected

Reselect decay products momenta isotropically in phase space.

Reselect decay products momenta isotropically in phase space. Does not redo secondary vertex position!

Multiplicity and mother mass and four-momentum.

Description of two-body decays as simple special case.

Products and product masses.

Fill four-momenta in mother rest frame and then boost to lab frame.

Done for two-body decay.

Description of three-body decays as semi-simple special case.

Products and product masses.

Kinematical limits for 2+3 mass. Maximum phase-space weight.

Pick an intermediate mass m23 flat in the allowed range.

Translate into relative momenta and find phase-space weight.

If rejected, try again with new invariant masses.

Set up m23 -> m2 + m3 isotropic in its rest frame.

Set up 0 -> 1 + 23 isotropic in its rest frame.

Boost 2 + 3 to the 0 rest frame.

Boost from rest frame to lab frame.

Done for three-body decay.

Do a multibody decay using the M-generator algorithm.

Set up masses and four-momenta in a vector, with mother in slot 0.

Sum of daughter masses.

Begin setup of intermediate invariant masses.

Calculate the maximum weight in the decay.

Begin loop to find the set of intermediate invariant masses.

Find and order random numbers in descending order.

Translate into intermediate masses and find weight.

If rejected, try again with new invariant masses.

Perform two-particle decays in the respective rest frame.

Boost decay products to the mother rest frame and on to lab frame.

Done for multibody decay.

virtual bool finalKin ( )
pure virtual

A pure virtual method, wherein the accepted event kinematics is to be constructed in the derived class.

Implemented in PhaseSpaceLHA, PhaseSpace2to3yyycyl, PhaseSpace2to3tauycyl, PhaseSpace2to2nondiffractive, PhaseSpace2to3diffractive, PhaseSpace2to2diffractive, PhaseSpace2to2elastic, PhaseSpace2to2tauyz, and PhaseSpace2to1tauy.

void init ( bool  isFirst,
SigmaProcessPtr  sigmaProcessPtrIn 
)

Perform simple initialization and store pointers.

Store input pointers for future use.

Some commonly used beam information.

Flag if lepton beams, and if non-resolved ones.

Flags also for unresolved photons.

Flag if photons from leptons.

Set flags for (un)resolved photons according to gammaModes.

Standard phase space cuts.

Optionally separate phase space cuts for second hard process.

Cutoff against divergences at pT -> 0.

Special cut on DIS Q2 = -tHat.

For photons from lepton beams match the cuts to gm+gm system cuts.

When to use Breit-Wigners.

Whether generation is with variable energy.

Flags for maximization information and violation handling.

Know whether a Z0 is pure Z0 or admixed with gamma*.

Flags if user should be allowed to reweight cross section.

Parameters for simplified reweighting of 2 -> 2 processes.

Default event-specific kinematics properties.

Default cross section information.

bool limitTau ( bool  is2,
bool  is3 
)
protected

Calculate kinematical limits for 2 -> 1/2/3.

Find range of allowed tau values.

Trivial reply for unresolved lepton beams.

Requirements from allowed mHat range and allowed Q2Min.

Requirements from allowed pT range and masses.

Check that there is an open range.

bool limitY ( )
protected

Find range of allowed y values.

Trivial reply for unresolved lepton beams.

Requirements from selected tau value. Trivial for one unresolved beam.

For lepton beams requirements from cutoff for f_e^e.

Check that there is an open range.

bool limitZ ( )
protected

Find range of allowed z = cos(theta) values.

Default limits.

Requirements from pTHat limits.

Check that there is an open range so far.

Define two individual ranges.

Optionally introduce Q2 = -tHat cut.

Check that there is an open range.

virtual void rescaleSigma ( double  )
inlinevirtual

Functions to rescale momenta and cross section for new sHat Currently implemented only for PhaseSpace2to2tauyz class.

Reimplemented in PhaseSpace2to2tauyz.

bool select3Body ( )
protected

Select three-body phase space according to a cylindrically based form that can be chosen to favour low pT based on the form of propagators.

Upper and lower limits of pT choice for 4 and 5.

Check that pT ranges not closed.

Select pT4S according to c0 + c1/(M^2 + pT^2) + c2/(M^2 + pT^2)^2.

Select pT5S according to c0 + c1/(M^2 + pT^2) + c2/(M^2 + pT^2)^2.

Select azimuthal angles and check that third pT in range.

Calculate transverse masses and check that phase space not closed.

Select rapidity for particle 3 and check that phase space not closed.

Find momentum transfers in the two mirror solutions (in 4-5 frame).

Construct relative mirror weights and make choice.

Construct four-vectors in rest frame of subprocess.

Total weight to associate with kinematics choice.

Cross section of form |M|^2/(2 sHat) dPS_3 so need 1/(2 sHat).

Done.

void selectTau ( int  iTau,
double  tauVal,
bool  is2 
)
protected

Select kinematical variable between defined limits for 2 -> 1/2/3.

Select tau according to a choice of shapes.

Trivial reply for unresolved lepton beams.

Contributions from s-channel resonances.

Contributions from 1 / (1 - tau) for lepton beams.

Select according to 1/tau or 1/tau^2.

Select according to 1 / (1 - tau) for lepton beams.

Select according to 1 / (tau * (tau + tauRes)) or 1 / ((tau - tauRes)^2 + widRes^2) for resonances A and B.

Phase space weight in tau.

Calculate sHat and absolute momentum of outgoing partons.

void selectY ( int  iY,
double  yVal 
)
protected

Select y according to a choice of shapes.

Trivial reply for two unresolved lepton beams.

Trivial replies for one unresolved lepton beam.

For lepton beams skip options 3&4 and go straight to 5&6.

Standard expressions used below.

1 / cosh(y).

y - y_min or mirrored y_max - y.

exp(y) or mirrored exp(-y).

1 / (1 - exp(y - y_max)) or mirrored 1 / (1 - exp(y_min - y)).

Mirror two cases.

Phase space integral in y.

Calculate x1 and x2.

void selectZ ( int  iZ,
double  zVal 
)
protected

Select z = cos(theta) according to a choice of shapes. The selection is split in the positive- and negative-z regions, since a pTmax cut can remove the region around z = 0. Furthermore, a Q2 (= -tHat) cut can make the two regions asymmetric.

Mass-dependent dampening of pT -> 0 limit.

Common expressions of unity - z and unity + z limits, protected from 0.

Evaluate integrals over negative and positive z ranges. Flat in z.

1 / (unity34 - z).

1 / (unity34 + z).

1 / (unity34 - z)^2.

1 / (unity34 + z)^2.

Pick z value according to alternatives. Flat in z.

1 / (unity34 - z).

1 / (unity34 + z).

1 / (unity34 - z)^2.

1 / (unity34 + z)^2.

Safety check for roundoff errors. Combinations with z.

Phase space integral in z.

Calculate tHat and uHat. Also gives pTHat.

void setup3Body ( )
protected

Much common code for normal 2 -> 1, 2 -> 2 and 2 -> 3 cases:

Determine how 3-body phase space should be sampled.

Determine how phase space should be sampled.

Check for massive t-channel propagator particles.

Find coefficients of different pT2 selection terms. Mirror choice.

void setupMass1 ( int  iM)
protected

Setup mass selection for one resonance at a time. Split in two parts.

Setup mass selection for one resonance at a time - part 1.

Identity for mass seletion; is 0 also for light quarks (not yet selected).

Masses and widths of resonances.

gmZmode == 1 means pure photon propagator; set at lower mass limit.

Mass and width combinations for Breit-Wigners.

Simple Breit-Wigner range, upper edge to be corrected subsequently.

void setupMass2 ( int  iM,
double  distToThresh 
)
protected

Setup mass selection for one resonance at a time - part 2.

Store reduced Breit-Wigner range.

Prepare to select m3 by BW + flat + 1/s_3. Determine relative coefficients by allowed mass range.

For gamma*/Z0: increase 1/s_i part and introduce 1/s_i^2 part.

Normalization integrals for the respective contribution.

virtual bool setupSampling ( )
pure virtual

A pure virtual method, wherein an optimization procedure is used to determine how phase space should be sampled.

Implemented in PhaseSpaceLHA, PhaseSpace2to3yyycyl, PhaseSpace2to3tauycyl, PhaseSpace2to2nondiffractive, PhaseSpace2to3diffractive, PhaseSpace2to2diffractive, PhaseSpace2to2elastic, PhaseSpace2to2tauyz, and PhaseSpace2to1tauy.

bool setupSampling123 ( bool  is2,
bool  is3 
)
protected

Determine how phase space should be sampled.

Optional printout.

Check that open range in tau (+ set tauMin, tauMax).

Reset coefficients and matrices of equation system to solve.

Number of used coefficients/points for each dimension: tau, y, c.

Identify if any resonances contribute in s-channel.

Check resonances have non-zero widths.

More sampling in tau if resonances in s-channel.

More sampling in tau (and different in y) if incoming lepton beams.

Special case when both resonances have same mass.

Default z value and weight required for 2 -> 1. Number of dimensions.

Initial values, to be modified later.

Step through grid in tau. Set limits on y and z generation.

Step through grids in y and z.

2 -> 1: calculate cross section, weighted by phase-space volume.

2 -> 2: calculate cross section, weighted by phase-space volume and Breit-Wigners for masses

2 -> 3: repeat internal 3-body phase space several times and keep maximal cross section, weighted by phase-space volume and Breit-Wigners for masses

Allow possibility for user to modify cross section. (3body??)

Check if current maximum exceeded.

Optional printout. Protect against negative cross sections.

Sum up tau cross-section pieces in points used.

Sum up y cross-section pieces in points used.

Integrals over z expressions at tauMax, to be used below.

Sum up z cross-section pieces in points used.

End of loops over phase space points.

Fail if no non-vanishing cross sections.

Solve respective equation system for better phase space coefficients.

Provide cumulative sum of coefficients.

The last element should be > 1 to be on safe side in selection below.

Begin find two most promising maxima among same points as before.

Scan same grid as before in tau, y, z.

2 -> 1: calculate cross section, weighted by phase-space volume.

2 -> 2: calculate cross section, weighted by phase-space volume and Breit-Wigners for masses

2 -> 3: repeat internal 3-body phase space several times and keep maximal cross section, weighted by phase-space volume and Breit-Wigners for masses

Allow possibility for user to modify cross section. (3body??)

Optional printout. Protect against negative cross section.

Check that point is not simply mirror of already found one.

Add to or insert in maximum list. Only first two count.

Found two most promising maxima.

Read out starting position for search.

Starting point and step size in parameter space.

Run through (possibly a subset of) tau, y and z.

Define new parameter-space point by step in one dimension.

Convert to relevant variables and find derived new limits.

Evaluate cross-section.

2 -> 1: calculate cross section, weighted by phase-space volume.

2 -> 2: calculate cross section, weighted by phase-space volume and Breit-Wigners for masses

2 -> 3: repeat internal 3-body phase space several times and keep maximal cross section, weighted by phase-space volume and Breit-Wigners for masses

Allow possibility for user to modify cross section.

Optional printout. Protect against negative cross section.

Save new maximum. Final maximum.

Optional printout.

Done.

void solveSys ( int  n,
int  bin[8],
double  vec[8],
double  mat[8][8],
double  coef[8] 
)
protected

Solve equation system for better phase space coefficients in 2 -> 1/2/3.

Solve linear equation system for better phase space coefficients.

Optional printout.

Local variables.

Check if equation system solvable.

Solve to find relative importance of cross-section pieces.

Share evenly if failure.

Normalize coefficients, with piece shared democratically.

Optional printout.

pair<double,double> tRange ( double  sIn,
double  s1In,
double  s2In,
double  s3In,
double  s4In 
)
inlineprotected

Standard methods to find t range of a 2 -> 2 process and to check whether a given t value is in that range.

virtual bool trialKin ( bool  inEvent = true,
bool  repeatSame = false 
)
pure virtual

A pure virtual method, wherein a trial event kinematics is to be selected in the derived class.

Implemented in PhaseSpaceLHA, PhaseSpace2to3yyycyl, PhaseSpace2to3tauycyl, PhaseSpace2to2nondiffractive, PhaseSpace2to3diffractive, PhaseSpace2to2diffractive, PhaseSpace2to2elastic, PhaseSpace2to2tauyz, and PhaseSpace2to1tauy.

bool trialKin123 ( bool  is2,
bool  is3,
bool  inEvent = true 
)
protected

Select a trial kinematics phase space point.

Select a trial kinematics phase space point. Note: by In is meant the integral over the quantity multiplying coefficient cn. The sum of cn is normalized to unity.

Allow for possibility that energy varies from event to event.

Find shifted tauRes values.

Choose tau according to h1(tau)/tau, where h1(tau) = c0/I0 + (c1/I1) * 1/tau

  • (c2/I2) / (tau + tauResA)
  • (c3/I3) * tau / ((tau - tauResA)^2 + widResA^2)
  • (c4/I4) / (tau + tauResB)
  • (c5/I5) * tau / ((tau - tauResB)^2 + widResB^2)
  • (c6/I6) * tau / (1 - tau).

Choose y according to h2(y), where h2(y) = (c0/I0) * 1/cosh(y)

  • (c1/I1) * (y-ymin) + (c2/I2) * (ymax-y)
  • (c3/I3) * exp(y) + (c4/i4) * exp(-y) (for hadron; for lepton instead)
  • (c5/I5) * 1 / (1 - exp(y-ymax)) + (c6/I6) * 1 / (1 - exp(ymin-y)).

Choose z = cos(thetaHat) according to h3(z), where h3(z) = c0/I0 + (c1/I1) * 1/(A - z) + (c2/I2) * 1/(A + z)

  • (c3/I3) * 1/(A - z)^2 + (c4/I4) * 1/(A + z)^2, where A = 1 + 2*(m3*m4/sH)^2 (= 1 for massless products).

2 -> 1: calculate cross section, weighted by phase-space volume.

2 -> 2: calculate cross section, weighted by phase-space volume and Breit-Wigners for masses

2 -> 3: also sample internal 3-body phase, weighted by 2 -> 1 phase-space volume and Breit-Wigners for masses

Allow possibility for user to modify cross section.

Check if maximum violated.

Violation strategy 1: increase maximum (always during initialization).

Violation strategy 2: weight event (done in ProcessContainer).

Check if negative cross section.

Optional printout of (all) violations.

Set event weight, where relevant.

Done.

void trialMass ( int  iM)
protected

Do mass selection and find the associated weight.

Select Breit-Wigner-distributed or fixed masses.

References to masses to be set.

Distribution for m_i is BW + flat(s) + 1/sqrt(s_i) + 1/s_i + 1/s_i^2.

Distribution for m_i is simple BW.

Else m_i is fixed at peak value.

double weightMass ( int  iM)
protected

Naively a fixed-width Breit-Wigner is used to pick the mass. Here come the correction factors for (i) preselection from BW + flat in s_i + 1/sqrt(s_i) + 1/s_i + 1/s_i^2, (ii) reduced allowed mass range, (iii) running width, i.e. m0*Gamma0 -> s*Gamma0/m0. In the end, the weighted distribution is a running-width BW.

Reference to mass and to Breit-Wigner weight to be set.

Default weight if no Breit-Wigner.

Weight of generated distribution.

Weight of distribution with running width in Breit-Wigner.

?? Alternative recipe, taking into account that decay channels close at different mass thresholds Needs refining, e.g. no doublecouting with openFrac and difference numerator/denominator. double mwRun = mSet * particleDataPtr->resWidthOpen(idMass[iM], mSet);

Done.

Member Data Documentation

const double LEPTONXMIN = 1e-10
staticprotected

Cutoff for f_e^e at x < 1 - 10^{-10} to be used in phase space selection. Note: the ...MIN quantities come from 1 - x_max or 1 - tau_max.

const int NMAXTRY = 2
staticprotected

Constants: could only be changed in the code itself.

Number of trial maxima around which maximum search is performed.

The PhaseSpace class. Base class for phase space generators. Constants: could be changed here if desired, but normally should not. These are of technical nature, as described for each.

const double WTCORRECTION
staticprotected
Initial value:
= { 1., 1., 1.,
2., 5., 15., 60., 250., 1250., 7000., 50000. }

These numbers are hardwired empirical parameters, intended to speed up the M-generator.


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