**********************************************************************
**********************************************************************
** **
** September 1991 **
** **
** A Manual to **
** **
** The Lund Monte Carlo for Hadronic Processes **
** **
** PYTHIA version 5.6 **
** **
** Hans-Uno Bengtsson **
** Department of Theoretical Physics, University of Uppsala **
** **
** Torbjorn Sjostrand **
** CERN/TH, CH-1211 Geneva 23 **
** BITNET/EARN address TORSJO@CERNVM **
** Tel. +22 - 767 28 20 **
** **
** Copyright Hans-Uno Bengtsson and Torbjorn Sjostrand **
** **
**********************************************************************
**********************************************************************
* *
* Table of Contents *
* *
* 1. Introductory Material *
* 1.1. Program Objective *
* 1.2. Update History *
* 1.3. Major Changes from PYTHIA 4.8 to 5.3 *
* 1.4. Major Changes from PYTHIA 5.3 to 5.4 *
* 1.5. Major Changes from PYTHIA 5.4 to 5.5 *
* 1.6. Major Changes from PYTHIA 5.5 to 5.6 *
* 1.7. Installation of Program *
* 1.8. Programming Philosophy *
* *
* 2. The Program Components *
* 2.1. The Main Subroutines *
* 2.2. The Physics Processes *
* 2.3. Comments on Physics Processes *
* 2.4. Switches for Event Type and Kinematics Selection *
* 2.5. The General Switches and Parameters *
* 2.6. General Event Information *
* 2.7. The Event Record *
* 2.8. Other Routines and Commonblocks *
* 2.9. The JETSET Routines *
* 2.10. On Cross-Sections *
* 2.11. Examples *
* *
* Acknowledgements *
* References *
* *
**********************************************************************
* Legend: *
* >= larger than or equal to *
* <= smaller than or equal to *
* /= not equal to *
* -> goes to (rightarrow) *
* ^ what follows next is to be read as an upper index *
* _ what follows next is to be read as a lower index *
* (D=..) default value for commonblock parameter *
* (R) commonblock variable which user may read but not change *
* (I) commonblock variable for purely internal use *
**********************************************************************
1. Introductory Material
PYTHIA is a program intended for the study of hadronic processes.
Historically, this has meant incoming hadron, but today also lepton
beams are included, so that e.g. e+e- and ep events can be generated,
in addition to pp ones. (Further combinations, involving antiprotons,
neutrons, pions, photons, muons, etc., are also possible.) For
completeness, it is also possible to generate a few non-hadronic
final states. In hadronic interactions, the emphasis is on high-pT
physics, but PYTHIA also covers the domain of low-pT interactions
as an integral part of the total cross-section. In its present form it
includes hard scattering matrix elements, structure functions and
initial and final state parton showers. Fragmentation is performed using
the ordinary Lund fragmentation model, JETSET version 7.3, but an
important task for PYTHIA is to set up the correct string configuration,
particularly nontrivial for the low-pT target remnants.
______________________________________________________________________
1.1. Program Objective
Leaving aside for the moment a discussion of to what extent we can
trust the various pieces that make up a full-fledged Monte Carlo
program like PYTHIA, i.e., assuming that PYTHIA indeed generates
events as they would occur if the underlying theory and assumptions
were true, three important tasks can be singled out as the main
objectives for PYTHIA:
1) to test the underlying theory by comparison with experiments at
present colliders;
2) to help designing techniques and strategies for probing the
physics at present and future colliders;
3) to help disentangle possible signals for new physics at colliders
by giving accurate estimates of standard model backgrounds and
some scenarios for physics beyond the standard model.
Over the years, PYTHIA has indeed been put to all of these uses; in
particular it has been extensively used to look at the possibilities
of detection of predicted standard model particles like the top quark
and the Higgs boson at LHC and SSC.
It should also be stated what PYTHIA is not. It is not a program for
high-precision tests of the electroweak model, like e.g. Z0 line
shape studies. The program includes an 'improved Born term' approach
(running alpha_em, s-dependent width, etc.) to resonances, and a
modified leading log description of initial and final state radiation,
which should ensure a reasonably accurate representation, but certainly
not without its limitations.
______________________________________________________________________
1.2. Update History
The PYTHIA program has undergone a contined rapid expansion, in
terms of the number of subprocesses included, and also in terms of
the physics aspects covered. This is a continuous process, with the
official numbered versions little more than snapshots of this process.
Some versions have never been given a larger distribution, but have
instead been handed to a few people for specific purposes. For the
record, we below give a list of the versions, with some brief notes.
no date publ. main new or improved features
1 Dec82 [Ben84] synthesis of predecessors COMPTON, HIGHPT and
KASSANDRA
2 -----
3.1 -----
3.2 -----
3.3 Feb84 [Ben84a] scale-breaking structure functions
3.4 Sep84 [Ben85] more efficient kinematics selection
4.1 Dec84 initial and final state parton showers, W and Z
4.2 Jun85 multiple interactions
4.3 Aug85 W W, W Z, Z Z, R
4.4 Nov85 gamma W, gamma Z, gamma gamma
4.5 Jan86 H0 production, diffractive and elastic events
4.6 May86 angular correlation in resonance pair decays
4.7 May86 Z'0, H+
4.8 Jan87 [Ben87] variable impact parameter in multiple interactions
4.9 May86 g H+
5.1 May86 massive matrix elements for heavy quarks
5.2 Jun87 intermediate boson scattering
5.3 Oct89 new particle and subprocess codes, new commonblock
structure, new kinematics selection, some
lepton-lepton and lepton-hadron interactions,
new subprocesses
5.4 Jun90 s-dependent widths, resonances not on mass-shell,
new processes, new structure functions
5.5 Jan91 improved e+e- and ep, several new processes
5.6 Sep91 [this] reorganized structure functions
Versions preceding 4.8 by now should be considered obsolete. Versions
4.9, 5.1 and 5.2 are all minor expansions on 4.8, and have been given
only limited distribution. Version 5.3 is therefore the first major
revision since 4.8. Since the current version is almost backwards
compatible with 5.3, and since several errors have been found in the
original distribution of 5.3, users are recommended to switch from
5.3 to 5.6.
______________________________________________________________________
1.3. Major Changes from PYTHIA 4.8 to 5.3
Version 5.3 of PYTHIA, the Lund Monte Carlo for Hadronic Processes,
represented a major break in continuity from a programming point
of view, with no backwards compatibility to PYTHIA 4.8, the former
standard version. Many major rewritings have been prompted by the
adoption of the new particle numbering scheme developed under the aegis
of the Particle Data Group [PDG88]. This scheme is also the standard
adopted for JETSET 7. The internal scheme for subprocess
classification has also been expanded to cover the future inclusion
of further new processes. Also a number of other changes have been
made, to systematize current features and leave space for new ones.
Here is a more extensive list of major changes in PYTHIA 5.3.
- New KF particle codes have been introduced and are consistently used
in the program.
- The program has, also in other respects, been updated and extended to
run with JETSET 7.2. In particular, this applies to the commonblock
LUJETS, which has been expanded to contain more event information.
- The ISUB subprocess classification scheme has been reorganized and
expanded.
- A number of new subprocesses have been included: heavy quark
production with massive matrix elements, gauge boson pair scattering,
etc.
- Many of the subprocesses can also be simulated in lepton-lepton
interactions (and are prepared for lepton-hadron ones).
- The possibility of event kinematics preselection has been added.
- Switches for event flavour preselection have been reorganized and
coordinated with JETSET.
- Generation of kinematical variables has been changed and systematized.
- It is possible to generate pileup events, where one single
bunch crossing gives rise to several hadron-hadron interactions.
______________________________________________________________________
1.4. Major Changes from PYTHIA 5.3 to 5.4
The updates from version 5.3 to 5.4 are all minor, and just about any
program that ran with version 5.3 will also work with PYTHIA 5.4.
The following changes might give compatibility problems from programs
based on PYTHIA 5.3:
- The routine PYTHIA has been renamed PYEVNT (to avoid confusion
between the program-as-a-whole and one specific subroutine);
a dummy routine PYTHIA which calls PYEVNT is provided, however.
- Subprocesses 142 and 143 have been moved to 143 and 144, respectively.
The related MSEL codes 22 and 23 have also been moved to 23 and 24,
respectively.
- The non-standard decay channel Z -> H+ + H- has been removed,
instead the decay Z' -> H+ + H- has been introduced to fill the
same function.
- The program should be run with JETSET version 7.3 rather than 7.2
(to take into account the new decay channels and the running of
alpha_em). Note that many decay channels appear with new numbers.
- Several errors have been corrected, where previously the program gave
wrong results.
In addition, the following changes have been made, where no
compatibility problems are involved:
- A subprocess 142 has been added; in general the Z'/W' sector has
been much expanded, and couplings to ordinary fermions and gauge
bosons appear as free parameters (defined in JETSET).
- Four new structure function parametrization by Morfin and Tung
[Mor90] come with the program, and two new by Gluck, Reya and Vogt
[GRV89].
- alpha_em is now taken to be running; if the user so decides it can
be frozen, however.
- s-channel resonances are generated with shat-dependent widths;
this should give an improved description of line shapes.
- In the decays H -> Z Z or W W, the contribution to the total
Higgs width is included also for Higgs so light that either or
both products are off mass-shell.
- The couplings of H to quarks is taken with running quark masses,
both for production q + qb -> H0 and decay H0 -> q + qb; see
MSTP(37) for possibility to use fix masses.
- In a number of 2 -> 2 processes it is now possible to select
specific mass ranges for resonances in the final state, and in
particular to generate events below the nominal threshold mass.
- Also for other decays of a resonance into two new resonances it
is possible to set mass ranges on the decay products. This feature
is not fully developed yet, however, so some limitations exist.
- Process 22, f + fb -> Z0 + Z0, has been expanded so that the two Z0
by default represent the correct mixture of gamma*/Z0, with the
possibility to switch to Z0 only or gamma* only, according to
MSTP(43) value.
- The process f + q -> f' + Q for single heavy flavour production by
W exchange in the t channel has been implemented.
- The matrix element for process 73, W_L Z_L -> W_L Z_L, has been
corrected according to formulae of [Ter90].
- Processes 71, 72, 73, 76 and 77 have been expanded so that they
alternatively can be used to simulate the strong interaction
Z_L/W_L models of Dobado, Herrero and Terron [DHT90].
- The process g + g -> Z + q + qbar has been introduced with massive
matrix elements, as calculated and programmed by Kleiss [Kle90];
the code used in PYTHIA is a only very slightly modified copy.
- Some processes may now be used also in leptoproduction. This in
particular includes ZZ, ZW and WW interactions, like in H production.
PYTEST has been expanded to cover this possibility.
- The table produced with PYSTAT(2) also gives the flavour (KF) codes
and current masses of resonances, and the decay channel numbers in
JETSET (IDC) for the channels listed. By default, fourth generation
channels are not listed.
- A number of internal administrative changes, which should not be
visible to the outside user.
______________________________________________________________________
1.5. Major Changes from PYTHIA 5.4 to 5.5
The updates from version 5.4 to 5.5 are all minor, and just about any
program that ran with version 5.4 will also work with PYTHIA 5.5.
The following changes might give compatibility problems from programs
based on PYTHIA 5.4:
- The routine PYTHIA has now been removed; all calls to this routine
should henceforward be replaced by calls to PYEVNT.
- The expansion of process 11 to also include electroweak interactions
has been split off and is now found as process 10.
- The size of the XPQ array in PYSTFU and PYSTFE has been expanded
from -6:6 to -25:25.
- Some COEF coefficients have been moved so as to make space for the
e+e- motivated new coefficients.
- Some MINT and VINT variables have been moved to leave space for new
variables.
In addition, the following changes have been made, where no
compatibility problems are involved:
- The use of symmetrization factors for identical final state
particles has been systematized; see section 2.10.
- It is possible to have a photon as incoming particle, with the
quark and gluon structure functions of Drees and Grassie [Dre85].
- Initial state gamma radiation in e+e- annihilation has been
included in the effective structure function approach.
In an approximate fashion, one can also this way include effects
of beamstrahlung. See MSTP(11).
- Parton (quark/gluon) structure functions for the electon have
been implemented by numerical convolution of the photon structure
function inside an electron with the parton structure functions
inside a photon. Therefore a number of parton-parton scattering
processes, otherwise available only for hadron-hadron collisions,
can also be simulated in e+e- and ep events. See MSTP(12).
- Decays Z' -> Z + H and W' -> W + H have been included in the framework
of left-right-symmetric models (see [Coc90]). The strength of these
couplings can be set with PARU(145) and PARU(146), respectively.
- Several processes for leptoquark production have been implemented,
ISUB numbers 145, 162, 163, and 164. The LQ -> q + l coupling
strength can be modifed in PARU(151).
- Process 103, gamma + gamma -> H0, has been added.
- A number of photon-induced subprocesses have been included:
f gamma -> f g (33), f gamma -> f gamma (34), g gamma -> q qbar (54),
gamma gamma -> f fbar (58), g gamma -> Q Qbar (massive) (84), and
gamma gamma -> F Fbar (massive) (85).
- Process 10 has been introduced to include neutral and/or charged
current interactions in the t-channel between partons, e.g. as in
leptoproduction or Bhabha scattering. For options, see MSTP(21).
- Processes 123 and 124, f f' -> f f' H and f f' -> f" f"' H, have
beend added. These provide the exact matrix elements for what is
simulated in the effective W approximation in processes 5 and 8.
- Processes 35, f gamma -> f Z, 36, f gamma -> f' W, 69,
gamma gamma -> W+ W-, and 70, gamma W -> Z0 W, have been implemented.
However, only in the first are Z decay angles correctly included,
while the other decays currently are isotropic.
- Processes 121, g g -> Q Q~ H, and 122, q q~ -> Q Q~ H, have been
implemented, using the formulae of [Kun84].
- Simulation of the two further neutral Higgs particles, H'0 and A0,
of the two Higgs doublet scenario. This includes a total of 18 new
processes, 151 - 153, 156 - 158, 171 - 174, 176 - 179, 181 - 182,
and 186 - 187, obtained by analogy with the corresponding processes
for the standard model Higgs. Couplings can either be set freely,
see PARU(161) - PARU(185), or in accordance with the minimal
supersymmetric scenario [GHK90], see MSTP(4). This work is still
in progress: Born level couplings should work fine, while loop
contributions still are incomplete.
- Possibility of anomalous W magnetic moment in process 20.
- Possibility to simulate contact interaction (quark and lepton
substructure) in processes 11, 12, 165 and 166.
- Production of excited d* and u* quarks, by q + g -> q*, in
processes 147 and 148.
- Possibility to have incoming hyperon beams: Lambda, Sigma (-,0,+),
Xi (-,0), and Omega-. Since no structure functions are available,
an average valence quark distribution is constructed from the
proton structure function set used.
______________________________________________________________________
1.6. Major Changes from PYTHIA 5.5 to 5.6
The updates from version 5.5 to 5.6 are all minor, and just about any
program that ran with version 5.5 will also work with PYTHIA 5.6.
The following changes might give compatibility problems from programs
based on PYTHIA 5.5:
- Structure functions have been reorganized, such that the user
interface PYSTFE has been replaced by in-code interfaces to the
PDFLIB [Plo91], PAKPDF [Cha91] and PHOPDF [Cha91a] structure
function packages. To enable the appropriate interfaces, look for
lines with C! in the first two columns of the FORTRAN source file.
Also, the MSTP(51) - MSTP(58) switches for structure functions have
been changed.
In addition, the following changes have been made, where no
compatibility problems are involved:
- Addition of a routine PYEVWT to allow users to modify cross-sections
and weight events; see section 2.1.
- QCD shower evolution for photon inside electron. Intrinsic k_T
in photon; see MSTP(93).
- Additional scale choices are possible for Deep Inelastic Scattering
processes; see MSTP(22) and MSTP(32).
- It is possible to have x and Q2 of the scattered electron agree with
the originally generated x and Q2 for Deep Inelastic Scattering
processes, also when including parton showers etc.; see MSTP(23).
______________________________________________________________________
1.7. Installation of program
The program is written entirely in Fortran 77, and should run on any
machine with such a compiler. The only external program required is
JETSET 7.3. Comments on compiler optimization level, random number
generators and machine precision problems are the same as given in the
JETSET 7.3 manual. (Those wishing to use the IBM AUTODBL complier option
should note that a dummy variable MSEDUM should be inserted after MSEL
in commonblock PYSUBS.)
SAVE statements have been included in accordance with the Fortran
standard. Since most ordinary machines take SAVE for granted, this
part is not particularly well tried out, however. Users on machines
without automatic SAVE are therefore warned to be on the lookout for
any variables which may have been missed.
A test program, PYTEST, is included in the PYTHIA package. It is
disguised as a subroutine, so that the user has to run a main program
CALL PYTEST(0)
END
This program will generate some events of different types. If PYTHIA
has not been properly installed, this program is likely to crash, or
at least generate erroneous events. If everything works properly, as
far as the program can judge, a final message to that effect is
produced, else various error messages may appear. If PYTEST(1) is
called instead of PYTEST(0), the same program is run through, but with
more complete initialization and cross-section information, and with
a listing of one event of each type generated. Finally, PYTEST(2) will
in addition give an extensive listing of the initial search for
cross-section coefficients and maxima; normally there is no reason to
use this latter option.
______________________________________________________________________
1.8. Programming Philosophy
The Monte Carlo program is built as a slave system, i.e. the user
supplies the main program, and from this the various subroutines are
called on to execute specific tasks, after which control is returned to
the main program. Some of these tasks may be very trivial, whereas the
"high-level" routines by themselves may make a large number of
subroutine calls.
It should be noted that, while the physics content is obviously at
the center of attention, the Lund Monte Carlo also contains a more
extensive setup of auxiliary service routines than any other physics
event generator. The hope is that this will provide a comfortable
working environment, where not only events are generated, but users
also linger on to perform a lot of the subsequent studies. (As for
the relatively small attention given to physics in this manual,
the reason is that the physics is documented separately in a series
of papers, but the program pieces only here.)
The general rule is that all routines have names six characters long,
beginning with PY. Also commonblocks have names starting with PY.
Before events can be generated, an initialization call is necessary,
unlike the case of JETSET. Default values and other data are stored
in BLOCK DATA PYDATA. Thus this subprogram, as well as BLOCK DATA LUDATA
in JETSET, must be linked, which does not occur automatically with all
loaders.
Apart from writing initialization information, printing error messages
if need be, and responding to explicit requests for listings, all tasks
of the program are performed silently. All output is directed to unit
MSTU(11), by default 6, and it is up to the user to see to it that this
unit is open for write.
The Lund Monte Carlo is extremely versatile, but the price to be paid
for this is a large number of adjustable parameters and switches for
alternative modes of operation. No single user is ever likely to have
need for more than a fraction of the options available.
Since all these parameters and switches are assigned sensible default
values, there is no reason to worry about them until the need arises.
A number of error checks are performed during execution. Serious errors,
in particular those that may be found already at initialization time,
lead to the program being aborted. Less serious errors may lead to the
treatment of a particular event being cut short; see the JETSET manual.
Consequences are unpredictable when using integer valued input variables
with values not defined, or real-valued variables outside the physically
sensible range. Users beware!
**********************************************************************
2. The Program Components
It is useful to distinguish three phases in a normal run with PYTHIA.
In the first phase, the initialization, the general character of the
run is determined. At a minimum, this requires the specification of the
incoming hadrons and the energies involved. At the discretion of the
user, it is also possible to select specific final states, and to make
a number of decisions about details in the subsequent generation.
This step is finished by a PYINIT call, at which time several variables
are initialized in accordance with the values set. The second phase
consists of the main loop over the number of events, with each new
event being generated by a PYEVNT call. This event may then be analyzed,
using information stored in some commonblocks, and the statistics
accumulated. In the final phase, results are presented. This may often
be done without the invocation of any PYTHIA routines. From PYSTAT,
however, it is possible to obtain a useful list of cross-sections for
the different subprocesses.
In the presentation in this section, the ordering above is not
strictly adhered to. In particular, information less important for
an efficient use of PYTHIA has been delegated closer to the end.
- In subsection 2.1 the subroutines to be called by the user are
introduced, particularly PYINIT, PYEVNT and PYSTAT.
- The ISUB classification code for subprocesses included in PYTHIA
is tabulated in 2.2. Some comments on the physics of these
subprocesses, and what special possibilities are open, are given
in 2.3.
- The following two subsections, 2.4 and 2.5, deal with variables
that can (with a few exceptions) be set only in the initialization
phase, i.e. before the PYINIT call, if the default values should not
be desirable. In 2.4 is collected the switches for the type of
process to generate and kinematical constraints. Subsection 2.5
covers all the switches and parameters which regulate the details
of the generation, such as choice of structure functions or Q^2
scale, on/off switches for parton showers or fragmentation, etc.
- The information available on each new event is covered in 2.6 and
2.7, in the former the general type of process and some kinematical
variables, in the latter the detailed list of all particles
generated.
- In 2.8 the further subroutines, functions and commonblocks in
PYTHIA are listed, with brief information about their purpose.
- The fragmentation routines of JETSET 7.3 are amply documented
elsewhere [Sjo90], but in 2.9 a brief reminder is given, with
special emphasis on aspects of relevance when running PYTHIA.
- A few comments on cross-section calculations in the program are
collected in 2.10.
- Finally, subsection 2.11 contains some examples on how to run the
program.
_______________________________________________________________________
2.1. The Main Subroutines
There are two routines that users must know: PYINIT for initialization
and PYEVNT for the subsequent generation of each new event. In addition,
the cross-section and other kinds of information available with PYSTAT
is frequently useful. The three other routines described here, PYFRAM,
PYKCUT, and PYEVWT, are of more specialized interest.
SUBROUTINE PYINIT(FRAME,BEAM,TARGET,WIN)
Purpose: to initialize the generation procedure.
FRAME : a character variable used to specify the frame of the
experiment. Uppercase and lowercase letters may be freely mixed.
= 'CMS' : colliding beam experiment in CM frame, with beam momentum
in +z direction and target momentum in -z direction.
= 'FIXT' : fixed target experiment, with beam particle momentum
pointing in +z direction.
= 'USER' : full freedom to specify frame by giving beam momentum in
P(1,1), P(1,2) and P(1,3) and target momentum in P(2,1), P(2,2)
and P(2,3) in common block LUJETS.
= 'NONE' : there will be no initialization of any processes, but
only of resonance widths and a few other process-independent
variables. Subsequent to such a call, PYEVNT can not be used
to generate events, so this option is mainly intended for those
who will want to construct their own events afterwards, but
still want to have access to some of the PYTHIA facilities.
In this option, the BEAM, TARGET and WIN arguments are dummy.
BEAM, TARGET : character variables to specify beam and target particles.
Uppercase and lowercase letters may be freely mixed. An
antiparticle may be denoted either by '~' or 'bar' at the end of
the name. It is also possible to leave out the underscore ('_')
directly after 'nu' in neutrino names, and the charge for
proton and neutron.
= 'p+' : proton.
= 'p~-' : antiproton.
= 'n0' : neutron.
= 'n~0' : antineutron.
= 'pi+' : positive pion.
= 'pi-' : negative pion.
= 'e-' : electron.
= 'e+' : positron.
= 'nu_e' : electron neutrino.
= 'nu_e~' : electron antineutrino.
= 'mu-' : muon.
= 'mu+' : antimuon.
= 'nu_mu' : muon neutrino.
= 'nu_mu~' : muon antineutrino.
= 'gamma' : photon (real, i.e. on mass-shell).
= 'Lambda0' : Lambda baryon.
= 'Sigma-' : Sigma- baryon.
= 'Sigma0' : Sigma0 baryon.
= 'Sigma+' : Sigma+ baryon.
= 'Xi-' : Xi- baryon.
= 'Xi0' : Xi0 baryon.
= 'Omega-' : Omega- baryon.
WIN : related to energy of system, exact meaning depends on FRAME.
FRAME='CMS' : total energy of system (in GeV).
FRAME='FIXT' : momentum of beam particle (in GeV/c).
FRAME='USER' : dummy (information is taken from the P vectors,
see above).
SUBROUTINE PYEVNT
Purpose: to generate one event of the type specified by the PYINIT call.
(This is the main routine, which calls a number of other routines
for specific tasks.)
Note: Previously this routine was called PYTHIA; it has now been renamed
PYEVNT to avoid confusion between a subroutine and the program-as-a-
whole.
SUBROUTINE PYSTAT(MSTAT)
Purpose: to print out cross-sections statistics, decay widths,
branching ratios, status codes and parameter values. PYSTAT may be
called at any time, after the PYINIT call, e.g. at the end of the
run, or not at all.
MSTAT : specification of information desired.
= 1 : prints a table of how many events of the different kinds that
have been generated and the corresponding cross-sections. All
numbers already include the effects of cuts required by the user
in PYKCUT.
= 2 : prints a table of the resonances defined in the program, with
their particle codes (KF), and all allowed decay channels.
(If the number of generations in MSTP(1) is 3, however, channels
involving fourth generation particles are not displayed.)
For each decay channel is shown the sequential channel number
(IDC) of the JETSET decay tables, the partial decay width,
branching ratio and effective branching ratio (in the event
some channels have been excluded by the user).
= 3 : prints a table with the allowed hard interaction flavours
KFIN(I,J) for beam and target particles.
= 4 : prints a table of the kinematical cuts CKIN(I) set by the user
in the current run.
= 5 : prints a table with all the values of the status codes
MSTP(I) and the parameters PARP(I) used in the current run.
SUBROUTINE PYFRAM(IFRAME)
Purpose: to transform event between different frames, if so desired.
IFRAME : specification of frame the event is to be boosted to.
= 1 : frame specified by user in the PYINIT call.
= 2 : CM frame of incoming particles.
SUBROUTINE PYKCUT(MCUT)
Purpose: to enable a user to reject a given set of kinematic variables
at an early stage of the generation procedure (before evaluation
cross-sections), so as not to spend unnecessary time on the
generation of events that are not wanted.
The routine will not be called unless the user requires is by
setting MSTP(141) = 1, and never if "minimum bias" type events
(class d of section 2.2) are to be generated as well.
A dummy routine PYKCUT is included in the program file, so as to
avoid unresolved external references when the routine is not used.
MCUT : flag to signal effect of user-defined cuts.
= 0 : event is to be retained and generated in full.
= 1 : event is to be rejected and a new one generated.
Remark : at the time of selection, several variables in the MINT and
VINT arrays in the PYINT1 commonblock (see section 2.8) contain
information that can be used to make the decision. The routine
provided in the program file explicitly reads the variables that
have been defined at the time PYKCUT is called, and also calculates
some derived quantities. The list of information given includes
subprocess type ISUB, E_CM, s-hat, t-hat, u-hat, p_T-hat,
x_1, x_2, x_F, tau, y*, tau', cos(theta-hat), and a few more.
Some of these may not be relevant for the process under study,
and are then set to zero.
SUBROUTINE PYEVWT(WTXS)
Purpose: to allow a user to reweight event cross-sections, by process
type and kinematics of the hard scattering. There exists two
separate modes of usage, described in the following.
For MSTP(142)=1, it is assumed that the cross-section of the
process is correctly given by default in PYTHIA, but that one
wishes to generate events biased to a specific region of phase
space. While the WTXS factor therefore multiplies the naive
cross-section in the choice of subprocess type and kinematics,
the produced event comes with a compensating weight
PARI(10)=1./WTXS, which should be used when filling histograms
etc. In the PYSTAT(1) table, the cross-sections are unchanged
(up to statistical errors) compared to the standard cross-sections,
but the relative composition of events may be changed and need
no longer be in proportion to relative cross-sections. A typical
example of this usage is if one wishes to enhance the production
of high-p_T events; then a weight like WTXS=(p_T/p_T_0)**2
(with p_T_0 some fixed number) might be appropriate.
For MSTP(142)=2, on the other hand, it is assumed that the true
cross-section is really to be modifed by the multiplicative factor
WTXS. The events generated therefore come with unit weight, just
as usual. This option is really equivalent to replacing the basic
cross-sections coded in PYTHIA, but allows more flexibility: no
need to recompile the whole of PYTHIA.
The routine will not be called unless MSTP(142) >= 1, and never if
"minimum bias" type events (class d of section 2.2) are to be
generated as well. Further, cross-sections for additional multiple
interactions or pileup events are never affected.
WTXS: multiplication factor to ordinary event cross-section; to be set
(by user) in PYEVWT call.
Remark : at the time of selection, several variables in the MINT and
VINT arrays in the PYINT1 commonblock (see section 2.8) contain
information that can be used to make the decision. The routine
provided in the program file explicitly reads the variables that
have been defined at the time PYEVWT is called, and also calculates
some derived quantities. The list of information given includes
subprocess type ISUB, E_CM, s-hat, t-hat, u-hat, p_T-hat,
x_1, x_2, x_F, tau, y*, tau', cos(theta-hat), and a few more.
Some of these may not be relevant for the process under study,
and are then set to zero.
WARNING: the weights only apply to the hard scattering subprocesses.
There is no way to reweight the shape of initial and final state
showers, fragmentation, or other aspects of the event.
______________________________________________________________________
2.2. The Physics Processes
A wide selection of fundamental 2 -> 1 and 2 -> 2 tree processes of the
standard model (electroweak and strong) has been included in PYTHIA,
and slots are provided for those not yet implemented. In addition,
a few "minimum bias" type processes (like elastic scattering), loop
graphs, box graphs, 2 -> 3 tree graphs and some non-standard model
processes are included. The classification is not always unique. A
process that proceeds only via an s-channel state is classified as
a 2 -> 1 process (e.g. q q~ -> gamma* -> e+ e-), but a 2 -> 2
cross-section may well have contributions from s-channel diagrams
(g g -> g g obtains contributions from g g -> g* -> g g). Also, in
the program, 2 -> 1 and 2 -> 2 graphs may sometimes be folded with
two 1 -> 2 splittings to form effective 2 -> 3 or 2 -> 4 processes
(W+ W- -> H0 is folded with q -> q' W to give q q' -> q" q'" H0).
It is possible to select a combination of subprocesses to simulate,
and also afterwards to know which subprocess was actually selected in
each event. For this purpose, all subprocesses are numbered according
to an ISUB code. The list of possible codes is given in this section,
while its detailed use will be made clear in sections 2.4 and 2.6.
Only processes marked with a + sign in the first column have been
implemented in the program to date. The others are still given here,
so that the user will easier understand how the classification works.
In the following f_i represents a fundamental fermion of flavour i,
i.e. either of d, u, s, c, b, t, l, h, e-, nu_e, mu-, nu_mu, tau-,
nu_tau, chi- or nu_chi. A corresponding antifermion is denoted f_i~.
In several cases, some classes of fermions are explicitly excluded,
since they do not couple to the g or gamma (no e+ e- -> g g, e.g.).
When processes have only been included for quarks, while leptons
might also have been possible, the notation q_i is used. In processes
where fermion masses are explicitly included in the matrix elements,
an F is used to denote an arbitrary fermion and a Q a quark.
Flavours appearing already in the initial state are denoted i and j,
whereas new flavours in the final state are denoted k and l.
Charge conjugate channels are always assumed included as well (where
separate), and processes involving a W+ also imply those involving
a W-. Wherever Z0 is written, it is understood that gamma* and
gamma*/Z0 interference should be included as well (with possibilities
to switch off either, if so desired). In practice, the full gamma*/Z0
structure is only included in subprocesses 1 and 22; for the other
processes currently a Z0 does not contain the gamma* piece.
Correspondingly, Z'0 denotes the complete set gamma*/Z0/Z'0 (or some
subset of it). Thus the notation gamma is only used for an
on-mass-shell photon.
In the last column below, references are given to works from which
formulae have been taken. Sometimes these references are to the
original work on the subject, sometimes only to the place where the
formulae are given in the most convenient or accessible form. In
several instances, errata have been obtained from the authors. Often
the formulae given in the literature have been generalized to include
trivial radiative corrections etc.
Subprocess References
a) 2 -> 1, tree
+ 1 f_i f_i~ -> Z0 EHL84
+ 2 f_i f_j~ -> W+ EHL84
+ 3 f_i f_i~ -> H0 EHL84
4 gamma W+ -> W+
+ 5 Z0 Z0 -> H0 EHL84,Cha85
6 Z0 W+ -> W+
7 W+ W- -> Z0
+ 8 W+ W- -> H0 EHL84,Cha85
b) 2 -> 2, tree
+ 10 f_i f_j(~) -> f_i f_j(~) (QFD) Ing87
+ 11 f_i f_j(~) -> f_i f_j(~) (QCD) Cut78,Ben84,Lan91
+ 12 f_i f_i~ -> f_k f_k~ Cut78,Ben84,Lan91
+ 13 f_i f_i~ -> g g Cut78,Ben84
+ 14 f_i f_i~ -> g gamma Hal78,Ben84
+ 15 f_i f_i~ -> g Z0 EHL84
+ 16 f_i f_j~ -> g W+ EHL84
17 f_i f_i~ -> g H0
+ 18 f_i f_i~ -> gamma gamma Ber84
+ 19 f_i f_i~ -> gamma Z0 EHL84
+ 20 f_i f_j~ -> gamma W+ EHL84,Sam90
21 f_i f_i~ -> gamma H0
+ 22 f_i f_i~ -> Z0 Z0 EHL84,Gun86
+ 23 f_i f_j~ -> Z0 W+ EHL84,Gun86
+ 24 f_i f_i~ -> Z0 H0 EHL84
+ 25 f_i f_i~ -> W+ W- EHL84,Gun86
+ 26 f_i f_j~ -> W+ H0 EHL84
27 f_i f_i~ -> H0 H0
+ 28 f_i g -> f_i g Cut78,Ben84
+ 29 f_i g -> f_i gamma Hal78,Ben84
+ 30 f_i g -> f_i Z0 EHL84
+ 31 f_i g -> f_k W+ EHL84
32 f_i g -> f_i H0
+ 33 f_i gamma -> f_i g Duk82
+ 34 f_i gamma -> f_i gamma Duk82
+ 35 f_i gamma -> f_i Z0 Gab86
+ 36 f_i gamma -> f_k W+ Gab86
37 f_i gamma -> f_i H0
38 f_i Z0 -> f_i g
39 f_i Z0 -> f_i gamma
40 f_i Z0 -> f_i Z0
41 f_i Z0 -> f_k W+
42 f_i Z0 -> f_i H0
43 f_i W+ -> f_k g
44 f_i W+ -> f_k gamma
45 f_i W+ -> f_k Z0
46 f_i W+ -> f_k W+
47 f_i W+ -> f_k H0
48 f_i H0 -> f_i g
49 f_i H0 -> f_i gamma
50 f_i H0 -> f_i Z0
51 f_i H0 -> f_k W+
52 f_i H0 -> f_i H0
+ 53 g g -> f_k f_k~ Cut78,Ben84
+ 54 g gamma -> f_k f_k~ Duk82
55 g Z0 -> f_k f_k~
56 g W+ -> f_k f_l~
57 g H0 -> f_k f_k~
+ 58 gamma gamma -> f_k f_k~ Bar90
59 gamma Z0 -> f_k f_k~
60 gamma W+ -> f_k f_l~
61 gamma H0 -> f_k f_k~
62 Z0 Z0 -> f_k f_k~
63 Z0 W+ -> f_k f_l~
64 Z0 H0 -> f_k f_k~
65 W+ W- -> f_k f_k~
66 W+ H0 -> f_k f_l~
67 H0 H0 -> f_k f_k~
+ 68 g g -> g g Cut78,Ben84
+ 69 gamma gamma -> W+ W- Kat83
+ 70 gamma W+ -> Z0 W+ Kun87
+ 71 Z0 Z0 -> Z0 Z0 (longitudinal)
+ 72 Z0 Z0 -> W+ W- "
+ 73 Z0 W+ -> Z0 W+ " Ter90
74 Z0 H0 -> Z0 H0
75 W+ W- -> gamma gamma
+ 76 W+ W- -> Z0 Z0 " BLS87
+ 77 W+ W+(-) -> W+ W+(-) " Dun86
78 W+ H0 -> W+ H0
79 H0 H0 -> H0 H0
c) 2 -> 2, tree, massive final quarks
+ 81 f_i f_i~ -> Q_i Q_i~ Com77
+ 82 g g -> Q_i Q_i~ Com77
+ 83 q_i f_j -> Q_k f_l Zer90
+ 84 g gamma -> Q_i Q_i~ Fon81
+ 85 gamma gamma -> F_i F_i~ Bar90
d) "minimum bias"
+ 91 elastic scattering Blo85
+ 92 single diffraction Gou83,Blo85
+ 93 double diffraction Gou83,Blo85
94 central diffraction
+ 95 low-pT production Sjo87
e) 2 -> 1, loop
101 g g -> Z0
+ 102 g g -> H0 EHL84
+ 103 gamma gamma -> H0 Dre89
f) 2 -> 2, box
+ 111 f_i f_i~ -> g H0 Ell88
+ 112 f_i g -> f_i H0 Ell88
+ 113 g g -> g H0 Ell88
+ 114 g g -> gamma gamma Con71,Ber84,Dic88
+ 115 g g -> g gamma Con71,Ber84,Dic88
116 g g -> gamma Z0
117 g g -> Z0 Z0
118 g g -> W+ W-
119 gamma gamma -> g g
g) 2 -> 3, tree
+ 121 g g -> Q_k Q_k~ H0 Kun84
+ 122 q + q~ -> Q_k Q_k~ H0 Kun84
+ 123 f_i f_j -> f_i f_j H0 Cah84
+ 124 f_i f_j -> f_k f_l H0 Cah84
+ 131 g g -> Z0 Q_k Q_k~ Kle90
h) non-standard model, 2 -> 1
+ 141 f_i f_i~ -> Z'0 Alt89
+ 142 f_i f_j~ -> W'+ Alt89
+ 143 f_i f_j~ -> H+
+ 144 f_i f_j~ -> R Ben85a
+ 145 q_i l_j -> LQ Wud86
+ 147 d g -> d* Bau90
+ 148 u g -> u* Bau90
+ 151 f_i f_i~ -> H'0 EHL84
+ 152 g g -> H'0 EHL84
+ 153 gamma gamma -> H'0 Dre89
+ 156 f_i f_i~ -> A0 EHL84
+ 157 g g -> A0 EHL84
+ 158 gamma gamma -> A0 Dre89
i) non-standard model, 2 -> 2 and 2 -> 3
+ 161 f_i g -> f_k H+
+ 162 q g -> l LQ Hew88
+ 163 g g -> LQ LQ~ Hew88,EHL84
+ 164 q q~ -> LQ LQ~ Hew88
+ 165 f_i f_i~ -> f_k f_k~ (gamma/Z) Lan91
+ 166 f_i f_j~ -> f_k f_l~ (W) Lan91
+ 171 f_i f_i~ -> Z0 H'0 EHL84
+ 172 f_i f_j~ -> W+ H'0 EHL84
+ 173 f_i f_j -> f_i f_j H'0 Cah84
+ 174 f_i f_j -> f_k f_l H'0 Cah84
+ 176 f_i f_i~ -> Z0 A0 EHL84
+ 177 f_i f_j~ -> W+ A0 EHL84
+ 178 f_i f_j -> f_i f_j A0 Cah84
+ 179 f_i f_j -> f_k f_l A0 Cah84
+ 181 g g -> Q_k Q_k~ H'0 Kun84
+ 182 q q~ -> Q_k Q_k~ H'0 Kun84
+ 186 g g -> Q_k Q_k~ A0 Kun84
+ 187 q q~ -> Q_k Q_k~ A0 Kun84
For many of the subprocesses, further notes and qualifications may
be of interest. These are given in the following section.
Also note that some groups of subprocesses are available with the
MSEL switch, see section 2.4.
----------------------------------------------------------------------
2.3. Comments on Physics Processes
This section contains some useful comments on the processes included
in the program, grouped by physics interest rather than sequentially
by ISUB or MSEL code. The different ISUB and MSEL codes that can be
used to simulate the different groups are given. ISUB codes within
brackets indicate the kind of processes that indirectly involve the
given physics topic, although only as part of a larger whole. Some
obvious examples, like the possibility to produce jets in just about
any process, are not spelled out in detail.
The text at times contains information on which special switches or
parameters are of particular interest to a given process. All these
switches are described in detail in the following sections, but are
alluded to here so as to provide a more complete picture of the
possibilities available for the different subprocesses. However, the
list of possibilities is certainly not exhausted by the text below.
____________________
2.3.1 QCD Jets
MSEL = 1,2
ISUB = 11,12,13,28,53,68
The basic cross-sections are taken from [Cut78]. However, a string-based
fragmentation scheme such as the Lund model needs cross-sections for the
different colour flows; these have been calculated in [Ben84] and differ
from the usual calculations by interference terms of the order 1/N^2. By
default, these interference terms are excluded; however, they can be
introduced by changing MSTP(34). In this case, the interference terms
are distributed on the various colour flows according to the pole
structure of the terms.
As an example, consider subprocess 28, q + g -> q + g. The total
cross-section for this process, obtained by summing and squaring the
Feynman s-, t-, and u-channel graphs, is [Cut78]:
2*(1 - u*s/t^2) - 4/9*(s/u + u/s) - 1.
(The hats on s, t and u have been suppressed, and an overall factor
alpha_S^2/s^2 ignored.) Using the identity of the Mandelstam variables
for the massless case:
s + t + u = 0,
this can be rewritten as
(s^2 + u^2)/t^2 - 4/9*(s/u - u/s).
On the other hand, the cross-sections for the two possible colour flows
of this subprocess are [Ben84]:
A: 4/9*(2*u^2/t^2 - u/s)
B: 4/9*(2*s^2/t^2 - s/u),
the sum of which is:
8/9*(s^2 + u^2)/t^2 - 4/9*(s/u - u/s).
The difference between this expression and that of [Cut78],
corresponding to the interference between the two colour flow
configurations, is then
1/9*(s^2 + u^2)/t^2,
and can be naturally divided between colour flow A and B:
A: 1/9*u^2/t^2
B: 1/9*s^2/t^2.
This procedure is followed also for the other QCD subprocesses.
All the matrix elements in this group are for massless quarks (although
final state quarks are of course put on mass-shell). As a consequence,
some kind of regularization for p_T -> 0 is required. Normally the user
is expected to set the desired pTmin value in CKIN(3).
The new flavour produced in the annihilation processes (ISUB = 12,53)
is determined by the flavours allowed for gluon splitting into
quark-antiquark; see switch MDME.
For production of massive quarks, see below.
____________________
2.3.2 Heavy Flavours
MSEL=4,5,6,7,8,35,36,37,38
ISUB = 81,82,83,84,85
The cross-sections are taken from [Com77] for ISUB = 81,82 and from
[Zer90] for ISUB = 83. The former two processes are pure QCD ones,
and normally dominate. The last process proceeds via t channel W
exchange, and is mainly of interest for the production of very heavy
flavours, where the possibility of producing just one heavy quark
is kinematically favoured over pair production.
The matrix elements in this group differ from the corresponding ones
in the group above in that they correctly take into account the quark
masses. As a consequence, the cross-sections are finite for p_T -> 0.
It is therefore not necessary to introduce any special cuts.
The flavour produced is determined by the heaviest flavour allowed for
gluon splitting into quark-antiquark; see switch MDME. When one of the
MSEL options is used, MDME is automatically set by the program. Note
that only one heavy flavour at a time is allowed; if more than one is
turned on, only the heaviest will be produced (as opposed to the case
for ISUB = 12,53 above, where more than one flavour is allowed
simultaneously). Process 85 is set by the 'decay' channels of the
gamma rather than the gluon, and can therefore be set separately.
Here also lepton production is possible and, since lepton KF codes
come after quark ones, they are counted as being 'heavier', and thus
take precedence if they have been allowed, since again only the
heaviest flavour is simulated.
The lowest order processes above just represent one source of heavy
flavour production. Heavy quarks can also be present in the structure
functions at the Q^2 scale of the hard interaction, leading to
processes like Q + g -> Q + g, so-called flavour excitation, or be
created by gluon splittings g -> Q + Q~ in initial or final state
shower evolution. In fact, as the CM energy is increased, these other
processes gain in importance relative to the lowest order production
graphs above. As as example, only 10% of the b production at LHC
energies come from the lowest order graphs. The figure is even smaller
for charm, while it is at or above 50% for top. At LHC/SSC energies,
the specialized treatment described in this subsection is therefore
only of interest for top (and potential fourth generation quarks) -
the higher order corrections can here be approximated by an effective
K factor, except possibly in some rare corners of phase space. For
charm and bottom, on the other hand, it is necessary to simulate the
full event sample (within the desired kinematics cuts), and then only
keep those events with b/c either from lowest order production, or
flavour excitation, or gluon splitting. Obviously this may be a
time-consuming enterprise - although the probability for a high-p_T
event to contain (at least) one charm or bottom pair is fairly large,
most of these heavy flavours are carrying a small fraction of the total
p_T flow of the jets, and therefore do not survive normal experimental
cuts.
As an aside, it is not only for the lowest order graphs that events
may be generated with a guaranteed heavy flavour content. One may also
generate the flavour excitation process by itself, in the massless
approximation, using ISUB = 28 and setting the KFIN array appropriately.
No trick exists to force the gluon splittings without introducing
undesirable biases, however.
The cross-section for a heavy flavour pair close to threshold can be
modified according to the formulae of [Fad89], see MSTP(35). This
affects the total rate and also kinematical distributions.
____________________
2.3.3 Minimum Bias
MSEL = 1,2
ISUB = 91,92,93,95
The total and elastic cross-sections are given by the parametrizations
in [Blo85]. Several different parametrizations are available; these
can be selected with MSTP(31). For the single and double diffraction
cross-sections, the ansatz in [Gou83] is used. The remaining part of
the cross-section is automatically assigned to low-p_T events.
The simulation of the that variable in elastic or diffractive scattering
is still fairly primitive. Combined with the imprecise kinematics
treatment of small scattering angles, this means that the program should
not be used for studies of the scattered (undiffracted) hadron. In
diffractive scattering, the structure of the hadronic system selected
may be regulated with MSTP(101). No high-p_T jet production in
diffractive events is included so far.
The subprocess 95, low-p_T events, is somewhat unique, in
that no meaningful physical borderline to high-p_T events can be
defined. Even if the QCD 2 -> 2 high-p_T processes are formally
switched off, some of the events generated will be classified as
belonging to this group, with a p_T spectrum of interactions to
match the "minimum bias" event sample. Only with the option
MSTP(82) = 0 will subprocess 95 yield strictly low-p_T events,
events which will then probably not be compatible with any
experimental event sample. A number of options exists for the detailed
structure of low-p_T events, see in particular MSTP(81) and MSTP(82).
Details of the model(s) used for these events may be found in
[Sjo87].
____________________
2.3.4 Prompt Photons
MSEL = 10
ISUB = 14,18,29,114,115
Processes ISUB = 14,29 give the main source of single gamma production,
with ISUB = 115 giving an additional contribution which in some
kinematics regions may become important. For gamma pair production,
the process ISUB = 18 is often overshadowed in importance by ISUB = 114.
Cross-sections for the Born term graphs 14, 18 and 29 are found e.g. in
[EHL84], while the box graphs 114 and 115 are given in [Con71,Ber84,
Dic88]. For the issue of normalization of single photon vs. diphoton,
see section 2.10.
Another source of photons is bremsstrahlung off incoming or outgoing
quarks. This has to be treated on an equal footing with QCD parton
showering. For timelike parton shower evolution, i.e. in the final
state showering and in the side branches of the initial state showering,
photon emission may be switched on with MSTJ(41) (see JETSET manual).
Photon radiation off the spacelike incoming quark legs is not yet
included, but should be of lesser importance for production at
reasonably large p_T values. Radiation off an incoming electron is
included in a leading log approximation.
WARNING: The cross-sections for the box graphs 114 and 115 become
very complicated, numerically unstable and slow when the
full quark mass dependence is included. For quark masses much
below the shat scale, the simplified massless expressions are
therefore used - a fairly accurate approximation. However, there
is another set of subtle numerical cancellations between different
terms in the massive matrix elements in the region of small-angle
scattering. The associated problems have not been sorted out yet.
There are therefore two possible solutions. One is to use the
massless formulae throughout. The program then becomes faster and
numerically stable, but does e.g. not give the characteristic dip
(due to destructive interference) at top threshold. This is the
current default procedure, with five flavours assumed, but this
number can be changed in MSTP(38). The other possibility is to
impose cuts on the scattering angle of the hard process, see
CKIN(27) and CKIN(28), since the numerically unstable regions
are when abs(cos(theta-hat)) is close to unity. It is then also
necessary to change MSTP(38) to 0.
____________________
2.3.5 Single Z/W Production
MSEL = 11,12,13,14,(21)
ISUB = 1,2,15,16,30,31,35,36,131,(141)
This group consists of 2 -> 1 processes, single resonance production,
and 2 -> 2 processes, resonance + jets. With initial state parton
showers turned on, the 2 -> 1 processes also generate resonance + jets;
in order to avoid double counting, the corresponding 2 -> 2 processes
should not be turned on simultaneously. The basic rule is to use the
2 -> 1 processes for inclusive generation of Z/W, whereas the 2 -> 2
processes should be used for study of the subsample with high
transverse momentum.
The Z0 of subprocess 1 includes the full interference structure
Z0/gamma*; via MSTP(43) the user can select to produce only gamma*,
only Z0, or the full Z0/gamma*. The same holds true for the Z'0 of
subprocess 141; via MSTP(44) any combination of gamma*, Z0 and Z'0
can be selected. Thus, subprocess 141 with MSTP(44) = 4 is essentially
equivalent to subprocess 1 with MSTP(43) = 3; however, the former
also includes the possibility of a decay into Higgses. The Z0 of
other processes than 1 and 141 are currently pure Z0 only.
For the 2 -> 1 processes, the Breit-Wigner includes an shat-dependent
width, which should provide an improved description of line shapes.
The process e+e- -> e+e- + Z0 can be simulated in two different ways.
One is to make use of e 'sea' quark structure function inside e, i.e.
have splittings e -> gamma -> e (currently in the collinear limit; and
with an e+e- system in the forward direction represented by a gamma).
This can be obtained, together with ordinary Z production, by using
subprocess 1, with MSTP(11)=1 and MSTP(12)=1. Then the contribution of
the type above is 5.0 pb for a 500 GeV e+e- collider, compared to the
correct 6.2 pb [Hag90]. Alternatively one may use process 35, with
MSTP(11)=1 and MSTP(12)=0. To catch the singularity in the forward
direction, regularized by the electron mass, it is necessary to
set CKIN(3)=CKIN(5)=0.01 - using lower values will only slow down
execution, not significantly increase cross-section. One then obtains
5.1 pb, i.e. 20% below the correct, and also generates a pT distribution
for the Z0; this is therefore to be preferred.
Process 36, f + gamma -> f' + W may have corresponding problems; except
that in e+e- the forward radiation amplitude for e + gamma -> nu W
is killed (radiation zero). It is therefore feasible to use the
default CKIN(3) and CKIN(5) values in e+e-, and one also comes closer
to the correct cross-section.
One single true 2 -> 3 process is included in this class as well;
namely g g -> Z q q~, with full massive matrix elements. The more
complicated phase space and the lengthy matrix element evaluations
makes this process extremely slow. With q = b, it may form an
important background to intermediate mass Higgs searches in the
multilepton channel. The quark flavour is stored in KFPR(131,2);
the default is 5 = b. The kinematics is set up in terms of a Z
recoiling against the q q~ system, and all ordinary kinematics cut
for a 2 -> 2 process can be used on this level, including CKIN(43)
and CKIN(44) to restrict the range of the q q~ invariant mass.
In addition, for this process alone, CKIN(51) - CKIN(54) can be used
to set the p_T range of the two quarks; as is to be expected that of
the Z is set by CKIN(3) - CKIN(4). Since the optimization procedure
is not set up to probe the full multidimensional phase space allowed
in this process, maximum violations may be quite large. It may then
be useful to make a preliminary run to find how big the violations
are in total, and then use the MSTP(121) = 1 option in the full run.
____________________
2.3.6 Neutral Higgs Production
MSEL = 16
ISUB = 3,5,8,24,26,(71,72,73,76,77),102,103,111,112,113,121,122,
123,124
This section deals with the production of a standard model Higgs
exclusively, as obtained by default for the processes above.
The issue of additional Higgses is addressed in section 2.3.9.
Many different processes can be involved, as seen from the list
above. The proper choice depends on the actual Higgs mass, and
(occasionally) on the desired search strategy.
For a Higgs which can still be handled not unreasonably well in a
narrow width approximation, i.e. with a mass below 700 GeV, say,
the production processes that are involved are ISUB = 3,102,103,123,
124, as obtained for MSEL = 16.
The subprocess t t~ -> H0 (a subset of the more general subprocess
3, but the only subset of importance for heavy Higgs production at
hadronic colliders) is by now known to overestimate the cross-section
for heavy Higgs production as compared to a more careful calculation
based on the subprocess g g -> t t~ H0 (121). This process, as well
as the less important q + q~ -> t t~ H0 (122), are now included,
but the complicated multibody phase space makes it run slower than
most other processes. One should therefore think twice before using
it. As usual, it would be doublecounting to include both 3 and 121.
The complications related to g g -> t t~ H are actually not
major, since heavy Higgs production is anyway dominated by
subprocesses 102,123 and 124. The two subprocesses 5 and 8 are
effective W approximations to the matrix elements now correctly
included in 123 and 124. For light Higgs masses, the effective W
approximation is not good (it overestimates cross-sections), and
subprocesses 5 and 8 should not be used there. For a heavier Higgs,
agreement is better, and either approach may be chosen.
The subprocesses 5 and 8, V V -> H0, which contribute to the
processes V V -> V' V' (V and V' intermediate vector bosons) show
bad high energy behaviour, which can be cured only by the inclusion
of all V V -> V' V' graphs, as is done in subprocesses 71, 72, 73,
76 and 77. In particular, subprocesses 5 and 8 give rise to a
fictitious high-mass tail of the Higgs. If this tail is thrown away,
however, the agreement between the s-channel graphs (subprocesses
5 and 8) and the full set of graphs (subprocesses 71 etc.) is very
good: for a Higgs of nominal mass 300 (800) GeV, a cut at 600 (1200)
GeV retains 95% (84%) of the total cross-section, and differs from
the exact calculation, cut at the same values, by only 2% (11%)
(numbers for SSC energies). With this prescription there is
therefore no need to use subprocesses 71 etc. rather than
subprocesses 5 and 8.
For subprocesses 71, 72, 76 and 77, an option is included (see
MSTP(46)) whereby the user can select only the s-channel Higgs
graph; this will then be essentially equivalent to running
subprocess 5 or 8 with the proper decay channels (i.e. Z0Z0 or
W+W-) set via MDME. The difference is that the Breit-Wigner in
subprocesses 5 and 8 contain an shat-dependent width, whereas the
width in subprocesses 71 - 77 is calculated at nominal Higgs mass;
also, higher order corrections to the widths are treated more
accurately in subprocesses 5 and 8. Further, processes 71 - 77 assume
on-mass-shell incoming W/Z, with associated kinematics factors,
while processes 5 and 8 have W/Z correctly spacelike. All this
leads to differences in the cross-sections by up to a factor 1.5.
For subprocesses 71 - 77, also read comments in next subsection.
A subprocess like 113, with a Higgs recoiling against a gluon jet, is
also effectively generated by initial state corrections to subprocess
102; thus, in order to avoid double counting, just as for the case
of Z0/W+ above, these subprocesses should not be switched on
simultaneously, but 3, 5, 8, and 102 be used for inclusive production
of Higgs, and 111 - 113 for the study of the Higgs subsample with high
transverse momentum.
The Higgs can also be produced in association with a Z0/W,
ISUB = 24,26. In pp collisions, these processes have a lower
cross-section than single Higgs production, but are still of
interest at the lower Higgs mass range because of a potentially
better signal to background ratio. Of course, the process
Z0 -> Z0 + H0 is particularly interesting for e+e- colliders. The
other processes possible here are 103 and, in particular, 123 and 124.
As above, we note that 123 and 124 give much more accurate
cross-sections than the equivalent processes 5 and 8 based on the
effective W approximation.
The branching ratios of the Higgs are very strongly dependent on the
actual mass.In principle, the program is set up to calculate these
correctly at initialization. However, higher order corrections may at
times be important and not fully unambiguous; see e.g. MSTP(37).
____________________
2.3.7 gamma/Z/W Pairs
MSEL = 15
ISUB = 19,20,22,23,25,69,70,71,72,73,76,77
This heading contains two different types of processes:
fermion-antifermion annihilation into gamma/Z/W pairs, and scattering
of intermediate vector bosons. Obviously other sources of gamma/Z/W
pairs also may exist, like the Higgs or new gauge bosons.
The first set contains the standard gamma/Z/W produced by f + fbar
(f + fbar') annihilation. We remind that the notation gamma means a
massless gamma, while Z in principle stands for the full gamma*/Z0
structure. In fact, currently the Z0 appearing in almost all
subprocesses are pure Z0; only the Z0 of process 22 contains the full
interference structure, modifiable with MSTP(43). Currently subprocess
22 does not contain all possible contributions, but only those where
both gamma*/Z0 are produced off the quark line, i.e. excluding e.g.
emission of a gamma* off a Z0 decay product. This approximation is
fully valid for two on-mass-shell Z0:s, but should give the dominant
contribution also in more general situations. We also remind the user
that the mass ranges of the two daughters may be set with the
CKIN(41) - CKIN(44) parameters; this is particularly convenient e.g.
to pick one Z almost on mass-shell and the other not.
The second set is of interest in that only the inclusion of the full
set of VV -> V'V' graphs will cure the bad high energy behaviour of
VV -> H0 (see 2.3.6). There is an option (see MSTP(46)) that selects
only the s-channel Higgs exchange from subprocesses 71,72,76 and 77;
with this option, these subprocesses will be essentially equivalent
to subprocesses 5 and 8 if the Higgs is allowed to decay into Z or W
(see, however, further comments in the Higgs subsection).
For subprocess 77, there is an option (see MSTP(45)) to select the
charge combination of the scattering W's: like-sign, opposite-sign
(relevant for Higgs), or both.
In the absence of a Higgs, the sector of longitudinal Z and W
scattering will become strongly interacting at energies above 1 TeV.
The models proposed by Dobado, Herrero and Terron [DHT90] to describe
this kind of physics have been included as alternative matrix
elements for subprocesses 71,72,73,76,77, selectable by MSTP(46).
(From the point of view of the general classification scheme for
subprocesses, this kind of models should appropriately be included
as separate subprocesses with numbers above 100, but the current
solution allows a more efficient reuse of existing code.)
By a proper choice of parameters, it is also here possible to
simulate the production of a techni-rho.
WARNING: For subprocess 77, the option for like-sign W scattering
presently gives a non-sensical cross-section and should not be
used; the default is set for W+W- scattering. The error does not
affect the Dobado, Herrero and Terron options of this process.
____________________
2.3.8 New Gauge Bosons
MSEL = 21,22,24
ISUB = 141,142,144
The Z' of subprocess 141 contains the full gamma*/Z0/Z'0 interference
structure; the choice between different combinations is made via
MSTP(44). The couplings of Z' to quarks and leptons can be set via
PARU(121) - PARU(128), the coupling to W via PARU(129) - PARU(130), and
to H+ via PARU(143). The coupling of the Z to H+ is set via PARU(142).
By a suitable setting of these parameters, it is possible to simulate
several different physics scenarios. The default values agree with the
'extended gauge model' of [Alt89]. Further description of the coupling
parameters are given in the JETSET 7.3 manual.
The W' of subprocess 142 so far does not contain interference with the
standard model W - in practice this should not be a major limitation.
The couplings of the W' to quarks and leptons are set via PARU(131) -
PARU(134), the coupling to Z + W via PARU(135) - PARU(136). Further
comments as for Z'; in particular, default couplings again agree with
the 'extended gauge model' of [Alt89].
The strenght of couplings Z' -> Z H and W' -> W H can be set via
PARU(145) and PARU(146), respectively. Default values are in agreement
with expectations in a left-right symmmetric model [Coc90].
The R boson (particle code 40) of subprocess 144 represents one
possible scenario for a horizontal gauge boson, i.e. a gauge boson
that couples between the generations, inducing processes like
s + dbar -> mu- + e+. Experimental limits on flavour changing neutral
currents forces such a boson to be fairly heavy. The model implemented
is the one described in [Ben85a].
____________________
2.3.9 Extended Higgs scenario
MSEL = 23
ISUB = 3,24,26,102,103,121,122,123,124,(141),143,151,152,153,156,
157,158,161,171,172,173,174,176,177,178,179
The particle content of a two Higgs doublet scenario is included:
two neutral scalar particles, 25 and 35, one pseudoscalar one, 36,
and a charged doublet, 37. (Of course, these particles may also be
associated with corresponding Higgs states in larger multiplets).
By convention we choose to call the lighter scalar Higgs H0 and
the heavier H'0 - this differs from the convention in the minimal
supersymmetric extension of the standard model (MSSM), where the
ligher is called h0 and the heavier H0, but allows us to call the
Higgs of the one-Higgs scenario H0. The pseudoscalar is called A0
and the charged H+ and H-.
A number of H0 processes have been duplicated for H'0 and A0.
The correspondence between ISUB numbers is as follows.
H0 H'0 A0
f_i f_i~ -> 3 151 156
g g -> 102 152 157
gamma gamma -> 103 153 158
f_i f_i~ -> Z + 24 171 176
f_i f_j~ -> W+ + 26 172 177
g g -> Q_k Q_k~ + 121 181 186
q_i + q_i~ -> Q_k Q_k~ + 122 182 187
f_i f_j -> f_i f_j + 123 173 178
f_i f_j -> f_k f_l + 124 174 179
Note that several of the processes above are not expected to take
place at all, due to vanishing Born term couplings. We have still
included them for flexibility in simulating arbitrary couplings at
Born or loop level.
A few standard model Higgs processes have no correspondence in the
scheme above. These include: 5 and 8, which anyway have been
superseded by 123 and 124; 71, 72, 73, 76, 77, which deal with
what happens if there is no light Higgs, and so is a scenario
complementary to the one above, where several light Higgses are
assumed; and 111, 112 and 113, which describe the high-p_T tail
of the Higgs production, and are less interesting for most Higgs
studies.
In processes 121,122,181,182,186,187 the recoiling heavy flavour
is assumed to be top, which is the only one of interest in the
standard model, and the one where the structure function approach
invoked in processes 3,151,156 is least reliable. However, it is
possible to change the quark flavour in 121 etc.; for each process
ISUB this flavour is given by KFPR(ISUB,2).
By default, the H0 has the couplings of the standard model Higgs,
while the H'0 and A0 have couplings set in PARU(171) - PARU(177)
and PARU(181) - PARU(187). The default values for the latter have
no deep physics motivation, but are set just so the program will
not crash due to the absence of any couplings whatsoever. The user
should therefore set the couplings above to his desired values if
he wants to simulate either H'0 or A0. Also the couplings of the
H0 particle can be modified, in PARU(161) - PARU(165), provided
that MSTP(4) is set to 1.
For MSTP(4) = 2, the mass of the H0 (in PMAS(25,1)) and the
tan(beta) value (in PARU(141)) are used to derive the masses of
the other Higgses, as well as all Higgs couplings. PMAS(35,1) -
PMAS(37,1) and PARU(161) - PARU(195) are overwritten accordingly.
The relations used are the ones of the Born level MSSM [GHK90].
Note that not all combinations of m_H and tan(beta) are allowed;
the requirement of a finite A mass imposes the constraint
m_H < m_Z * (tan^2(beta) - 1)/(tan^2(beta) + 1),
or, equivalently,
tan^2(beta) > (m_Z + m_H)/(m_Z - m_H).
If this condition is not fulfilled, the program will crash.
The basic subprocess for charged Higgs production is ISUB = 143.
However, this process is dominated by t + b~ -> H+, and so
depends on the choice of t structure function. A better
representation is provided by subprocess 161, f + g -> f' + H+-;
i.e. actually b + g -> t + H-. It is therefore recommended to use
161 and not 143. The tan(beta) parameter, which is relevant also
for charged Higgs couplings, is set via PARU(141).
Note in particular, that subprocess 141 is one possibility for
H+ production, in addition to subprocesses 143 and 161. Note also,
that it is only via subprocess 141 that a Z can be made to decay
into charged Higgses. The coupling of the Z to H+ + H- is regulated
by PARU(142), and that of the Z' by PARU(143).
Process 141 can also be used to simulate Z0 -> H0 A0 and
Z0 -> H'0 A0 for associated neutral Higgs production. Finally,
heavier Higgses may decay into lighter ones, if kinematically
allowed, in processes like A -> Z0 H0 or H+ -> W+ H0.
____________________
2.3.10 Leptoquarks
MSEL = 25
ISUB = 145,162,163,164
Several processes that can generate a leptoquark have been included.
Currently only one leptoquark has been implemented, as particle 39.
The leptoquark is assumed to carry specific quark and lepton quantum
numbers, by default u quark + electron. These flavour numbers are
conserved, i.e. a process like u + e -> LQ -> d + nu is not allowed.
This may be a bit restrictive, but represents one of many leptoquark
possibilities.
Although only one leptoquark is implemented, its flavours may be
changed arbitrarily to study the different possibilities. The flavours
of the leptoquark are defined by the quark and lepton flavours in the
decay mode list. Since only one decay channel is allowed, this means
that the quark flavour is stored in KFDP(MDCY(39,2),1) and the lepton
one in KFDP(MDCY(39,2),2). The former must always be a quark, while
the latter could be a lepton or an antilepton; a charge conjugate
partner is automatically defined by the program. At initialization,
the charge is recalculated as a function of the flavours defined; also
the leptoquark name is redefined to be of the type 'LQ_(q)(l)',
where actual quark (q) and lepton (l) flavours are displayed.
The LQ -> q + l vertex contains an undetermined Yukawa coupling
strength, which affects both the width of the leptoquark and the
cross-section for many of the production graphs. This strenght may
be changed in PARU(151). The definition of PARU(151) corresponds to
the k factor of [Hew88], i.e. to lambda**2/(4*pi*alpha_em), where
lambda is the Yukawa coupling strenght of [Wud86]. Note that PARU(151)
thus is quadratic in the coupling.
The leptoquark is likely to be fairly long-lived, in which case it
has time to fragment into a mesonic or baryonic type states, which
would decay later on. This is a bit tedious to handle, so therefore
the leptoquark is always assumed to decay before fragmentation has
to be considered. This may give some imperfections in the event
generation, but should not be off by much in the final analysis.
Inside the program, the leptoquark is treated as a resonance, actually
the only coloured such. This requires some extra care, and in particular
it is not allowed to put the leptoquark stable, by modifying either
MDCY(39,1) or MSTP(41): then the leptoquark would be handed undecayed to
JETSET, which would try to fragment it (as it does with any other
coloured object), and most likely crash.
____________________
2.3.11 Compositeness and Anomalous Couplings
ISUB = 11,12,20,165,166
Some processes have been set up to allow anomalous coupling to be
introduced, in addition to the standard model ones. These can be
switched on by MSTP(5)>=1; the default MSTP(5)=0 corresponds to the
standard model behaviour.
In process 11 and 12, quark substructure is included in the left-left
isoscalar model [Lan91] for MSTP(5)=1; with compositeness scale Lambda
given in PARU(155) (default 1000 GeV) and sign eta of interference
term in PARU(156) (default +1; only other alternative -1). The model
above assumes that only u and d quarks are composite (at least at the
scale studied); with MSTP(5)=2 compositeness terms are included in the
interactions between all quarks.
The processes 165 and 166 are basically equivalent to 1 and 2, i.e.
gamma/Z0 and W+- exchange, respectively, but a bit less fancy (no
s-dependent width etc.). Only one final state flavour at the time is
generated; this flavour should be set in KFPR(165,1) and KFPR(166,1),
respectively. For process 166 one gives the down-type flavour, and the
program will associate the up-type flavour of the same generation.
Defaults are 11 in both cases, i.e. e+e- and e+nu_e/e-nu_e~ final
states. While MSTP(5)=0 gives the standard model results, MSTP(5)=1
contains the left-left model isoscalar model (which does not affect
process 166), and MSTP(5)=3 the helicity-nonconserving model (which
affects both) [Lan91]. Both the models above assume that only u and d
quarks are composite; with MSTP(5)=2 or 4, respectively, contact terms
are included for all quarks in the initial state. Parameters are
PARU(155) and PARU(156), as above.
Note that processes 165 and 166 are bookkept as 2 -> 2 processes,
while 1 and 2 are 2 -> 1 ones. This means that the default Q2 scale in
structure functions is p_T^2 for the former and shat for the latter.
To make contact between the two, it is recommended to set MSTP(32)=4,
so as to use shat as scale also for processes 165 and 166.
In process 20, for W + gamma pair production, it is possible to set
an anomalous magnetic moment for the W in PARU(153) (= eta = kappa-1;
where kappa = 1 is the standard model value). The production process
is affected according to the formulae of [Sam90]; W decay currently
remains unaffected. It is necessary to set MSTP(5)=1 to enable this
extension.
____________________
2.3.12 Excited fermions
ISUB = 147,148
Compositeness scenarios may also give rise to sharp resonances of
excited quarks and leptons. If MSTP(6)=1, then at initialization
the standard fourth generation of fermions will be overwritten,
and made to correspond to an excited copy of the first generation,
consisting of spin 1/2 particles d* (code 7), u* (8), e* (17) and
nu*_e (18). Since the original fourth generation information is lost,
it is then not in the same run possible to generate fourth
generation particles.
The current implementation only contains gauge interactions, i.e.
no contact interactions. The couplings f, f' and f_s to the SU(2),
U(1) and SU(3) groups are stored in PARU(157) - PARU(159), the
scale parameter Lambda in PARU(155); users are also expected to
change f* masses in accordance with what is desired - see [Bau90]
for details on conventions. This means that decay processes are of
the types q* -> q + g, q + gamma, q + Z0 or q' + W+-. Production
is currently only by quark-gluon fusion.
______________________________________________________________________
2.4. Switches for Event Type and Kinematics Selection
By default, only QCD 2 -> 2 processes are generated by PYTHIA,
composed of hard interactions above p_T-hat_min = PARP(81),
with low-p_T processes added on so as to give the full (parametrized)
inelastic, non-diffractive cross-section. With the help of the
commonblock PYSUBS, it is possible to select for the generation of
another process, or a combination of processes. It is also allowed
to restrict the generation to specific incoming partons/particles
at the hard interaction. This often automatically also restricts
final state flavours but, in processes like resonance production
(Z0, W+, H0, Z'0, H+ or R0) or QCD/QED production of new flavours
(ISUB = 12, 53, 54, 81, 82, 83, 84, 85), switches in the JETSET
program may be used to this end; see section 2.9.
The CKIN array may be used to impose specific kinematics cuts.
The user should here be warned that, if kinematical variables are
too strongly restricted, the generation time per event may become
very long. In extreme cases, where the cuts effectively close the
full phase space, the event generation may run into an infinite
loop. The generation of 2 -> 1 resonance production is performed in
terms of the m-hat and y* variables, and so the ranges CKIN(1) -
CKIN(2) and CKIN(7) - CKIN(8) may be arbitrarily restricted without
a significant loss of speed. For 2 -> 2 processes, cos(theta-hat) is
added as a third generation variable, and so additionally the range
CKIN(27) - CKIN(28) may be restricted without any danger.
In addition to the variables found in PYSUBS, also those in the
PYPARS commonblock may be used to select exactly what one wants
to have simulated. These possibilities will be described in the
following subsection.
The notation used above and in the following is that 'hat' denotes
internal variables in the hard scattering subsystem, while '*' is
for variables in the CM frame of the event-as-a-whole. Effects from
initial and final state radiation are not included, since they are not
known at the time the kinematics at the hard interaction is selected.
The sharp kinematical cutoffs that can be imposed on the generation
process are therefore smeared, both by QCD radiation and by
fragmentation. In a study of events within a given window of
experimentally defined variables, it is up to the user to leave
such liberal margins that no events are missed. In other words, cuts
have to be chosen such that a negligible fraction of events migrate
from outside the simulated region to inside the interesting region.
Often this may lead to low efficiency in terms of what fraction of
the generated events are actually of interest to the user.
COMMON/PYSUBS/MSEL,MSUB(200),KFIN(2,-40:40),CKIN(200)
Purpose: to allow the user to run the program with any desired subset
of processes, or restrict flavours or kinematics.
MSEL (D=1) a switch to select between full user control and some
preprogrammed alternatives.
= 0 : desired subprocesses have to be switched on in MSUB, i.e.
full user control.
= 1 : depending on incoming particles, three different alternatives
are used.
Lepton-lepton: Z or W production (ISUB = 1 or 2).
Lepton-hadron: deep inelastic scattering (ISUB = 10).
Hadron-hadron: QCD high-p_T processes (ISUB = 11, 12, 13, 28,
53, 68); additionally low-p_T production if CKIN(3) <
PARP(81) or PARP(82), depending on MSTP(82) (ISUB = 95).
If low-p_T is switched on, the other CKIN cuts are not
used.
= 2 : as MSEL = 1 for lepton-lepton and lepton-hadron. For
hadron-hadron all QCD processes, including low-pT, single and
double diffractive and elastic scattering, are included (ISUB =
11, 12, 13, 28, 53, 68, 91, 92, 93, 95). The CKIN cuts are
here not used.
= 4 : charm (cc~) production with massive matrix elements
(ISUB = 81, 82, 84, 85).
= 5 : bottom (bb~) production with massive matrix elements
(ISUB = 81, 82, 84, 85).
= 6 : top (tt~) production with massive matrix elements
(ISUB = 81, 82, 84, 85).
= 7 : low (ll~) production with massive matrix elements
(ISUB = 81, 82, 84, 85).
= 8 : high (hh~) production with massive matrix elements
(ISUB = 81, 82, 84, 85).
= 10 : prompt photons (ISUB = 14, 18, 29).
= 11 : Z0 production (ISUB = 1).
= 12 : W+/- production (ISUB = 2).
= 13 : Z0 + jet production (ISUB = 15, 30).
= 14 : W+/- + jet production (ISUB = 16, 31).
= 15 : pair production of different combinations of gamma, Z0 and
W+/- (except gamma + gamma; see MSEL = 10) (ISUB = 19, 20,
22, 23, 25).
= 16 : H0 production (ISUB = 3, 102, 103, 123, 124).
= 17 : H0 + Z0 or H0 + W+/- (ISUB = 24, 26).
= 18 : H0 production, combination relevant for e+e- annihilation
(ISUB = 24, 103, 123, 124).
= 19 : H0, H'0 and A0 production, excepting pair production
(ISUB = 24, 103, 123, 124, 153, 158, 171, 173, 174, 176,
178, 179).
= 21 : Z'0 production (ISUB = 141).
= 22 : W'+/- production (ISUB = 142).
= 23 : H+/- production (ISUB = 143).
= 24 : R0 production (ISUB = 144).
= 25 : LQ (leptoquark) production (ISUB = 145, 162, 163, 164).
= 35: single bottom production by W exchange (ISUB = 83).
= 36: single top production by W exchange (ISUB = 83).
= 37: single low production by W exchange (ISUB = 83).
= 38: single high production by W exchange (ISUB = 83).
MSUB : (D=200*0) array to be set when MSEL = 0 (for MSEL >= 1 relevant
entries are set in PYINIT) to choose which subset of subprocesses to
include in the generation. The ordering follows the ISUB code given
in subsection 2.2 (with comments as given there).
MSUB(ISUB) = 0 : the subprocess is excluded.
MSUB(ISUB) = 1 : the subprocess is included.
Note: when MSEL >= 1 the MSUB values set by the user are never
changed by PYTHIA. If the user wants to combine several
different 'subruns', each with its own PYINIT call, into one
single run, it is up to him to remember not only to switch on
the new processes before each new PYINIT call, but also to
switch off the old ones that are no longer desired.
KFIN(I,J) : provides an option selectively to switch on and off
contributions to the cross-sections from the different incoming
partons/particles at the hard interaction. In combination with the
JETSET resonance decay switches, this also allows the user to set
restrictions on flavours appearing in the final state.
I : is 1 for beam side of event and 2 for target side.
J : enumerates flavours according to the KF code; see section 2.7,
or the JETSET manual.
KFIN(I,J) = 0 : the parton/particle is forbidden.
KFIN(I,J) = 1 : the parton/particle is allowed.
By default, everything is on, except for J = 0, which does not
have a physical meaning.
CKIN : kinematics cuts that can be set by the user before the PYINIT
call, and that affect the region of phase space within which
events are generated. Some cuts are 'hardwired' while most are
'softwired'. The hardwired ones are directly related to the
kinematical variables used in the event selection procedure,
and therefore have negligible effects on program efficiency.
The most important of these are CKIN(1) - CKIN(8), CKIN(27) -
CKIN(28), and CKIN(31) - CKIN(32). The softwired ones are most
of the remaining ones, that can not fully be taken into account
in the kinematical variable selection, so that generation in
constrained regions of phase space may be slow. In extreme
cases the phase space may be so small that the maximization
procedure fails to find any allowed points at all (though some
small region might still exist somewhere), and therefore switches
off some subprocesses, or aborts altogether.
CKIN(1), CKIN(2) : (D=2.,-1.) range of allowed m-hat = sqrt(s-hat)
values. IF CKIN(2) < 0., the upper limit is inactive.
CKIN(3), CKIN(4) : (D=0.,-1.) range of allowed p_T-hat values for
hard 2 -> 2 processes, with transverse momentum p_T-hat defined
in the rest frame of the hard interaction. If CKIN(4) < 0.,
the upper limit is inactive. For processes which are singular
in the limit p_T-hat -> 0 (see CKIN(6)), CKIN(5) provides
an additional constraint. The CKIN(3) and CKIN(4) limits can
also be used in 2 -> 1 -> 2 processes. Here, however, the product
masses are not known and hence assumed vanishing in the event
selection. The actual p_T range for massive products is thus
shifted downwards with respect to the nominal one.
CKIN(5) : (D=1.) lower cutoff on p_T-hat values, in addition to the
CKIN(3) cut above, for processes which are singular in the
limit p_T-hat -> 0 (see CKIN(6)).
CKIN(6) : (D=1.) hard 2 -> 2 processes, which do not proceed only
via an intermediate resonance (i.e. are 2 -> 1 -> 2 processes),
are classified as singular in the limit p_T-hat -> 0 if either
or both of the two final state products has a mass m < CKIN(6).
CKIN(7), CKIN(8) : (D=-10.,10.) range of allowed scattering
subsystem rapidities y* in the CM frame of the event,
where y* = (1/2) * ln(x_1/x_2).
CKIN(9), CKIN(10) : (D=-10.,10.) range of allowed (true) rapidities
for the product with largest rapidity in a 2 -> 2 or a
2 -> 1 -> 2 process, defined in the CM frame of the event,
i.e. max(y*_3, y*_4).
CKIN(11), CKIN(12) : (D=-10.,10.) range of allowed (true) rapidities
for the product with smallest rapidity in a 2 -> 2 or a
2 -> 1 -> 2 process, defined in the CM frame of the event,
i.e. min(y*_3, y*_4). Consistency thus requires
CKIN(11) <= CKIN(9) and CKIN(12) <= CKIN(10).
CKIN(13), CKIN(14) : (D=-10.,10.) range of allowed pseudorapidities
for the product with largest pseudorapidity in a 2 -> 2 or a
2 -> 1 -> 2 process, defined in the CM frame of the event,
i.e. max(eta*_3, eta*_4).
CKIN(15), CKIN(16) : (D=-10.,10.) range of allowed pseudorapidities
for the product with smallest pseudorapidity in a 2 -> 2 or a
2 -> 1 -> 2 process, defined in the CM frame of the event,
i.e. min(eta*_3, eta*_4). Consistency thus requires
CKIN(15) <= CKIN(13) and CKIN(16) <= CKIN(14).
CKIN(17), CKIN(18) : (D=-1.,1.) range of allowed cos(theta*) values
for the product with largest cos(theta*) value in a 2 -> 2 or a
2 -> 1 -> 2 process, defined in the CM frame of the event,
i.e. max(cos(theta*_3),cos(theta*_4)).
CKIN(19), CKIN(20) : (D=-1.,1.) range of allowed cos(theta*) values
for the product with smallest cos(theta*) value in a 2 -> 2 or a
2 -> 1 -> 2 process, defined in the CM frame of the event,
i.e. min(cos(theta*_3),cos(theta*_4)). Consistency thus requires
CKIN(19) <= CKIN(17) and CKIN(20) <= CKIN(18).
CKIN(21), CKIN(22) : (D=0.,1.) range of allowed x_1 values for
the parton on side 1 that enters the hard interaction.
CKIN(23), CKIN(24) : (D=0.,1.) range of allowed x_2 values for
the parton on side 2 that enters the hard interaction.
CKIN(25), CKIN(26) : (D=-1.,1.) range of allowed Feynman-x values,
where x_F = x_1 - x_2.
CKIN(27), CKIN(28) : (D=-1.,1.) range of allowed cos(theta-hat) values
in a hard 2 -> 2 scattering, where theta-hat is the scattering
angle in the rest frame of the hard interaction.
CKIN(31), CKIN(32) : (D=2.,-1.) range of allowed m'-hat values, where
m'-hat is the mass of the complete three- or four-body final state
in 2 -> 3 or 2 -> 4 processes (while m-hat, constrained in CKIN(1)
and CKIN(2), here corresponds to the one- or two-body central
system). If CKIN(32) < 0., the upper limit is inactive.
CKIN(35), CKIN(36) : (D=0.,-1.) range of allowed abs(t-hat) (= -that)
values in 2 -> 2 processes. Note that for deep inelastic scattering
this is nothing but the Q2 scale, in the limit that initial and
final state radiation is neglected. If CKIN(36) < 0., the upper
limit is inactive.
CKIN(37), CKIN(38) : (D=0.,-1.) range of allowed abs(u-hat) (= -uhat)
values in 2 -> 2 processes. If CKIN(38) < 0., the upper limit is
inactive.
CKIN(41) - CKIN(44) : (D=12.,-1.,12.,-1.) range of allowed mass values
of the two (or one) resonances produced in a 'true' 2 -> 2
process, i.e. one not (only) proceeding through a single s-channel
resonance (2 -> 1 -> 2). (These are the ones listed as 2 -> 2 in the
table in section 2.2.) Only particles with a width above PARP(41)
are considered as bona fide resonances and tested against the CKIN
limits; particles with a smaller width are put on mass-shell
without applying any cuts. The exact interpretation of the CKIN
variables depends on the flavours of the produced two resonances.
For two resonances like Z + W (produced from f + f' -> Z + W),
which are not identical and which are not each other's
antiparticles, one has
CKIN(41) < m1 < CKIN(42), and
CKIN(43) < m2 < CKIN(44),
where m1 and m2 are the actually generated masses of the two
resonances, and 1 and 2 are defined by the order in which they
given in the production process specification, see section 2.2.
For two resonances like Z + Z, which are identical, or W+ + W-,
which are each other's antiparticles, one instead has
CKIN(41) < min(m1,m2) < CKIN(42), and
CKIN(43) < max(m1,m2) < CKIN(44).
In addition, whatever limits are set on CKIN(1) and, in particular,
CKIN(2) obviously affects the masses actually selected.
Note 1: If MSTP(42) = 0, so that no mass smearing is allowed,
the CKIN values have not effect (the same as for particles
with too narrow a width).
Note 2: If CKIN(42) < CKIN(41) or CKIN(44) < CKIN(43) it means
that the CKIN(42) or CKIN(44) limit is inactive.
Note 3: If limits are active and the resonances are identical, it
is up to the user to ensure that CKIN(41) <= CKIN(43) and
CKIN(42) <= CKIN(44).
Note 4: For identical resonances, it is not possible to preselect
which of the resonances is the lighter one; if e.g. one Z is
to decay to leptons and the other to quarks, there is no
mechanism to guarantee that the lepton pair has a smaller mass
than the quark one.
Note 5: The CKIN values are applied to all relevant 2 -> 2
processes equally, which may not be what one desires if several
processes are generated simultaneously. Some caution is
therefore urged in the use of the CKIN(41) - CKIN(44) values.
Also in other respects, users are recommended to take proper
care - if a Z is only allowed to decay into b + bbar, e.g.,
setting its mass range to be 2 - 8 GeV is obviously not a
good idea.
Note 6: In principle, the machinery should work for any 2 -> 2
process with resonances in the final state, but so far it has
only been checked for processes 22 - 26, so also from this
point some caution is urged.
CKIN(45) - CKIN(48) : (D=12.,-1.,12.,-1.) range of allowed mass values
of the two (or one) secondary resonances produced in 2 -> 1 -> 2
process (like g + g -> H0 -> Z0 + Z0) or even a 2 -> 2 -> 4 (or 3)
process (like q + qbar -> Z0 + H0 -> Z0 + W+ + W-). Note that these
CKIN values only affect the secondary resonances; the primary ones
are constrained by CKIN(1), CKIN(2) and CKIN(41) - CKIN(44).
(indirectly, of course, the choice of primary resonance masses
affects the allowed mass range for the secondary ones).
What is considered to be a resonance is defined by PARP(41);
particles with a width smaller than this are automatically put on
mass-shell. The description closely parallels the one given for
CKIN(41) - CKIN(44). Thus, for two resonances which are not
identical or each other's antiparticles, one has
CKIN(45) < m1 < CKIN(46), and
CKIN(47) < m2 < CKIN(48),
where m1 and m2 are the actually generated masses of the two
resonances, and 1 and 2 are defined by the order in which they
given in the decay channel specification in the program (see e.g.
output from PYSTAT(2) or LULIST(12)). For two resonances which
are identical or each other's antiparticles, one instead has
CKIN(45) < min(m1,m2) < CKIN(46), and
CKIN(47) < max(m1,m2) < CKIN(48).
Notes 1 - 5: as for CKIN(41) - CKIN(44), with trivial modifications.
Note 6: Setting limits on secondary resonance masses is possible
in any of the channels of the allowed types (see above).
However, so far only H -> Z0 + Z0 and H -> W+ + W- have been
fully implemented, such that an arbitrary mass range below
the naive mass threshold may be picked. For other possible
resonances, any restrictions made on the allowed mass range
are not reflected in the cross-section; and further it is not
recommendable to pick mass windows that makes an on-mass-shell
decay impossible. These limitations will be relaxed in future
versions.
CKIN(51) - CKIN(56) : (D=0.,-1.,0.,-1.,0.,-1.) range of allowed
transverse momenta in a true 2 -> 3 process. Currently two
different alternatives are around. For subprocess 131, the p_T
of the first product (the Z) is set by CKIN(3) and CKIN(4),
while for the quark and antiquark p_T:s one has
CKIN(51) < min(p_T_q,p_T_q~) < CKIN(52)
CKIN(53) < max(p_T_q,p_T_q~) < CKIN(54)
Negative CKIN(52) and CKIN(54) values means that the corresponding
limits are inactive.
For subprocesses 121-124, and their H'0 and A0 equivalents (173,
174,178,179,181,182,186,187), CKIN(51) - CKIN(54) again corresponds
to p_T ranges for scattered partons, but in order of appearance,
i.e. CKIN(51) - CKIN(52) for the parton scattered off the beam and
CKIN(53) - CKIN(54) for the one scattered off the target. CKIN(55)
and CKIN(56) here sets p_T limits for the third product, the H0,
i.e. the CKIN(3) and CKIN(4) values have no effect for this
process. Since the p_T of the Higgs is not one of the primary
variables selected, any constraints here may mean reduced Monte
Carlo efficiency, while for these processes CKIN(51) - CKIN(54)
are 'hardwired' and therefore do not cost anything.
______________________________________________________________________
2.5. The General Switches and Parameters
In addition to the event information described in section 2.6, the
PYPARS commonblock contains the status code and parameters which
regulate the performance of the program. All of them are provided with
sensible default values, so that a novice user can neglect them, and
only gradually explore the full range of possibilities.
COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
Purpose: to give access to status code and parameters which regulate the
performance of the program. If the default values, below denoted by
(D=...), are not satisfactory, they must in general be changed
before the PYINIT call. Exceptions, i.e. variables which can be
changed for each new event, are denoted by (C).
MSTP(1) : (D=3) maximum number of generations. Automatically set <= 4.
MSTP(2) : (D=1) calculation of alpha_strong at hard interaction, in the
routine ULALPS.
= 0 : alpha_strong is fixed at value PARU(111).
= 1 : first order running alpha_strong.
= 2 : second order running alpha_strong.
MSTP(3) : (D=2) selection of Lambda value in alpha_strong for
MSTP(2) >= 1.
= 1 : Lambda is given by PARP(1) for hard interactions, by PARP(61)
for spacelike showers and by PARJ(81) for timelike ones. This
Lambda is assumed valid for 5 flavours; for the hard interaction
the number of flavours assumed can be changed by MSTU(112).
= 2 : Lambda value is chosen according to the structure function
parametrizations, i.e. Lambda = 0.20 GeV for EHLQ1, = 0.29 GeV
for EHLQ2, = 0.20 GeV for DO1, = 0.40 GeV for DO2, = 0.187 GeV
for MT1, = 0.212 GeV for MT2, = 0.191 GeV for MT3, = 0.155 GeV
for MT4, = 0.22 GeV for GRV1, = 0.16 GeV for GRV2, etc. for
structure functions in the PDFLIB and PAKPDF libraries (cf.
(MSTP(51), MSTP(52)). The choice is always based on the proton
structure function set selected, i.e. is unaffected by pion and
photon structure function selection. All the Lambda values above
are assumed to refer to 4 flavours, and MSTU(112) is set
accordingly. This Lambda value is used both for the hard
scattering and the initial and final state radiation. The
ambiguity in the choice of Q^2 argument still remains (see
MSTP(32), MSTP(64) and MSTJ(44)). This Lambda value is used
also for MSTP(57) = 0, but the sensible choice here would be
to use MSTP(2) = 0 and have no initial or final state radiation.
MSTP(4) : (D=0) treatment of the Higgs sector, predominantly the
neutral one.
= 0 : the H is given the standard model Higgs couplings, while
H'0 and A0 couplings should be set by the user in PARU(171) -
PARU(175) and PARU(181) - PARU(185), respectively.
= 1 : the user should set couplings for all three Higgses, for
the H0 in PARU(161) - PARU(165), and for the H'0 and A0 as
above.
= 2 : the mass of H0 in PMAS(25,1) and the tan(beta) value in
PARU(141) are used to derive H'0, A0 and H+- masses, and H0,
H'0, A0 and H+- couplings, using the relations of the minimal
supersymmetric extension to the standard model at Born level
[GHK90]. Existing mass and couplings are overwritten by the
derived values. See section 2.3.9 for discussion in parameter
constraints.
= 3: as =2, but using relations at the one loop level. This option
is not yet implemented.
MSTP(5) : (D=0) presence of anomalous couplings in processes.
= 0 : absent.
>= 1 : present, wherever implemented. See section 2.3.11 for further
details.
MSTP(6) : (D=0) usage of the fourth generation fermions to simulate
other fermions.
= 0 : none, i.e. can be used as a standard forurht generation.
= 1 : excited fermions, as present in compositeness scenarios;
see section 2.3.12.
MSTP(11) : (D=0) use of electron structure function in e+e- and ep
interactions.
= 0 : no, i.e. electron carries whole beam energy.
= 1 : yes, i.e. electron only carries fraction of beam energy in
agreement with next-to-leading electron structure function
in [Kle89] (p. 34, eqs. (58), (60)), thereby including the
effects of initial state bremsstrahlung.
= 2 : include the effects of beamstrahlung in terms of a simple,
approximate effective electron structure function, of the
type D(x) = N * x^a (1-x)^b, where a = PARP(11), b = PARP(12),
and N a normalization factor so as to make the integral of
D(x) equal to one. The mean x value is = (1+a)/(2+a+b).
= 3: include, approximately, the effects of both beamstrahlung and
bremsstrahlung, in terms of a parametrization as in = 2, but
with b = PARP(12) + 0.5 * beta, where
beta = (2.*alpha_em/pi) * (ln(Q^2/m_e^2) - 1)
is the traditional factor associated with bremsstrahlung.
MSTP(12) : (D=0) use of e- ('sea', i.e. from e -> gamma -> e), e+,
quark and gluon structure functions inside an electron.
= 0 : off.
= 1 : on, provided that MSTP(11) >= 1. The e- and e+ structure
functions are given by [Che75], divided by two (reduction
due to Q2 ordering). The other structure functions
are obtained by numerical convolution of the photon content
inside an electron (as given by the bremsstrahlung spectrum of
MSTP(11) = 1) with the quark and gluon content inside a photon
(as given by the Drees and Grassie structure functions).
Required numerical precision is set by PARP(14).
Since the need for numerical integration makes this option
somewhat more time-consuming than ordinary structure function
evaluation, one should only use it when studying processes
where it is needed.
MSTP(13) : (D=1) choice of Q^2 range over which electrons are assumed
to radiate photons; affects normalization of e- (sea), e+, photon,
quark and gluon structure functions inside electron.
= 1 : range set by Q^2 argument of structure function call, i.e.
by Q^2 scale of hard interaction. Therefore structure functions
are proportional to ln(Q^2/m_e^2). This is normally most
appropriate for e+e- annihilation.
= 2 : range set by user determined Q_max^2, given in PARP(13).
Structure functions are assumed proportional to
ln((Q_max^2/m_e^2)*(1-x)/x^2). This is normally most appropriate
for photoproduction, where the electron is supposed to go
undetected, i.e. scatter less than Q_max^2.
MSTP(14) : (D=1) structure of incoming photon beam or target (does not
affect photon inside electron, only photons appearing as argument
in the PYINIT call).
= 0 : photon is assumed pointlike, i.e. can only interact in
processes which contain an incoming photon, such as ISUB =
33 - 36.
= 1 : photon is assumed resolved; i.e. can only interact through
its constituent quarks and gluons, in processes such as
ISUB = 11 - 31.
Note: currently it is not possible to mix the two kinds of
photon-induced processes in one run.
MSTP(21) : (D=1) nature of fermion-fermion scatterings simulated in
process 10 by t-channel exchange.
= 0 : all off.
= 1 : full mixture of gamma/Z neutral current and W charged current.
= 2 : gamma neutral current only.
= 3 : Z neutral current only.
= 4 : gamma/Z neutral current only.
= 5 : W charged current only.
MSTP(22) : (D=0) special override of normal Q^2 definition used for
maximum of parton shower evolution. This option only affects
processes 10 and 83 (Deep Inelastic Scattering) and only in
lepton-hadron events.
= 0 : use the scale as given in MSTP(32).
= 1 : use the DIS Q^2 scale, i.e. -that.
= 2 : use the DIS W^2 scale, i.e. (1-x)/x * (-that).
= 3 : use the DIS Q*W scale, i.e. sqrt((1-x)/x) * (-that).
Note: in all of these alternatives, a multiplicative factor is
introduced by PARP(67) and PARP(71).
MSTP(23) : (D=0) for Deep Inelastic Scattering processes (10 and 83)
this option allows the x and Q^2 of the original hard scattering
to be retained by the final electron.
= 0 : no correction procedure, i.e. x and Q^2 of the scattered
electron differ from the originally generated x and Q^2.
= 1 : post facto correction, i.e. the change of electron momentum
by initial and final QCD radiation, primordial k_T and beam
remnant treatment is corrected for by a shuffling of momentum
between the electron and hadron side in the final state.
Note: the correction procedure will fail for a fraction of the
events, which thus are thrown (and new ones generated in their
place). The correction option is not unambiguous, and should
not be taken too seriously. For very small Q^2 values, the x
is not exactly preserved even after this procedure.
MSTP(31) : (D=1) parametrization of total and elastic cross-sections,
nuclear slope parameter B and curvature C [Blo85].
= 1 : Block-Cahn fit 1 for cross-section, fit 1 for slope parameter.
= 2 : Block-Cahn fit 2 for cross-section, fit 1 for slope parameter.
= 3 : Block-Cahn fit 3 for cross-section, fit 1 for slope parameter.
= 4 : Block-Cahn fit 6 for cross-section, fit 2 for slope parameter.
= 5 : Block-Cahn fit 8 for cross-section, fit 2 for slope parameter.
Note: sets 1-3 for cross-section and set 1 for slope parameter are
fits excluding recent measurements from Spp~S, whereas sets 4-5
for cross-section and set 2 for slope parameter are fits
including the Spp~S measurements.
MSTP(32) : (D=2) Q^2-definition in hard scattering for 2 -> 2 processes;
for resonance production Q^2 is always chosen to be m_res^2, where
m_res is the mass of the resonance.
= 1 : Q^2 = 2*shat*that*uhat/(shat^2 + that^2 + uhat^2).
= 2 : Q^2 = (m_T1^2 + m_T2^2)/2.
= 3 : Q^2 = min(-that, -uhat).
= 4 : Q^2 = shat.
= 5 : Q^2 = -that.
MSTP(33) : (D=0) inclusion of K factors in hard cross-sections for
parton-parton interactions (i.e. for incoming quarks and gluons).
= 0 : none, i.e. K = 1.
= 1 : a common K factor is used, as stored in PARP(31).
= 2 : separate factors are used for ordinary (PARP(31)) and colour
annihilation graphs (PARP(32)).
= 3 : A K factor is introduced by a shift in the alpha_strong Q^2
argument, alpha_strong = alpha_strong(PARP(33)*Q^2).
MSTP(34) : (D=0) use of interference term in matrix elements for QCD
processes.
= 0 : excluded (i.e. string-inspired matrix elements).
= 1 : included (i.e. conventional QCD matrix elements).
Note: for the option MSTP(34)=1, i.e., interference terms included,
these are divided between the different possible colour
configurations according to the pole structure of the
(string-inspired) matrix elements for the different colour
configurations.
MSTP(35) : (D=0) threshold behaviour for heavy flavour production, i.e.
ISUB = 81, 82, 84, 85, and also for Z and Z' decays. The nonstandard
options are mainly intended for top, but can be used, with less
theoretical reliability, also for charm and bottom (for Z and Z'
only top and heavier flavours are affected). See [Fad89] for
details.
= 0 : naive lowest order matrix element behaviour.
= 1 : enhancement or suppression close to threshold, according to
the colour structure of the process, see [Fad89] for the
expressions actually used. The alpha_strong value appearing
in the threshold factor (which is not the same as the
alpha_strong of the lowest order 2 -> 2 process) is taken to
be fixed at the value given in PARP(35). The threshold factor
used in an event is stored in PARI(81).
= 2 : as =1, but the alpha_strong value appearing in the threshold
factor is taken to be running, with argument
Q^2 = m_Q sqrt( (m-hat - 2m_Q)^2 + Gamma_Q^2).
Here m_Q is the nominal heavy quark mass, Gamma_Q is the width
of the heavy quark mass distribution, and m-hat is the invariant
mass of the heavy quark pair. The Gamma_Q value has to be stored
by the user in PARP(35). The regularization of alpha_strong at
low Q^2 is given by MSTP(36).
MSTP(36) : (D=2) regularization of alpha_strong in the limit Q^2 -> 0
for the threshold factor obtainable in the MSTP(35) = 2 option;
see MSTU(115) for a list of the possibilities.
MSTP(37) : (D=1) inclusion of running quark masses in Higgs production
(q + qb -> H0) and decay (H0 -> q + qb) couplings.
= 0 : not included, i.e. fix quark masses are used according to
the values in the PMAS array.
= 1 : included, with running starting from the value given in
the PMAS array, at a Q (i.e. Higgs mass) scale given by
PARP(37) times the quark mass itself.
MSTP(38) : (D=5) handling of quark loop masses in the box graphs
g + g -> gamma + gamma and g + g -> g + gamma.
= 0 : the program will for each flavour automatically choose
the massless approximation for light quarks and the full
massive formulae for heavy quarks, with the dividing line
between light and heavy quarks dependent on the actual
shat scale.
>= 1 : the program will use the massless approximation throughout,
assuming the presence of MSTP(38) effectively massless quark
species (however, at most 8). Normally one would use = 5 for
the inclusion of all quarks up to bottom, and = 6 to include
top as well.
WARNING: for =0, numerical instabilities may arise for scattering
at small angles. Users are therefore recommended in this case
to set CKIN(27) and CKIN(28) so as to exclude the range of
scattering angles that are not of interest anyway.
MSTP(41) : (D=1) master switch for all resonance decays (Z0, W+-, H0,
Z'0, W'+-, H'0, A0, H+-, LQ, R).
= 0 : off.
= 1 : on.
Note: also for MSTP(41) = 1 it is possible to switch off the decays
of specific resonances by using the MDCY(KC,1) switches in
JETSET. However, since the MDCY values are overwritten in the
PYINIT call, individual resonances should be switched off after
the PYINIT call.
WARNING: leptoquark decays must not be switched off if one later
on intends to let it decay (with LUEXEC); see 2.3.10.
MSTP(42) : (D=1) mass treatment in 2 -> 2 processes, where the final
state resonances have finite width (see PARP(41)).
(Does not apply for the production of a single s-channel resonance,
where the mass appears explicitly in the cross-section of the
process, and thus is always selected with width.)
= 0 : particles are put on mass-shell.
= 1 : mass generated according to a Breit-Wigner.
MSTP(43) : (D=3) treatment of Z0/gamma* interference in matrix elements.
So far only implemented in subprocesses 1 and 22; in other processes
what is called a Z0 is really a Z0 only, without the gamma* piece.
= 1 : only gamma* included.
= 2 : only Z0 included.
= 3 : complete Z0/gamma* structure (with interference) included.
MSTP(44) : (D=7) treatment of Z'0/Z0/gamma* interference in matrix
elements.
= 1 : only gamma* included.
= 2 : only Z0 included.
= 3 : only Z'0 included.
= 4 : only Z0/gamma* (with interference) included.
= 5 : only Z'0/gamma* (with interference) included.
= 6 : only Z'0/Z0 (with interference) included.
= 7 : complete Z'0/Z0/gamma* structure (with interference) included.
MSTP(45) : (D=2) treatment of WW -> WW structure (ISUB = 77).
= 1 : only W+W+ -> W+W+ and W-W- -> W-W- included.
= 2 : only W+W- -> W+W- included.
= 3 : all charge combinations WW -> WW included.
MSTP(46) : (D=1) treatment of VV -> V'V' structures (ISUB = 71 - 77).
= 0 : only s-channel Higgs exchange included (where existing). With
this option, subprocesses 71 - 72 and 76 - 77 will essentially
be equivalent to subprocesses 5 and 8, respectively, with the
proper decay channels (i.e. only Z0Z0 or W+W-) set via MDME.
The description obtained for subprocesses 5 and 8 in this case
is more sophisticated, however; see section 2.3.6.
= 1 : all graphs contributing to VV -> V'V' processes are included.
= 2 : only graphs not involving Higgs exchange (either in s, t or u
channel) are included; this option then gives the naive
behaviour if no Higgs exists, including unphysical unitarity
violations at high energies.
= 3 : the strongly interacting Higgs-like model of Dobado, Herrero
and Terron [DHT90] with Pade unitarization. Note that to use
this option it is necessary to set the Higgs mass to a large
number like 20 TeV (i.e. PMAS(25,1)=20000.). The parameter nu
is stored in PARP(44), but should not have to be changed.
= 4 : as =3, but with K-matrix unitarization.
= 5 : the strongly interacting QCD-like model of Dobado, Herrero
and Terron [DHT90] with Pade unitarization. The parameter nu
is stored in PARP(44), but should not have to be changed.
The effective techni-rho mass in this model is stored in
PARP(45); by default it is 2054 GeV, which is the expected
value for three technicolours, based on scaling up the ordinary
rho mass appropriately.
= 6 : as =5, but with K-matrix unitarization. While PARP(45) still
is a parameter of the model, this type of unitarization does
not give rise to a resonance at a mass of PARP(45).
MSTP(47) : (D=1) (C) angular orientation of decay products of resonances
(Z, W, H, Z', W', R), either when produced singly or in pairs (also
from a H0 decay), or in combination with a single quark, gluon or
photon.
= 0 : independent decay of each resonance, isotropic in CM frame of
the resonance.
= 1 : correlated decay angular distributions according to proper
matrix elements, to the extent these are known.
MSTP(51) : (D=1) choice of proton structure functions; see also
MSTP(52).
= 1 : EHLQ set 1 (1986 updated version).
= 2 : EHLQ set 2 (1986 updated version).
= 3 : Duke-Owens set 1.
= 4 : Duke-Owens set 2.
= 5 : Morfin-Tung set 1.
= 6 : Morfin-Tung set 2.
= 7 : Morfin-Tung set 3.
= 8 : Morfin-Tung set 4.
= 9 : Gluck-Reya-Vogt LO set.
= 10 : Gluck-Reya-Vogt HO set.
Note: since all parametrizations have some region of applicability,
the structure functions are assumed frozed below the lowest Q^2
and above the highest Q^2 covered by the parametrizations. The
extrapolation to low x is covered by PARP(51).
MSTP(52) : (D=1) choice of proton structure function library.
= 1 : the internal PYTHIA one, with structure functions according
to MSTP(51) above.
= 2 : the PDFLIB one [Plo91], with the PDFLIB MODE (set) number
to be given in MSTP(51).
= 3 : the PAKPDF one [Cha91], with the PAKPDF 100*IPAR+ISET number
to be given in MSTP(51).
NOTE: to make use of options 2 and 3, it is necessary to enable
the calls to PDFLIB and PAKPDF, respectively, in subroutines
PYINIT and PYSTFU. What to do is described inline in the
Fortran file; look for lines with C! in the first two columns.
WARNING: For external structure function libraries, PYTHIA does not
check whether MSTP(51) corresponds to a valid code, neither if
special x and Q2 restrictions exist for a given set, such that
crazy values could be returned. This puts an extra
responsibility on the user.
MSTP(53) : (D=1) choice of pion structure functions; see also
MSTP(54).
= 1 : Owens set 1.
= 2 : Owens set 2.
MSTP(54) : (D=1) choice of pion structure function library.
= 1 : the internal PYTHIA one, with structure functions according
to MSTP(53) above.
= 2 : the PDFLIB one [Plo91], with the PDFLIB MODE (set) number
to be given in MSTP(53).
NOTE: to make use of option 2, it is necessary to enable the calls
to PDFLIB in subroutines PYINIT and PYSTFU. What to do is
described inline in the Fortran file; look for lines with
C! in the first two columns.
WARNING: For external structure function libraries, PYTHIA does not
check whether MSTP(53) corresponds to a valid code, neither if
special x and Q2 restrictions exist for a given set, such that
crazy values could be returned. This puts an extra
responsibility on the user.
MSTP(55) : (D=1) choice of photon structure functions; see also
MSTP(52). Currently there are no alternatives available for
this switch.
= 1 : Drees-Grassie.
MSTP(56) : (D=1) choice of photon structure function library.
= 1 : the internal PYTHIA one, with structure functions according
to MSTP(53) above.
= 3 : the PHOPDF one [Cha91a], with the PHOPDF 100*IPAR+ISET number
to be given in MSTP(55).
NOTE: to make use of option 3, it is necessary to enable the calls
to PHOPDF, in subroutines PYINIT, PYSTFU and PYSTEL. What to
do is described inline in the Fortran file; look for lines with
C! in the first two columns.
WARNING: For external structure function libraries, PYTHIA does not
check whether MSTP(53) corresponds to a valid code, neither if
special x and Q2 restrictions exist for a given set, such that
crazy values could be returned. This puts an extra
responsibility on the user.
MSTP(57) : (D=1) choice of Q^2 dependence in structure functions.
This option only applies to the proton, pion and photon structure
functions that come with PYTHIA, and not to externally linked ones.
= 0 : structure functions are evaluated at nominal lower cutoff
value Q_0^2, i.e. are made Q^2-independent.
= 1 : the parametrized Q^2 dependence is used.
MSTP(58) : (D=min(6,2*MSTP(1))) maximum number of quark flavours used
in structure functions, and thus also for initial state spacelike
showers. If some distributions (notably t) are absent in the
parametrization selected in MSTP(51), these are obviously
automatically excluded.
MSTP(61) : (D=1) (C) master switch for initial state QCD and QED
radiation.
= 0 : off.
= 1 : on.
MSTP(62) : (D=2) (C) level of coherence imposed on the spacelike parton
shower evolution.
= 1 : none, i.e. neither Q^2 values nor angles need be ordered.
= 2 : Q^2 values at branches are strictly ordered, increasing
towards the hard interaction.
= 3 : Q^2 values and opening angles of emitted (on-mass-shell or
timelike) partons are both strictly ordered, increasing
towards the hard interaction.
MSTP(63) : (D=2) (C) structure of associated timelike showers, i.e.
showers initiated by emission off the incoming spacelike partons.
= 0 : no associated showers are allowed, i.e. emitted partons are
put on mass-shell.
= 1 : a shower may evolve, with maximum allowed timelike virtuality
set by the phase space only.
= 2 : a shower may evolve, with maximum allowed timelike virtuality
set by phase space or by PARP(71) times the Q^2 value of the
spacelike parton created in the same vertex, whichever is the
stronger constraint.
= 3 : a shower may evolve, with maximum allowed timelike virtuality
set by phase space or by the requirement of angular ordering,
whichever is the stronger constraint. (Not yet implemented!)
MSTP(64) : (D=2) (C) choice of alpha_strong and Q^2 scale in spacelike
parton showers.
= 0 : alpha-strong is taken to be fix at the value PARU(111).
= 1 : first order running alpha_strong with argument PARP(63)*Q^2.
= 2 : first order running alpha_strong with argument
PARP(64)*k_T^2 = PARP(64)*(1-z)*Q^2.
MSTP(65) : (D=1) (C) treatment of soft gluon emission in spacelike
parton shower evolution.
= 0 : soft gluons are entirely neglected.
= 1 : soft gluon emission is resummed and included together
with the hard radiation as an effective z shift.
MSTP(71) : (D=1) (C) master switch for final state QCD and QED
radiation.
= 0 : off.
= 1 : on.
Note: additional switches (e.g. for conventional/coherent showers)
are available in JETSET, see MSTJ(41) - MSTJ(45), PARJ(81) and
PARJ(82). In particular, to allow photon emission off quarks
and leptons, it is necessary to put MSTJ(41)=2.
MSTP(81) : (D=1) master switch for multiple interactions.
= 0 : off.
= 1 : on.
MSTP(82) : (D=1) structure of multiple interactions. For QCD processes,
used down to p_T values below p_Tmin, it also affects the choice
of structure for the one hard/semihard interaction.
= 0 : simple two-string model without any hard interactions.
= 1 : multiple interactions assuming the same probability in all
events, with an abrupt p_Tmin cutoff at PARP(81).
= 2 : multiple interactions assuming the same probability in all
events, with a continuous turnoff of the cross-section at
p_T0 = PARP(82).
= 3 : multiple interactions assuming a varying impact parameter
and a hadronic matter overlap consistent with a Gaussian matter
distribution, with a continuous turnoff of the cross-section at
p_T0 = PARP(82).
= 4 : multiple interactions assuming a varying impact parameter
and a hadronic matter overlap consistent with a double Gaussian
matter distribution given by PARP(83) and PARP(84), with a
continuous turnoff of the cross-section at p_T0 = PARP(82).
Note 1: For MSTP(82) >= 2 and CKIN(3) > PARP(82), cross-sections
given with PYSTAT(1) may be somewhat too large, since (for
reasons of efficiency) the probability factor that the hard
interaction is indeed the hardest in the event is not
included in the cross-sections. It is included in the event
selection, however, so the events generated are correctly
distributed. For CKIN(3) values a couple of times larger than
PARP(82) this ceases to be a problem.
Note 2: The PARP(81) and, in particular, PARP(82) values are
sensitive to the choice of structure functions, Lambda_QCD,
etc., in the sense that a change in the latter variables
leads to a net change in the multiple interaction rate, which
has to be compensated by a retuning of PARP(81) or PARP(82)
if one wants to keep the net multiple interaction structure the
same. The default PARP(81) value is consistent with the other
default values give, i.e. EHLQ set 1 structure functions etc.
When options MSTP(82) = 2 - 4 are used, the default PARP(82)
value is to be used in conjunction with MSTP(2) = 2 and
MSTP(33) = 3. These switches should be set by the user.
MSTP(83) : (D=100) number of Monte Carlo generated phase space points
per bin (whereof there are 20) in the initialization (in PYINMU)
of multiple interactions for MSTP(82) >= 2.
MSTP(91) : (D=1) (C) primordial k_T distribution in hadron. See MSTP(93)
for photon.
= 0 : no primordial k_T.
= 1 : Gaussian, width given in PARP(91), upper cutoff in PARP(93).
= 2 : exponential, width given in PARP(92), upper cutoff in
PARP(93).
MSTP(92) : (D=4) (C) energy partitioning in hadron remnant. The energy
fraction chi taken by one of the two objects, with conventions as
described for PARP(94) - PARP(98), is chosen according to the
different distributions below. Here c_min = 2*m_quark/sqrt(s) =
0.6 GeV/sqrt(s).
= 1 : 1 for meson, 2*(1-chi) for baryon, i.e. simple counting rules.
= 2 : (k+1)*(1-chi)^k, with k as given in PARP(94) - PARP(97).
= 3 : as =2 for remnant splitting into hadron plus jet, but
proportional to (1-chi)^k/(chi^2+c_min^2)^0.25 for remnant
splitting into two jets, with k given by PARP(94) or PARP(96).
= 4 : as =2 for remnant splitting into hadron plus jet, but
proportional to (1-chi)^k/(chi^2+c_min^2)^0.5 for remnant
splitting into two jets, with k given by PARP(94) or PARP(96).
= 5 : as =2 for remnant splitting into hadron plus jet, but
proportional to (1-chi)^k/(chi^2+c_min^2)^(PARP(98)/2) for
remnant splitting into two jets, with k given by PARP(94) or
PARP(96).
MSTP(93) : (D=1) (C) primordial k_T distribution in photon, either it
is one of the incoming particles or inside an electron.
= 0 : no primordial k_T.
= 1 : Gaussian, width given in PARP(99), upper cutoff in PARP(100).
= 2 : exponential, width given in PARP(99), upper cutoff in
PARP(100).
= 3 : powerlike of the type d(k_T^2)/(PARP(99)^2+k_T^2)^2, with
upper k_T cutoff in PARP(100).
= 4 : powerlike of the type d(k_T^2)/(PARP(99)^2+k_T^2), with
upper k_T cutoff in PARP(100).
= 5 : powerlike of the type d(k_T^2)/(PARP(99)^2+k_T^2), with
upper k_T cutoff given by the p_T of the hard process or by
PARP(100), whichever is smaller.
NOTE: for options 1 and 2 the PARP(100) value is of minor
importance, once PARP(100) >> PARP(99). However, options 3 and
4 correspond to distributions with infinite if the
k_T spectrum is not cut off, and therefore the PARP(100) value
is as important for the overall distribution as is PARP(99).
MSTP(101) : (D=1) (C) structure of diffractive system.
= 1 : forward moving diquark + interacting quark.
= 2 : forward moving diquark + quark joined via interacting gluon
('hairpin' configuration).
MSTP(111) : (D=1) (C) master switch for fragmentation and decay, as
obtained with a LUEXEC call.
= 0 : off.
= 1 : on.
= -1 : only choose kinematical variables for hard scattering,
i.e. no jets are defined. This is useful e.g. to calculate
cross-sections (by Monte Carlo integration) without wanting
to simulate events; information obtained with PYSTAT(1)
will be correct.
MSTP(112) : (D=1) (C) cuts on partonic events; only affects an
exceedingly tiny fraction of events.
= 0 : no cuts (can be used only with independent fragmentation,
at least in principle).
= 1 : string cuts (as normally required for fragmentation).
MSTP(113) : (D=1) (C) recalculation of energies of partons from their
momenta and masses, to be done immediately before fragmentation, to
compensate in parts for some numerical problems appearing at high
energies.
= 0 : not performed.
= 1 : performed.
MSTP(121) : (D=0) calculation of kinematics selection coefficients and
differential cross-section maxima for subprocesses included
(by user or default).
= 0 : not known; to be calculated at initialization.
= 1 : not known; to be calculated at initialization; however, the
maximum value then obtained is to be multiplied by PARP(121)
(this may be useful if a violation factor has been observed
in a previous run of the same kind).
= 2 : known; kinematics selection coefficients stored by user in
COEF(ISUB,J) (J = 1 - 20) in /PYINT2/ and maximum of the
corresponding differential cross-section times Jacobians in
XSEC(ISUB,1) in /PYINT5/. This is to be done for each included
subprocess ISUB before initialization, with the sum of all
XSEC(ISUB,1) values, except for ISUB = 95, stored in XSEC(0,1).
MSTP(122) : (D=1) initialization and differential cross-section
maximization printout (see also MSTP(127)).
= 0 : none.
= 1 : short message.
= 2 : detailed message, including full maximization.
MSTP(123) : (D=2) reaction to violation of maximum differential cross-
section.
= 0 : stop generation, print message.
= 1 : continue generation, print message for each subsequently
larger violation.
= 2 : as 1, but also increase value of maximum.
MSTP(124) : (D=1) (C) frame for presentation of event.
= 1 : as specified in PYINIT.
= 2 : CM frame of incoming particles.
MSTP(125) : (D=1) (C) documentation of partonic process, see section
'The event record' for details.
= 0 : only list ultimate string/particle configuration.
= 1 : additionally list short summary of the hard process.
= 2 : list complete documentation of intermediate steps of parton
shower evolution.
MSTP(126) : (D=20) number of lines in the beginning of event record
that are reserved for event history information; see section 2.7.
This value should never be reduced, but may be increased at a
later date if more complicated processes are included.
MSTP(127) : (D=1) writing of header (version number and last date of
change) on output file.
= 0 : not done.
= 1 : header is written at first PYINIT call, at which time
MSTP(127) is set =0.
MSTP(128) : (D=0) storing of copy of resonance decay products in
the documentation section of the event record, and mother
pointer (K(I,3)) relation of the actual resonance decay
products (stored in the main section of the event record)
to the documentation copy.
= 0 : products are stored also in the documentation section,
and each product stored in the main section points back
to the corresponding entry in the documentation section.
= 1 : products are stored also in the documentation section,
but the products stored in the main section point back to
the decaying resonance copy in the main section.
= 2 : products are not stored in the documentation section,
the products stored in the main section point back to the
the decaying resonance copy in the main section.
MSTP(129) : (D=10) for the maximization of 2 -> 3 processes
(ISET(ISUB)=5) each phase space point in tau, y* and tau' is
tested MSTP(129) times in the other dimensions (at randomly
selected points) to determine the effective maximum in the
(tau, y*, tau') point.
MSTP(131) : (D=0) master switch for pileup events, i.e. several
independent hadron-hadron interactions generated in the same
bunch-bunch crossing, with the events following one after the
other in the event record.
= 0 : off, i.e. only one event is generated at a time.
= 1 : on, i.e. several events are allowed in the same event record.
Information on the processes generated may be found in MSTI(41)
- MSTI(50).
MSTP(132) : (D=4) the processes that are switched on for pileup events.
The first event may be set up completely arbitrarily, using the
switches in the PYSUBS commonblock, while all the subsequent events
have to be of one of the 'inclusive' processes which dominate the
cross-section, according to the options below. It is thus not
possible to generate two rare events in the pileup option.
= 1 : low-p_T processes (ISUB = 95) only. The low-pT model actually
used, both in the hard event and in the pileup events, is the
one set by MSTP(81) etc. This means that implicitly also
high-pT jets can be generated in the pileup events.
= 2 : low-p_T + double diffractive processes (ISUB = 95 and 93).
= 3 : low-p_T + double diffractive + single diffractive processes
(ISUB = 95, 93 and 92).
= 4 : low-p_T + double diffractive + single diffractive + elastic
processes, together corresponding to the full hadron-hadron
cross-section (ISUB = 95, 93, 92 and 91).
MSTP(133) : (D=0) multiplicity distribution of pileup events.
= 0 : selected by user, before each PYEVNT call, by giving the
MSTP(134) value.
= 1 : a Poissonian multiplicity distribution in the total number
of pileup events. This is the relevant distribution if the
switches set for the first event in PYSUBS give the same
subprocesses as are implied by MSTP(132). In that case the
mean number of events per beam crossing is nbar =
(sum of cross-section for allowed processes) * PARP(131).
Since bunch crossing which do not give any events at all
(probability exp(-nbar)) are not simulated, the actual average
number per PYEVNT call is = nbar/(1-exp(-nbar)).
= 2 : a biased distribution, as is relevant when one of the
events to be generated is assumed to belong to an event class
with a cross-section much smaller than the total hadronic
cross-section. If sigma_hard is the cross-section for this
rare process (or the sum of the cross-sections of several rare
processes) and sigma_soft the cross-section for the processes
allowed by MSTP(132), then define nbar = sigma_soft * PARP(131)
and f = sigma_hard/sigma_soft. The probability that a bunch
crossing will give i events is then
Prob(i) = exp(-nbar) * nbar^i/i! * f * i
i.e. the naive Poissonian is suppressed by a factor f since
one of the events will be rare rather than frequent, but
enhanced by a factor i since any of the i events may be the
rare one. Only beam crossings which give at least one event
of the required rare type are simulated, and the distribution
above normalized accordingly.
Note: for practical reasons, it is required that nbar < 120, i.e.
that an average beam crossing does not contain more than 120
pileup events. The multiplicity distribution is truncated
above 200, or when the probability for a multiplicity has fallen
below 10^(-6), whichever occurs sooner. Also low multiplicities
with probabilities below 10^(-6) are truncated. See also
PARI(91) - PARI(93).
MSTP(134) : (D=1) a user selected multiplicity, i.e. total number
of pileup events, to be generated in the next PYEVNT call.
May be reset for each new event, but must be in the range
1 <= MSTP(134) <= 200.
MSTP(141) : (D=0) calling of PYKCUT in event generation chain, for
inclusion of user-specified cuts.
= 0 : not called.
= 1 : called.
MSTP(142) : (D=0) calling of PYEVWT in the event generation chain,
either to give weighted events or to modify standard cross-sections.
See PYEVWT description in section 2.1 for further details.
= 0 : not called.
= 1 : called; the distribution of events among subprocesses and
in kinematics variables is modified by the factor WTXS, set
by the user in the PYEVWT call, but events come with a
compensating weight PARI(10)=1./WTXS, such that total
cross-sections are unchanged.
= 2 : called; the cross-section itself is modified by the factor
WTXS, set by the user in the PYEVWT call.
MSTP(181) : (R) PYTHIA version number.
MSTP(182) : (R) PYTHIA subversion number.
MSTP(183) : (R) last year of change for PYTHIA.
MSTP(184) : (R) last month of change for PYTHIA.
MSTP(185) : (R) last day of change for PYTHIA.
PARP(1) : (D=0.25 GeV) nominal Lambda-QCD used in running alpha-strong
for hard scattering (see MSTP(3)).
PARP(2) : (D=10. GeV) lowest CM energy for the event-as-a-whole that
the program will accept to simulate.
PARP(11), PARP(12) : (D=0.,-0.67) coefficients of a beamstrahlung-
induced effective electron structure function, see MSTP(11) = 2.
To ensure proper normalization, it is necessary to have
PARP(11) and PARP(12) both > -1; in practice, values smaller
than roughly -0.95 are not recommended.
PARP(13) : (D=25. GeV^2) Q_max^2 scale, to be set by user for defining
maximum scale allowed for photoproduction when using the option
MSTP(13)=2.
PARP(14) : (D=0.01) in the numerical integration of quark and gluon
structure functions inside an electron, the successive halvings
of evaluation point spacing is interrupted when two values agree
in relative size, abs(new-old)/(new+old), to better than PARP(14).
There are hardwired lower and upper limits of 2 and 8 halvings,
respectively.
PARP(31) : (D=1.5) common K factor multiplying the differential
cross-section for hard parton-parton processes when MSTP(33) = 1
or 2, with the exception of colour annihilation graphs in the
latter case.
PARP(32) : (D=2.0) special K-factor multiplying the differential
cross-section in hard colour annihilation graphs, including
resonance production, when MSTP(33) = 2.
PARP(33) : (D=0.075) this factor is used to multiply the ordinary Q^2
scale in alpha_strong at the hard interaction for MSTP(33) = 3.
The effective K-factor thus obtained is in accordance with the
results in [Ell86].
PARP(35) : (D=0.20) fix alpha_strong value used in heavy flavour
threshold factor when MSTP(35) = 1.
PARP(36) : (D=0. GeV) the width Gamma_Q for the heavy flavour studied
in processes ISUB = 81 or 82; to be used for the threshold factor
when MSTP(35) = 2.
PARP(37) : (D=2.) for MSTP(37)=1 this regulates the point at which the
naive fix quark mass in Higgs couplings is assumed defined;
specifically the running quark mass is assumed to coincide with
fix one at PARP(37) times the fix quark mass, i.e.
m_running(PARP(37)*m_fix) = m_fix.
PARP(41) : (D=0.010 GeV) in the process of generating mass for
resonances, and optionally to force that mass to be in a given
range, only resonances with a total width in excess of PARP(41)
are generated according to a Breit-Wigner shape (if allowed by
MSTP(42)), while more narrow resonances are put on mass-shell.
PARP(42) : (D=2. GeV) minimum mass of resonances assumed allowed when
evaluating total width of H0 to Z0 + Z0 or W+ + W- for cases
when the H0 is so light that (at least) one Z/W is forced to be
off mass-shell. Also generally used as safety check on minimum
mass of resonance.
PARP(43) : (D=0.10) precision parameter used in numerical integration
of width into channel with at least one daughter off mass-shell.
PARP(44) : (D=1000.) the nu papameter of the stronly interacting Z/W
model of Dobado, Herrero and Terron [DHT90].
PARP(45) : (D=2054. GeV) the effective techni-rho mass parameter of the
strongly interacting model of Dobado, Herrero and Terron [DHT90];
see MSTP(46)=5. On physical grounds it should not be chosen smaller
than about 1 TeV or larger than the about default value.
PARP(51) : (D=1.) if structure functions for light flavours have to be
extrapolated to x values lower than covered by the parametrizations,
an x^(-PARP(51)) behaviour is assumed in that region. This option
only applies for the EHLQ and GRV proton structure functions that
are internal to PYTHIA.
PARP(61) : (D=0.25 GeV) (C) Lambda value used in spacelike parton
shower (see MSTP(64)). This value may be overwritten, see
MSTP(3).
PARP(62) : (D=1. GeV) (C) effective cutoff Q or k_T value (see
MSTP(64)), below which spacelike parton showers are not evolved.
PARP(63) : (D=0.25) (C) in spacelike shower evolution the virtuality
Q^2 of a parton is multiplied by PARP(63) for use as a scale in
alpha_strong and structure functions when MSTP(64) = 1.
PARP(64) : (D=1.) C) in spacelike parton shower evolution the transverse
momentum-squared evolution scale k_T^2 is multiplied by PARP(64)
for use as a scale in alpha_strong and structure functions when
MSTP(64) = 2.
PARP(65) : (D=2. GeV) (C) effective minimum energy (in CM frame) of
timelike or on-shell parton emitted in spacelike shower; see also
PARP(66).
PARP(66) : (D=0.001) (C) effective lower cutoff on 1-z in spacelike
showers, in addition to the cut implied by PARP(65).
PARP(67) : (D=4.) (C) the Q^2 scale of the hard scattering (see
MSTP(32)) is multiplied by PARP(67) to define the maximum parton
virtuality allowed in spacelike showers. This does not apply to
s-channel resonances, where the maximum virtuality is set by m^2.
PARP(68) : (D=1E-6) lower Q^2 cutoff for QED spacelike showers.
PARP(71) : (D=4.) (C) the Q^2 scale of the hard scattering (see
MSTP(32)) is multiplied by PARP(71) to define the maximum parton
virtuality allowed in timelike showers. This does not apply to
s-channel resonances, where the maximum virtuality is set by m^2.
PARP(81) : (D=1.45 GeV/c) effective minimum transverse momentum p^T_min
for multiple interactions with MSTP(82) = 1.
Note: The default value is changed compared to the one used up till
version 5.3 (1.60 GeV/c), not because of any physics change,
but simply because the treatment of alpha_strong and Lambda
at flavour thresholds has been improved, leading to a smaller
alpha_strong at small p_T values, which has to be compensated.
PARP(82) : (D=1.70 GeV) regularization scale p^T_0 of the transverse
momentum spectrum for multiple interactions with MSTP(82) >= 2.
Note: The default value is changed compared to the one used up till
version 5.3 (2.00 GeV/c), not because of any physics change,
but simply because the treatment of alpha_strong and Lambda
at flavour thresholds has been improved, leading to a smaller
alpha_strong at small p_T values, which has to be compensated.
PARP(83), PARP(84) : (D=0.5, 0.2) parameters of an assumed double
Gaussian matter distribution inside the colliding hadrons for
MSTP(82) = 4, of the form
(1-PARP(83))*exp(-r^2) + (PARP(83)/PARP(84)^3)*exp(-r^2/PARP(84)^2)
i.e. with a core of radius PARP(84) of the main radius and
containing a fraction PARP(83) of the total hadronic matter.
PARP(85) : (D=0.33) probability that an additional interaction in the
multiple interaction formalism gives two gluons, with colour
connections to 'nearest neighbours' in momentum space.
PARP(86) : (D=0.66) probability that an additional interaction in the
multiple interaction formalism gives two gluons, either as
described in PARP(85) or as a closed gluon loop. Remaining fraction
is supposed to consist of quark-antiquark pairs.
PARP(87), PARP(88) : (D=0.7, 0.5) in order to account for an assumed
dominance of valence quarks at low transverse momentum scales, a
probability is introduced that a gg-scattering according to naive
cross-section is replaced by a qq~ one; this is used only for
MSTP(82) >= 2. The probability is parametrized as
prob = PARP(87)*(1.-(p_T^2/(p_T^2+PARP(88)*(PARP(82))^2))^2).
PARP(91) : (D=0.44 GeV/c) (C) width of Gaussian primordial k_T
distribution inside hadron for MSTP(91) = 1.
PARP(92) : (D=0.44 GeV/c) (C) width of exponential primordial k_T
distribution inside hadron for MSTP(91) = 2.
PARP(93) : (D=2. GeV/c) (C) upper cutoff for primordial k_T
distribution inside hadron.
PARP(94) : (D=1.) (C) for MSTP(92) >= 2 this gives the value of the
parameter k for the case when a pion remnant is split into two
fragments (which is which is chosen at random).
PARP(95) : (D=0.) (C) for MSTP(92) >= 2 this gives the value of the
parameter k for the case when a pion remnant is split into a meson
and a spectator fragment jet, with chi giving the energy fraction
taken by the meson.
PARP(96) : (D=3.) (C) for MSTP(92) >= 2 this gives the value of the
parameter k for the case when a nucleon remnant is split into a
diquark and a quark fragment, with chi giving the energy fraction
taken by the quark jet.
PARP(97) : (D=1.) (C) for MSTP(92) >= 2 this gives the value of the
parameter k for the case when a nucleon remnant is split into a
baryon and a quark jet or a meson and a diquark jet, with chi giving
the energy fraction taken by the quark jet or meson, respectively.
PARP(98) : (D=0.75) (C) for MSTP(92) = 5 this gives the power of an
assumed basic 1/chi^PARP(98) behaviour in the splitting
distribution.
PARP(99) : (D=0.44 GeV/c) (C) width of primordial k_T distribution
inside photon; exact meaning depends on MSTP(93) value chosen.
PARP(100) : (D=2. GeV/c) (C) upper cutoff for primordial k_T
distribution inside photon.
PARP(101) : (D=-0.02 GeV^2) minimum value of t-hat in (diffractive and)
elastic scattering.
PARP(111) : (D=2. GeV) used to define the minimum invariant mass of
the remnant hadronic system (i.e. when interacting partons have been
taken away), together with original hadron masses and extra parton
masses.
PARP(121) : (D=1.) the maxima obtained at initial maximization are
multiplied by this factor if PARP(121)=1; typically PARP(121) would
be given as the product of the violation factors observed (i.e the
ratio of final maximum value to the initial maximum value) for the
given process(es).
PARP(122) : (D=0.4) fraction of total probability that is shared
democratically between the COEF coefficients open for the given
variable, with remaining fraction distributed according to the
optimization results of PYMAXI.
PARP(131) : (D=0.01 mb^(-1)) in the pileup events scenario (see
MSTP(131) - MSTP(133)), PARP(131) gives the assumed luminosity per
bunch-bunch crossing, i.e. if a subprocess has a cross-section
sigma, the average number of events of this type per bunch-bunch
crossing is nbar = sigma * PARP(131). PARP(131) may be obtained
by dividing the integrated luminosity over a given time (1 s, say)
by the number of bunch-bunch crossings that this corresponds to.
Since the program will not generate more than 200 pileup events,
the initialization procedure will crash if nbar is above 120.
______________________________________________________________________
2.6. General Event Information
When an event is generated with PYEVNT, some information on this event
is stored in the MSTI and PARI arrays of the LUDATA commonblock (often
copied directly from the internal MINT and VINT variables). Further
information is stored in the complete event record; see section 2.7.
Part of the information is only relevant for some subprocesses; by
default everything irrelevant is set 0. Kindly note that, like the CKIN
constraints described in section 2.4, kinematical variables normally
(where it is not explicitly stated otherwise) refer to the naive hard
scattering, before initial and final state radiation effects have been
included.
COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
Purpose: to provide information on latest event generated or, in a few
cases, on the accumulated statistics during the run.
MSTI(1) : specifies the general type of subprocess that has occured,
according to the ISUB code given in section 2.2.
MSTI(2) : whenever MSTI(1) (together with MSTI(15) and MSTI(16)) are
not enough to specify the type of process uniquely, MSTI(2)
provides an ordering of the different possibilities. This is
particularly relevant for the different colour flow topologies
possible in QCD 2 -> 2 processes. With i = MSTI(15), j = MSTI(16)
and k = MSTI(2), the QCD possibilities are, in the classification
scheme of [Ben84]:
ISUB = 11, i = j, q_i q_i -> q_i q_i;
k = 1 : colour configuration A.
k = 2 : colour configuration B.
ISUB = 11, i /= j, q_i q_j -> q_i q_j;
k = 1 : only possibility.
ISUB = 12, q_i q~_i -> q_l q~_l;
k = 1 : only possibility.
ISUB = 13, q_i q~_i -> g g;
k = 1 : colour configuration A.
k = 2 : colour configuration B.
ISUB = 28, q_i g -> q_i g;
k = 1 : colour configuration A.
k = 2 : colour configuration B.
ISUB = 53, g g -> q_l q~_l;
k = 1 : colour configuration A.
k = 2 : colour configuration B.
ISUB = 68, g g -> g g;
k = 1 : colour configuration A.
k = 2 : colour configuration B.
k = 3 : colour configuration C.
ISUB = 83, f q -> f' Q (by t-channel W exchange; does not
distingusih colour flows but result of user selection);
k = 1 : heavy flavour Q is produced on side 1.
k = 2 : heavy flavour Q is produced on side 2.
MSTI(3) : number of partons produced in the hard interactions, i.e.
the number n of the 2 -> n matrix elements used; is sometimes 3 or
4 when a basic 2 -> 1 or 2 -> 2 process has been folded with
two 1 -> 2 initial branchings (like q q' -> q" q'" H0).
MSTI(4) : number of documentation lines in beginning of common block
LUJETS that are given with K(I,1) = 21; 0 for MSTP(125) = 0.
MSTI(5) : number of events generated to date in current run.
MSTI(6) : current frame of event, cf. MSTP(124).
MSTI(7), MSTI(8) : line number for documentation of outgoing partons/
particles from hard scattering for 2 -> 2 or 2 -> 1 -> 2 processes
(else = 0).
MSTI(10) : is 1 if cross-section maximum was violated in current
event, and 0 if not.
MSTI(11) : KF flavour code for beam (side 1) particle.
MSTI(12) : KF flavour code for target (side 2) particle.
MSTI(13), MSTI(14) : KF flavour codes for side 1 and side 2 initial
state shower initiators.
MSTI(15), MSTI(16) : KF flavour codes for side 1 and side 2 incoming
partons to the hard interaction.
MSTI(17), MSTI(18) : flag to signal if particle on side 1 or side 2
has been scattered diffractively; 0 if no, 1 if yes.
MSTI(21) - MSTI(24) : KF flavour codes for outgoing partons from the
hard interaction. The number of positions actually used is
process-dependent, see MSTI(3); trailing positions not used are
set = 0.
MSTI(25), MSTI(26) : KF flavour codes of the products in the decay
of a single s-channel resonance formed in the hard interaction.
Are thus only used when MSTI(3) = 1 and the resonance is allowed
to decay.
MSTI(31) : number of hard or semihard scatterings that occured in
current event in the multiple interaction scenario; is = 0 for a
low-p_T event.
MSTI(41) : the number of pileup events generated in latest PYEVNT
call (including the first, 'hard' event).
MSTI(42) - MSTI(50) : ISUB codes for the events 2 - 10 generated in
the pileup events scenario. The first event ISUB code is stored
in MSTI(1). If MSTI(41) is less than 10, only as many positions are
filled as there are pileup events. If MSTI(41) is above 10, some
ISUB codes will not appear anywhere.
PARI(1) : total integrated cross-section for the processes under study,
in mb. This number is obtained as a by-product of the selection of
hard process kinematics, and is thus known with better accuracy
when more event have been generated. The value stored here is based
on all events up till the latest one generated.
PARI(2) : is the ratio PARI(1)/MSTI(5), i.e. the ratio of total
integrated cross-section and number of events generated.
Histograms generated with unit weight for events have to be
multiplied by this factor, at the end of the run, to convert
results to mb.
For MSTP(142)=1, MSTI(5) is replaced by the sum of PARI(10)
values. Histograms are then filled with weight PARI(10) for
each event and multiplied by PARI(2) at the end.
PARI(9) : is weight WTXS returned from PYEVWT call when MSTP(142) >= 1,
else is 1.
PARI(10) : is compensating weight 1./WTXS that should be associated
to events when MSTP(142) = 1, else is 1.
PARI(11) : E_CM, i.e. total CM energy.
PARI(12) : s, i.e. total CM energy-squared.
PARI(13) : m-hat = sqrt(s-hat), i.e. mass of the hard scattering
subsystem.
PARI(14) : s-hat of the hard subprocess (2 -> 2 or 2 -> 1).
PARI(15) : t-hat of the hard subprocess (2 -> 2 or 2 -> 1 -> 2).
PARI(16) : u-hat of the hard subprocess (2 -> 2 or 2 -> 1 -> 2).
PARI(17) : p_T-hat of the hard subprocess (2 -> 2 or 2 -> 1 -> 2),
evaluated in the rest frame of the hard interaction.
PARI(18) : p_T-hat^2 of the hard subprocess; see PARI(17).
PARI(19) : m'-hat, the mass of the complete three- or four-body final
state in 2 -> 3 or 2 -> 4 processes (while m-hat, given in PARI(13),
here corresponds to the one- or two-body central system).
Kinematically m-hat <= m'-hat <= E_CM.
PARI(20) : s'-hat = m'-hat^2; see PARI(19).
PARI(21) : Q of the hard subprocess. The exact definition is
process-dependent, see MSTP(32).
PARI(22) : Q^2 of the hard subprocess; see PARI(21).
PARI(23) : Q of the outer hard scattering subprocess. Agrees with
PARI(21) for a 2 -> 1 or 2 -> 2 process. For a 2 -> 3 or 2 -> 4
W/Z fusion process, it is set by the W/Z mass scale, and for
subprocesses 121 and 122 by the heavy quark mass.
PARI(24) : Q^2 of the outer hard scattering subprocess; see PARI(23).
PARI(25) : Q scale used as maximum virtuality in parton showers.
Is equal to PARI(23), except for Deep Inelastic Scattering
processes when MSTP(22) > 0.
PARI(26) : Q^2 scale in parton showers; see PARI(25).
PARI(31), PARI(32) : the momentum fractions x of the initial
state parton shower initiators on side 1 and 2, respectively.
PARI(33), PARI(34) : the momentum fractions x taken by the partons
at the hard interaction, as used e.g. in the structure functions.
PARI(35) : Feynman-x, x_F = x_1 - x_2 = PARI(33) - PARI(34).
PARI(36) : tau = s-hat/s = x_1*x_2 = PARI(33)*PARI(34).
PARI(37) : y* = (1/2) * log(x_1/x_2), i.e. rapidity of the hard
interaction subsystem in the CM frame of the event-as-a-whole.
PARI(38) : tau' = s'-hat/s = PARI(20)/PARI(12).
PARI(39), PARI(40) : the primordial k_T values selected in the two
beam remnants.
PARI(41) : cos(theta-hat), where theta-hat is the scattering angle of
a 2 -> 2 (or 2 -> 1 -> 2) interaction, defined in the rest frame
of the hard scattering subsystem.
PARI(42) : x_T, i.e. scaled transverse momentum of the hard
scattering subprocess, x_T = 2 * p_T-hat/E_CM.
PARI(43), PARI(44) : x_L3 and x_L4, i.e. longitudinal momentum
fractions of the two scattered partons, in the range -1 < x_L < 1,
in the CM frame of the event-as-a-whole.
PARI(45), PARI(46) : x_3 and x_4, i.e. scaled energy fractions of the
two scattered partons, in the CM frame of the event-as-a-whole.
PARI(47), PARI(48) : y*_3 and y*_4, i.e. rapidities of the two
scattered partons in the CM frame of the event-as-a-whole.
PARI(49), PARI(50) : eta*_3 and eta*_4, i.e. pseudorapidities of the
two scattered partons in the CM frame of the event-as-a-whole.
PARI(51), PARI(52) : cos(theta*_3) and cos(theta*_4), i.e. cosines of
the polar angles of the two scattered partons in the CM frame of
the event-as-a-whole.
PARI(53), PARI(54) : theta*_3 and theta*_4, i.e. polar angles of the
two scattered partons, defined in the range 0 < theta* < pi, in
the CM frame of the event-as-a-whole.
PARI(55), PARI(56) : azimuthal angles phi*_3 and phi*_4 of the two
scattered partons, defined in the range -pi < phi* < pi, in the
CM frame of the event-as-a-whole.
PARI(61) : multiple interaction enhancement factor for current event.
A large value corresponds to a central collision and a small value
to a peripheral one.
PARI(65) : sum of the transverse momenta of partons generated at
the hardest interaction of the event, excluding initial and
final state radiation, i.e. 2 * PARI(17).
PARI(66) : sum of the transverse momenta of all partons generated at
the hardest interaction, including initial and final state
radiation, resonance decay products, and primordial k_T.
PARI(67) : sum of transverse momenta of partons generated at hard
interactions, excluding the hardest one (see PARI(65)), and also
excluding initial and final state radiation. Is nonvanishing only
in the multiple interaction scenario.
PARI(68) : sum of transverse momenta of all partons generated at hard
interactions, excluding the hardest one (see PARI(66)), but
including initial and final state radiation. Is nonvanishing only
in the multiple interaction scenario.
PARI(69) : sum of transverse momenta of all partons generated in hard
interactions (PARI(66) + PARI(68)) and, additionally, of all beam
remnant partons.
PARI(71), PARI(72) : sum of the momentum fractions x taken by initial
state parton shower initiators on side 1 and and side 2, excluding
those of the hardest interaction. Is nonvanishing only in the
multiple interaction scenario.
PARI(73), PARI(74) : sum of the momentum fractions x taken by the
partons at the hard interaction on side 1 and side 2, excluding
those of the hardest interaction. Is nonvanishing only in the
multiple interaction scenario.
PARI(75), PARI(76) : the x value of a photon that branches into
quarks or gluons, i.e. x at interface between initial state QED
and QCD cascades.
PARI(77), PARI(78) : the chi values selected for beam remnants that
are split into two objects, describing how the energy is shared
(see MSTP(92) etc.); is 0. if no splitting is needed.
PARI(81) : size of the threshold factor (enhancement or suppression)
in the latest event with heavy flavour production; see MSTP(35).
PARI(91) : average multiplicity nbar of pileup events, defined as
(cross-section for processes allowed by MSTP(132)) * PARP(131).
Only relevant for MSTP(133) = 1 or 2.
PARI(92) : average multiplicity of pileup events as actually
simulated, i.e. with multiplicity = 0 events removed and the
high-end tail truncated. Only relevant for MSTP(133) = 1 or 2.
PARI(93) : for MSTP(133) = 1 it is the probability that a beam crossing
will produce a pileup event at all, i.e. that there will be at
least one hadron-hadron interaction; for MSTP(133) = 2 the
probability that a beam crossing will produce a pileup event with
one hadron-hadron interaction of the desired rare type.
______________________________________________________________________
2.7. The Event Record
When an event is generated by a PYEVNT call, it is stored in the
common block LUJETS. Here each parton and particle is represented by
one line of information, giving a status code, parton/particle code,
line of mother, lines of first and last daughter (or colour
flow information in hard interaction and parton showers),
momentum, energy, mass, production vertex and time, and invariant
lifetime (if unstable). The LUJETS common block is used extensively
both by the PYTHIA and the JETSET routines; indeed it provides the
bridge that allows the general utility routines in JETSET to be used
also for PYTHIA events. A detailed description of LUJETS is found in
the JETSET 7.3 manual, a companion file to the program. The PYTHIA
event listing begins (optionally) with a few lines of event summary,
specific to the hard process simulated and thus not described in the
JETSET manual. Therefore we here give a brief summary on the general
structure of the LUJETS common block, followed by details specific to
PYTHIA.
COMMON/LUJETS/N,K(4000,5),P(4000,5),V(4000,5)
N : number of entries, i.e. number of lines used in the K, P and V
matrices, in the current event.
K(I,1) : status code KS for the parton/particle stored in the line.
= 0 : empty line.
= 1 - 10 : undecayed particle or unfragmented parton.
= 11 - 20 : decayed particle or fragmented or showered parton.
= 21 : documentation lines for compressed event summary (see below).
K(I,2) : KF flavour code for partons and particles, following the 1988
Particle Data Group numbering conventions [PDG88]. Some of the most
frequently used ones are
1 d 11 e- 21 g 31
2 u 12 nu_e 22 gamma 32 Z'0
3 s 13 mu- 23 Z0 33
4 c 14 nu_mu 24 W+ 34 W'+
5 b 15 tau- 25 H0 35 H'0
6 t 16 nu_tau 26 36 A0
7 l 17 chi- 27 37 H+
8 h 18 nu_chi 28 38
9 19 29 39 LQ
10 20 30 40 R0
1103 dd_1 2103 ud_1 3101 sd_0 3101 su_0
2101 ud_0 2203 uu_1 3103 sd_1 3103 su_1
211 pi+ 213 rho+ 2112 n 1114 Delta-
311 K0 313 K*0 2212 p 2114 Delta0
321 K+ 323 K*+ 3122 Lambda0 2214 Delta+
411 D+ 413 D*+ 3112 Sigma- 2224 Delta++
421 D0 423 D*0 3212 Sigma0 3114 Sigma*-
431 D_s+ 433 D*_s+ 3222 Sigma+ 3214 Sigma*0
130 K_L0 3312 Xi- 3224 Sigma*+
310 K_S0 3322 Xi0 3334 Omega-
A negative KF code, where existing, always corresponds to the
antiparticle of the one listed above.
To designate diffractive states, three non-standard codes are made
available by JETSET.
210 pi_diffr+ 2110 n_diffr 2210 p_diffr
K(I,3) : line number of parent parton or jet, where known, else 0.
K(I,4) : normally line number of the first daughter; for K(I,1) = 3, 13
or 14 instead special colour flow information.
K(I,5) : normally line number of the last daughter; for K(I,1) = 3, 13
or 14 instead special colour flow information.
P(I,1) : p_x, momentum in the x direction (in GeV/c).
P(I,2) : p_y, momentum in the y direction (in GeV/c).
P(I,3) : p_z, momentum in the z direction (in GeV/c).
P(I,4) : E, energy (in GeV).
P(I,5) : m, mass (in GeV/c^2). Spacelike partons are
given with P(I,5) = -sqrt(-m^2).
V(I,1) : x position of production vertex (in mm).
V(I,2) : y position of production vertex (in mm).
V(I,3) : z position of production vertex (in mm).
V(I,4) : time of production (in mm/c = 3.33*10^-12 s).
V(I,5) : proper lifetime of the particle/parton (in mm/c = 3.33*10^-12
s); is 0 for a stable entry.
In most instances, only the actual partons and particles produced are of
interest. For MSTP(125) = 0, the event record starts off with the parton
configuration existing after hard interaction, initial and final state
radiation, multiple interactions and beam remnants have been considered.
The partons are arranged in colour singlet clusters, ordered as required
for string fragmentation. Also photons and leptons produced as part of
the hard interaction (e.g. from q q~ -> g gamma or u u~ -> Z0 -> e+ e-)
appear in this part of the event record. These original entries
appear with pointer K(I,3) = 0, whereas the products of the subsequent
fragmentation and decay have K(I,3) numbers pointing back to the line
of the parent.
The standard documentation, obtained with MSTP(125) = 1, includes a few
lines in the beginning of the event record, which contain a brief
summary of the process that has taken place. The number of lines used
depends on the nature of the hard process, and is stored in MSTI(4) for
the current event. These lines all have K(I,1) = 21. For all
processes, lines 1 and 2 give the two incoming hadrons. When listed with
LULIST, these two lines will be separated from subsequent lines by a
sequence of ====== signs, to improve readability. For diffractive and
elastic events, the two outgoing states in lines 3 and 4 completes the
list. Otherwise, lines 3 and 4 contain the two partons that initiate
the two initial state parton showers, and 5 and 6 the endproducts of
these showers, i.e. the partons that enter the hard interaction. With
initial state radiation switched off, lines 3 and 5 and lines 4 and 6
coincide. For a simple 2 -> 2 hard scattering, lines 7 and 8 give the
two outgoing partons/particles from the hard interaction, before any
final state radiation. For 2 -> 2 processes proceeding via an
intermediate resonance such as Z0/gamma*, W+/-, H0 or R, the resonance
is found in 7 and the two outgoing partons/particles in 8 and 9. In some
cases one or both of these may be a resonance in its own right, so that
further pairs of lines are added for subsequent decays. If the decay of
a given resonance has been switched off, then no decay products are
listed either in this initial summary or in the subsequent ordinary
listing. Whenever partons are listed, they are assumed on mass shell for
simplicity. The fact that effective masses may be generated by initial
and final state radiation is taken into account in the actual parton
configuration that is allowed to fragment, however. A special case is
provided by W+ W- or Z0 Z0 fusion to a H0. Then the virtual W:s or Z:s
are shown in lines 7 and 8, the H0 in line 9 and the two recoiling
quarks (that emitted the bosons) in 10 and 11, followed by the Higgs
decay products. Since the W:s and Z:s are spacelike, what is actually
listed as the mass for them is -sqrt(-m^2). The listing of the event
documentation closes with another line made up of ====== signs.
A few examples may help clarify the picture. For a single diffractive
event p + p~ -> p_diffr + p~, the event record will start with
I K(I,1) K(I,2) K(I,3) comment
1 21 2212 0 incoming p
2 21 -2212 0 incoming p~
========================= not part of record; appears in listings
3 21 27 1 outgoing p_diffr
4 21 -2212 2 outgoing p~
=========================
The typical QCD 2 -> 2 process would be
I K(I,1) K(I,2) K(I,3) comment
1 21 2212 0 incoming p
2 21 -2212 0 incoming p~
=========================
3 21 2 1 u picked from incoming p
4 21 -1 2 d~ picked from incoming p~
5 21 21 3 u evolved to g at hard scattering
6 21 -1 4 still d~ at hard scattering
7 21 21 0 outgoing g from hard scattering
8 21 -1 0 outgoing d~ from hard scattering
=========================
Note that, where well defined, the K(I,3) code does contain information
on which side the different partons come from, e.g. above the gluon in
line 5 points back to the u in line 3, which points back to the proton
in line 1. In the example above it would have been possible to
associate the scattered g in line 7 with the incoming one in line 5, but
this is not possible in the general case, consider e.g. g g -> g g.
As a final example, W+ W- fusion to a H0 might look like
I K(I,1) K(I,2) K(I,3) comment
1 21 2212 0 first incoming p
2 21 2212 0 second incoming p
=========================
3 21 2 1 u picked from first p
4 21 21 2 g picked from second p
5 21 2 3 still u after initial state radiation
6 21 -4 4 g evolved to c~
7 21 24 5 spacelike W+ emitted by u quark
8 21 -24 6 spacelike W- emitted by c~ quark
9 21 25 0 Higgs produced by W+ W- fusion
10 21 1 5 u turned into d by emission of W+
11 21 -3 6 c~ turned into s~ by emission of W-
12 21 23 9 first Z0 coming from decay of H0
13 21 23 9 second Z0 coming from decay of H0
14 21 12 12 electron neutrino from first Z0 decay
15 21 -12 12 electron antineutrino from first Z0 decay
16 21 5 13 b quark from second Z0 decay
17 21 -5 13 b~ quark from second Z0 decay
=========================
After these lines with initial information, the event record looks the
same as for MSTP(125) = 0, i.e. first comes the parton configuration
to be fragmented and, after another separator line ====== in the output
(but not the event record), the products of subsequent fragmentation
and decay chains. The K(I,3) pointers for the partons, as well as
leptons and photons produced in the hard interaction, are now pointing
towards the documentation lines above, however. In particular, beam
remnants point to 1 or 2, depending on which side they belong
to, and partons emitted in the initial state parton showers point to 3
or 4. In the second example above, the partons produced by final state
radiation will be pointing back to 7 and 8; as usual, it should be
remembered that a specific assignment to 7 or 8 need not be unique. For
the third example, final state radiation partons will come both from
partons 10 and 11 and from partons 16 and 17, and additionally there
will be a neutrino-antineutrino pair pointing to 14 and 15. The extra
pairs of partons that are generated by multiple interactions do not
point back to anything, i.e. they have K(I,3) = 0.
There exists a third documentation option, MSTP(125) = 2. Here the
history of initial and final state parton branchings may be traced,
including all details on colour flow. This information has not been
optimized for user-friendliness, and can not be recommended for
general usage. With this option, the initial documentation lines
are the same. They are followed by blank lines, K(I,1) = 0, up to
line 20 (can be changed in MSTP(126)). From line 21 and onwards
each parton with K(I,1) = 3, 13 or 14 appear with special colour flow
information in the K(I,4) and K(I,5) positions; see the JETSET 7.3
manual. For an ordinary 2 -> 2 scattering, the two incoming partons
at the hard scattering are stored in line 21 and 22, and the two
outgoing in 23 and 24. The colour flow between these partons has
to be chosen according to the proper relative probabilities in
cases when many alternatives are possible, see [Ben84]. If there
is initial state radiation, the two partons in line 21 and 22 are
copied down to line 25 and 26, from which the initial state showers are
reconstructed backwards step by step. The branching history may be read
by noting that, for a branching a -> b + c, the K(I,3) codes of b and c
point towards the line number of a. Since the showers are reconstructed
backwards, this actually means that parton b would appear in the listing
before parton a and c, and hence have a pointer to a position below
itself in the list. Associated timelike partons c may initiate timelike
showers, as may the partons of the hard scattering. Again a showering
parton or pair of partons will be copied down towards the end of the
list and allowed to undergo successive branchings c -> d + e, with d and
e pointing towards c. The mass of timelike partons is properly stored in
P(I,5), for spacelike partons instead -sqrt(-m~2) is stored. After this
section containing all the branchings comes the final parton
configuration, properly arranged in colour, followed by all
subsequent fragmentation and decay products, as usual.
______________________________________________________________________
2.8. Other Routines and Commonblocks
The subroutines and commonblocks that a user will come in direct
contact with have already been described. A number of other routines
and commonblocks exist, and are here briefly listed for the sake of
completeness.
SUBROUTINE PYINKI
Purpose: to initialize the kinematics given by the two incoming
hadrons/leptons.
SUBROUTINE PYINRE
Purpose: to initialize the widths and effective widths of resonances.
SUBROUTINE PYXTOT
Purpose: to give the parametrized total, double diffractive, single
diffractive and elastic cross-sections for different energies and
colliding hadrons.
SUBROUTINE PYMAXI
Purpose: to find optimal coefficients COEF for the selection of
kinematical variables, and to find the related maxima for
the differential cross-section times Jacobian factors,
for each of the subprocesses included.
SUBROUTINE PYPILE
Purpose: to determine the number of pileup events, i.e. events
appearing in the same beam-beam crossing.
SUBROUTINE PYRAND
Purpose: to generate the quantities characterizing a hard scattering
on the parton level, according to the relevant matrix elements.
SUBROUTINE PYSCAT
Purpose: to find outgoing flavours and to set up the kinematics and
colour flow of the hard scattering.
SUBROUTINE PYSSPA(IP1,IP2)
Purpose: to generate the spacelike showers of the initial state
radiation.
SUBROUTINE PYRESD
Purpose: to allow Z0, W+/-, H0, Z'0, H+/- and R0 resonances to decay,
including chains of successive decays and parton showers.
SUBROUTINE PYMULT
Purpose: to generate semihard interactions according to the
multiple interaction formalism.
SUBROUTINE PYREMN(IPU1,IPU2)
Purpose: to add on target remnants and include primordial k_T.
SUBROUTINE PYDIFF
Purpose: to handle diffractive and elastic scattering events.
SUBROUTINE PYDOCU
Purpose: to compute cross-sections of processes, based on current
Monte Carlo statistics, and to store event information in the
MSTI and PARI arrays.
SUBROUTINE PYWIDT
Purpose: to calculate widths and effective widths of resonances.
SUBROUTINE PYKLIM
Purpose: to calculate allowed kinematical limits.
SUBROUTINE PYKMAP
Purpose: to calculate the value of a kinematical variable when this
is selected according to one of the simple pieces.
Purpose: to give the differential cross-section (multiplied by the
relevant Jacobians) for a given subprocess and kinematical
setup.
SUBROUTINE PYSTFU(KF,X,Q2,XPQ)
SUBROUTINE PYSTEL(X,Q2,XPEL)
SUBROUTINE PYSTGA(X,Q2,XPGA)
SUBROUTINE PYSTPI(X,Q2,XPPI)
SUBROUTINE PYSTPR(X,Q2,XPPR)
Purpose: to give gluon, quark and antiquark structure functions for
given x and Q^2 values for e+-, gamma, p, p~, n, n~, pi+-, and
hyperons.
SUBROUTINE PYSPLI
Purpose: to give hadron remnant or remnants left when the reacting
parton is kicked out.
FUNCTION PYGAMM(X)
Purpose: to give the value of the ordinary gamma function Gamma(x)
(used in some structure function parametrizations).
SUBROUTINE PYWAUX(IAUX,EPS,WRE,WIM)
SUBROUTINE PYI3AU(EPS,RAT,Y3RE,Y3IM)
Purpose: to calculate the value of some auxiliary functions appearing
in the cross-section expressions in PYSIGH.
FUNCTION PYSPEN(XREIN,XIMIN,IREIM)
Purpose: to calculate the real and imaginary part of the Spence
function.
BLOCK DATA PYDATA
Purpose: to give sensible default values to all status codes and
parameters.
COMMON/PYINT1/MINT(400),VINT(400)
Purpose: to collect a host of integer and real valued variables used
internally in the program during the initialization and/or event
generation stage. These variables must not be changed by the user.
MINT(1) : specifies the general type of subprocess that has occured,
according to the ISUB code given in section 2.2.
MINT(2) : whenever MINT(1) (together with MINT(15) and MINT(16)) are
not enough to specify the type of process uniquely, MINT(2)
provides an ordering of the different possibilities, see
MSTI(2).
MINT(3) : number of partons produced in the hard interactions, i.e.
the number n of the 2 -> n matrix elements used; is sometimes 3 or
4 when a basic 2 -> 1 or 2 -> 2 process has been folded with
two 1 -> 2 initial branchings (like q q' -> q" q'" H0).
MINT(4) : number of documentation lines in beginning of common block
LUJETS that are given with K(I,1) = 21; 0 for MSTP(125) = 0.
MINT(5) : number of events generated to date in current run.
MINT(6) : current frame of event, cf. MSTP(124).
MINT(7), MINT(8) : line number for documentation of outgoing partons/
particles from hard scattering for 2 -> 2 or 2 -> 1 -> 2 processes
(else = 0).
MINT(10) : is 1 if cross-section maximum was violated in current
event, and 0 if not.
MINT(11) : KF flavour code for beam (side 1) particle.
MINT(12) : KF flavour code for target (side 2) particle.
MINT(13), MINT(14) : KF flavour codes for side 1 and side 2 initial
state shower initiators.
MINT(15), MINT(16) : KF flavour codes for side 1 and side 2 incoming
partons to the hard interaction.
MINT(17), MINT(18) : flag to signal if particle on side 1 or side 2
has been scattered diffractively; 0 if no, 1 if yes.
MINT(19), MINT(20) : flag to signal initial state structure with
parton inside photon inside electron on side 1 or side 2;
0 if no, 1 if yes.
MINT(21) - MINT(24) : KF flavour codes for outgoing partons from the
hard interaction. The number of positions actually used is
process-dependent, see MINT(3); trailing positions not used are
set = 0.
MINT(25), MINT(26) : KF flavour codes of the products in the decay
of a single s-channel resonance formed in the hard interaction.
Are thus only used when MINT(3) = 1 and the resonance is allowed
to decay.
MINT(31) : number of hard or semihard scatterings that occured in
current event in the multiple interaction scenario; is = 0 for a
low-p_T event.
MINT(35) : in a true 2 -> 3 process, where one particle is a resonance
with decay channel selected already before the PYRESD call, the
decay channel number (in the /LUDAT3/ numbering) is stored here.
MINT(41), MINT(42) : type of incoming beam or target particle;
1 for lepton and 2 for hadron. A photon counts as a lepton if it
is not resolved (MSTP(14)=0) and as a hadron if it is resolved
(MSTP(14)=1).
MINT(43) : combination of incoming beam and target particles. A photon
counts as a hadron.
= 1 : lepton on lepton.
= 2 : lepton on hadron.
= 3 : hadron on lepton.
= 4 : hadron on hadron.
MINT(44) : as MINT(43), but a photon counts as a lepton.
MINT(45), MINT(46) : structure of incoming beam and target particles.
= 1 : no internal structure, i.e. electron or photon which carries
full beam energy.
= 2 : defined with structure functions which are not peaked at
x = 1, i.e. hadrons and a resolved photon.
= 3 : defined with structure functions which are peaked at x = 1,
i.e. the electron.
MINT(47) : combination of incoming beam and target particle structure
function types.
= 1 : no structure function either for beam or target.
= 2 : structure functions for target but not for beam.
= 3 : structure functions for beam but not for target.
= 4 : structure functions both for beam and target, but not
both peaked at x = 1.
= 5 : structure functions both for beam and target, with both
peaked at x = 1.
MINT(48) : total number of subprocesses switched on.
MINT(49) : number of subprocesses that are switched on, apart from
elastic scattering and single, double and central diffractive.
MINT(51) : internal flag that event failed cuts.
= 0 : no problem.
= 1 : event failed; new one to be generated.
MINT(52) : internal counter for number of lines used (in /LUJETS/)
before multiple interactions are considered.
MINT(53) : internal counter for number of lines used (in /LUJETS/)
before beam remnants are considered.
MINT(55) : the heaviest new flavour switched on for QCD processes,
specifically the flavour to be generated for ISUB = 81, 82 or
83 or 84.
MINT(56) : the heaviest new flavour switched on for QED processes,
specifically for ISUB = 85. Note that, unlike MINT(55), the
heaviest flavour may here be a lepton, and that heavy means the
one with largest KF code
MINT(61) : internal switch for the mode of operation of resonance
width calculations in PYWIDT for gamma*/Z0 or gamma*/Z0/Z'0.
= 0 : without reference to initial state flavours.
= 1 : with reference to given initial state flavours.
= 2 : for given final state flavours.
MINT(62) : internal switch for use at initialization of H width.
= 0 : use widths into Z + Z* or W + W* calculated before.
= 1 : evaluate widths into Z + Z* or W + W* for current Higgs mass.
MINT(65) : internal switch to indicate initialization without
specified reaction.
= 0 : normal initialization.
= 1 : initialization with argument 'none' in PYINIT call.
MINT(71) : switch whether current process is singular for p_T -> 0
or not.
= 0 : non-singular process, i.e. proceeding via an s-channel
resonance or with both products having a mass above CKIN(6).
= 1 : singular process.
MINT(72) : number of s-channel resonances which may contribute to
the cross-section.
MINT(73) : KF code of first s-channel resonance; 0 if there is none.
MINT(74) : KF code of second s-channel resonance; 0 if there is none.
MINT(81) : number of pileup events selected.
MINT(82) : sequence number of currently considered pileup event.
MINT(83) : number of lines in the event record already filled by
previously considered pileup events.
MINT(84) : MINT(83) + MSTP(126), i.e. number of lines already filled
by previously considered events plus number of lines to be kept
free for event documentation.
VINT(1) : ECM, CM energy.
VINT(2) : s (=ECM^2) mass-square of complete system.
VINT(3) : mass of beam particle.
VINT(4) : mass of target particle.
VINT(5) : momentum of beam (and target) particle in CM frame.
VINT(6) - VINT(10) : theta, phi and beta for rotation and boost
from CM frame to user-specified frame.
VINT(11) : tau_min.
VINT(12) : y*_min.
VINT(13) : cos(theta-hat)_min for cos(theta-hat) <= 0.
VINT(14) : cos(theta-hat)_min for cos(theta-hat) >= 0.
VINT(15) : x_T^2_min.
VINT(16) : tau'_min.
VINT(21) : tau.
VINT(22) : y*.
VINT(23) : cos(theta-hat).
VINT(24) : phi* (azimuthal angle).
VINT(25) : x_T^2.
VINT(26) : tau'.
VINT(31) : tau_max.
VINT(32) : y*_max.
VINT(33) : cos(theta-hat)_max for cos(theta-hat) <= 0.
VINT(34) : cos(theta-hat)_max for cos(theta-hat) >= 0.
VINT(35) : x_T^2_max.
VINT(36) : tau'_max.
VINT(41), VINT(42) : the momentum fractions x taken by the partons
at the hard interaction, as used e.g. in the structure functions.
VINT(43) : m-hat = sqrt(s-hat), mass of hard scattering subsystem.
VINT(44) : s-hat of the hard subprocess (2 -> 2 or 2 -> 1).
VINT(45) : t-hat of the hard subprocess (2 -> 2 or 2 -> 1 -> 2).
VINT(46) : u-hat of the hard subprocess (2 -> 2 or 2 -> 1 -> 2).
VINT(47) : p_T-hat of the hard subprocess (2 -> 2 or 2 -> 1 -> 2),
i.e. transverse momentum evaluated in the rest frame of the
scattering.
VINT(48) : p_T-hat^2 of the hard subprocess; see VINT(47).
VINT(49) : m'-hat, the mass of the complete three- or four-body final
state in 2 -> 3 or 2 -> 4 processes.
VINT(50) : s'-hat = m'-hat^2; see VINT(49).
VINT(51) : Q of the hard subprocess. The exact definition is
process-dependent, see MSTP(32).
VINT(52) : Q^2 of the hard subprocess; see VINT(51).
VINT(53) : Q of the outer hard scattering subprocess. Agrees with
VINT(51) for a 2 -> 1 or 2 -> 2 process. For a 2 -> 3 or 2 -> 4
W/Z fusion process, it is set by the W/Z mass scale, and for
subprocesses 121 and 122 by the heavy quark mass.
VINT(54) : Q^2 of the outer hard scattering subprocess; see VINT(53).
VINT(55) : Q scale used as maximum virtuality in parton showers.
Is equal to VINT(53), except for Deep Inelastic Scattering
processes when MSTP(22) > 0.
VINT(56) : Q^2 scale in parton showers; see VINT(55).
VINT(57) : alpha_em value of hard process.
VINT(58) : alpha_strong value of hard process.
VINT(59) : sin(theta-hat) (cf. VINT(23)); used for improved numerical
precision in elastic and diffractive scattering.
VINT(61), VINT(62) : nominal m^2 values, i.e. without initial state
radiation effects, for the two partons entering the hard
interaction.
VINT(63), VINT(64) : nominal m^2 values, i.e. without final state
radiation effects, for the two (or one) partons/particles leaving
the hard interaction.
VINT(65) : p-hat_init, i.e. common nominal absolute momentum of the two
partons entering the hard interaction, in their rest frame.
VINT(66) : p-hat_fin, i.e. common nominal absolute momentum of the two
partons leaving the hard interaction, in their rest frame.
VINT(71) : p_T_min of process, i.e. CKIN(3) or CKIN(5), depending on
which is larger, and whether the process is singular in p_T -> 0
or not.
VINT(73) : tau = m^2/s value of first resonance, if any; see MINT(73).
VINT(74) : m*Gamma/s value of first resonance, if any; see MINT(73).
VINT(75) : tau = m^2/s value of second resonance, if any; see MINT(74).
VINT(76) : m*Gamma/s value of second resonance, if any; see MINT(74).
VINT(80) : correction factor (evaluated in PYOFSH) for the cross-section
of resonances produced in 2 -> 2 processes, if only some mass range
of the full Breit-Wigner shape is allowed by user-set mass cuts
(CKIN(2), CKIN(45) - CKIN(48)).
VINT(81) - VINT(84) : the cos(theta) and phi variables of a true 2 -> 3
process, where one product is a resonance, effectively giving
2 -> 4. The first two are cos(theta) and phi for the resonance
decay, the other two ditto for the effective system formed by the
other two particles.
VINT(85), VINT(86) : transverse momenta in a true 2 -> 3 process; one
is stored in VINT(47) (that of the Z in g g -> Z Q Q), while the
smaller of the two others is stored in VINT(85) and the larger
in VINT(86).
VINT(91), VINT(92) : gives a dimensionless suppression factor, to take
into account reduction in cross-section due to the allowed channels
for a W+W+ or W-W- pair, respectively, in the same sense as
WIDS(24,1) gives it for a W+W- pair.
VINT(98) : is sum of VINT(100) values for current run.
VINT(99) : is weight WTXS returned from PYEVWT call when MSTP(142) >= 1,
else is 1.
VINT(100) : is compensating weight 1./WTXS that should be associated
to events when MSTP(142) = 1, else is 1.
VINT(101) : total cross-section.
VINT(102) : elastic cross-section.
VINT(103) : single diffractive cross-section.
VINT(104) : double diffractive cross-section.
VINT(105) : central diffractive cross-section.
VINT(106) : total non-diffractive, inelastic cross-section.
VINT(108) : ratio of maximum differential cross-section observed
to maximum differential cross-section assumed in generation;
cf. MSTP(123).
VINT(109) : ratio of minimal (negative!) cross-section observed
to maximum differential cross-section assumed in generation;
could only become negative if cross-sections are incorrectly
included.
VINT(111) - VINT(116) : for MINT(61) = 1 gives kinematical factors
for the different pieces contributing to gamma*/Z0 or gamma*/Z0/Z'0
production, for MINT(61) = 2 gives sum of final state weights
for the same; coefficients are given in the order pure gamma*,
gamma*-Z0 interference, gamma*-Z'0 interference, pure Z0, Z0-Z'0
interference and pure Z'0.
VINT(117) : width of Z0; needed in gamma*/Z0/Z'0 production.
VINT(121) : nuclear slope parameter B in t-hat distribution for
(diffractive and) elastic scattering.
VINT(122) : curvature parameter C in t-hat distribution for (diffractive
and) elastic scattering.
VINT(131) : total cross-section (in mb) for subprocesses allowed in
the pileup events scenario according to the MSTP(132) value.
VINT(132) : nbar = VINT(131) * PARP(131), cf. PARI(91).
VINT(133) : = (sum_i Prob_i * i)/(sum_i Prob_i) as actually
simulated, i.e. 1 <= i <= 100 (or smaller), cf. PARI(92).
VINT(134) : is exp(-nbar) * sum_i nbar^i/i! for MSTP(133) = 1 and
exp(-nbar) * sum_i nbar^i/(i-1)! for MSTP(133) = 2, cf. PARI(93).
VINT(138) : size of the threshold factor (enhancement or suppression)
in the latest event with heavy flavour production; see MSTP(35).
VINT(141), VINT(142) : x values for the parton shower initiators of the
hardest interaction; used to find what is left for multiple
interactions.
VINT(143), VINT(144) : 1 - (sum of x-values) for all scatterings; used
for rescaling each new x-value in the multiple interaction structure
function evaluation.
VINT(145) : estimate of total parton-parton cross-section for multiple
interactions; used for MSTP(82) >= 2.
VINT(146) : common correction factor f_c in the multiple interaction
probability; used for MSTP(82) >= 2.
VINT(147) : average hadronic matter overlap; used for MSTP(82) >= 2.
VINT(148) : enhancement factor for current event in multiple interaction
probability, defined as the actual overlap divided by the average
one; used for MSTP(82) >= 2.
VINT(149) : x_T^2 cutoff or turnoff for multiple interactions.
For MSTP(82) <= 1 it is 4*p_Tmin^2/s, for MSTP(82) >= 2 it is
4*p_T0^2/s.
VINT(150) : probability to keep given event in multiple interaction
scenario, as given by the "Sudakov" form factor.
VINT(151), VINT(152) : sum of x values for all multiple interaction
partons.
VINT(153) : current differential cross-section value obtained from
PYSIGH; used in multiple interactions only.
VINT(155), VINT(156) : the x value of a photon that branches into
quarks or gluons, i.e. x at interface between initial state QED
and QCD cascades.
VINT(157), VINT(158) : the primordial k_T values selected in the two
beam remnants.
VINT(159), VINT(160) : the chi values selected for beam remnants that
are split into two objects, describing how the energy is shared
(see MSTP(92) etc.); is 0. if no splitting is needed.
VINT(161) - VINT(200) : sum of Cabibbo-Kobayashi-Maskawa matrix
element-squared that a given flavour is allowed to couple to.
Results are stored in format VINT(180+KF) for quark and lepton
flavours and antiflavours (which need not be the same; see
MDME(IDC,2). For leptons, these factors are normally unity.
VINT(201) - VINT(220) : additional variables needed in phase space
selection for 2 -> 3 processes with ISET(ISUB) = 5. Below indices
1, 2 and 3 refer to scattered partons 1, 2 and 3, except that
q variables are q_1 + q_2 -> q_1' q_2' q_3' for four-momenta.
All kinematical variables refer to the internal kinematics of
the 3-body final state - the kinematics of the system as a whole
is described by tau' and y*, and the mass distribution of
particle 3 (a resonance) by tau.
VINT(201) : m_1.
VINT(202) : p_T1^2.
VINT(203) : phi_1.
VINT(204) : M_1 (mass of propagator).
VINT(205) : weight for the p_T1^2 choice.
VINT(206) : m_2.
VINT(207) : p_T2^2.
VINT(208) : phi_2.
VINT(209) : M_2 (mass of propagator).
VINT(210) : weight for the p_T2^2 choice.
VINT(211) : y_3.
VINT(212) : y_3max.
VINT(213) : epsilon = +-1; choice between two mirror solutions 1 <-> 2.
VINT(214) : weight associated to epsilon-choice.
VINT(215) : t_1 = (q_1 - q_1')^2.
VINT(216) : t_2 = (q_2 - q_2')^2.
VINT(217) : q_1*q_2'.
VINT(218) : q_2*q_1'.
VINT(219) : q_1'*q_2'.
VINT(220) : sqrt((m_T12^2-m_T1^2-m_T2^2)^2-4.*m_T1^2*m_T2^2).
COMMON/PYINT2/ISET(200),KFPR(200,2),COEF(200,20),ICOL(40,4,2)
Purpose: to store information necessary for efficient generation of the
different subprocesses, specifically type of generation scheme and
coefficients of the Jacobian. Also to store cumulative
cross-sections needed for multiple interaction generation for
MSTP(82) >= 2. These variables must not be changed by the user.
ISET(ISUB) : gives the type of kinematical variable selection scheme
used for subprocess ISUB.
= 0 : elastic, diffractive and low-p_T processes.
= 1 : 2 -> 1 processes (irrespective of subsequent decays).
= 2 : 2 -> 2 processes (i.e. the bulk of processes).
= 3 : 2 -> 3 processes (like q q' -> q" q'" H0).
= 4 : 2 -> 4 processes (like q q' -> q" q'" W+ W-).
= 5 : 'true' 2 -> 3 processes, one method.
= 6 : 'true' 2 -> 3 processes, another method;
currently only g + g -> Z + Q + Q~.
= 9 : 2 -> 2 in multiple interactions (pT as kinematics variable).
= -1 : legitimate process which has not yet been implemented.
= -2 : ISUB is an undefined process code.
KFPR(ISUB,J) : give the KF flavour codes for the products produced in
subprocess ISUB. If there is only one product, the J = 2 position
is left blank. Also, quarks and leptons assumed massless in the
matrix elements are denoted by 0. The main application is thus to
identify resonances produced (Z0, W+, H0, etc.).
COEF(ISUB,J) : factors used in the Jacobians in order to speed up the
selection of kinematical variables. More precisely, the shape of
the cross-section is given as the sum of terms with different
behaviour, where the integral over the allowed phase space is
unity for each term. COEF gives the relative strength of these
terms, normalized so that the sum of coefficients for each variable
used is unity. Note that which coefficients are indeed used is
process-dependent.
ISUB : standard subprocess code.
J = 1 : tau selected according 1/tau.
J = 2 : tau selected according to 1/tau^2.
J = 3 : tau selected according to 1/(tau*(tau+tau_r)), where
tau_r = m^2/s is tau value of resonance; only used for
resonance production.
J = 4 : tau selected according to Breit-Wigner of form
1/((tau-tau_r)^2+gam_r^2), where tau_r = m^2/s is tau value
of resonance and gam_r = m*Gamma/s is its scaled mass times
width; only used for resonance production.
J = 5 : tau selected according to 1/(tau*(tau+tau'_r)), where
tau'_r = m'^2/s is tau value of second resonance; only
used for simultaneous production of two resonances.
J = 6 : tau selected according to second Breit-Wigner of form
1/((tau-tau'_r)^2+gam'_r^2), where tau'_r = m'^2/s is tau value
of second resonance and gam'_r = m'*Gamma'/s is its scaled
mass times width; is used only for simultaneous production
of two resonances, like gamma*/Z0/Z'0.
J = 7 : tau selected according to 1/(1-tau); only used when both
structure functions are peaked at x = 1.
J = 8 : y* selected according to y* - y*min;
J = 9 : y* selected according to y*max - y*;
J = 10 : y* selected according to 1/cosh(y*).
J = 11 : y* selected according to 1/(1-exp(y*-y*max)); only used
when beam structure function is peaked close to x = 1.
J = 12 : y* selected according to 1/(1-exp(y*min-y*)); only used
when target structure function is peaked close to x = 1.
J = 13 : z = cos(theta-hat) selected evenly between limits.
J = 14 : z = cos(theta-hat) selected according to 1/(a-z),
where a = 1 + 2*m_3^2*m_4^2/shat^2, m_3 and m_4 being the
masses of the two final state particles.
J = 15 : z = cos(theta-hat) selected according to 1/(a+z),
with a as above.
J = 16 : z = cos(theta-hat) selected according to 1/(a-z)^2,
with a as above.
J = 17 : z = cos(theta-hat) selected according to 1/(a+z)^2,
with a as above.
J = 18 : tau' selected according to 1/tau'.
J = 19 : tau' selected according to (1 - tau/tau')^3/tau'^2.
J = 20 : tau' selected according to 1/(1-tau'); only used when
both structure functions are peaked close to x = 1.
ICOL : contains information on different colour flow topologies in
hard 2 -> 2 processes.
COMMON/PYINT3/XSFX(2,-40:40),ISIG(1000,3),SIGH(1000)
Purpose: to store information on structure functions, subprocess cross-
sections and different final state relative weights. These variables
must not be changed by the user.
XSFX : current values of structure functions (multiplied by x) on
beam and target side.
ISIG(ICHN,1) : incoming parton/particle on the beam side to the hard
interaction for allowed channel no. ICHN. The number of channels
filled with relevant information is given by NCHN, one of the
arguments returned in a PYSIGH call. Thus only 1 <= ICHN <= NCHN
is filled with relevant information.
ISIG(ICHN,2) : incoming parton/particle on the target side to the hard
interaction for allowed channel no. ICHN. See also comment above.
ISIG(ICHN,3) : colour flow type for allowed channel no. ICHN; see
MSTI(2) list. See also comment above. For 'subprocess' 96 uniquely,
ISIG(ICHN,3) is also used to translate information on what is the
correct subprocess number (11, 12, 13, 28, 53 or 68); this is
used for reassigning subprocess 96 to either of these.
SIGH(ICHN) : evaluated differential cross-section for allowed channel
no. ICHN, i.e. hard matrix element value times structure
functions, for current kinematical setup (in addition, Jacobian
factors are included in the figures, as used to speed up
generation). See also comment for ISIG(ICHN,1).
COMMON/PYINT4/WIDP(21:40,0:40),WIDE(21:40,0:40),WIDS(21:40,3)
Purpose: to store partial and effective decay widths for the different
resonances. These variables must not be changed by the user.
WIDP(KF,J) : gives partial decay widths of resonances into different
channels (in GeV), given that all physically allowed final
states are included.
KF : standard KF code for resonance considered.
J : enumerates the different decay channels possible for
resonance KF, as stored in the JETSET LUDAT3 commonblock,
with the first channel in J = 1, etc.
WIDE(KF,J) : gives effective decay widths of resonances inte different
channels (in GeV), given the decay modes actually left open in
the current run. The on/off status of decay modes is set by the
MDME switches in JETSET; see section 2.9.
KF : standard KF code for resonance considered.
J : enumerates the different decay channels possible for
resonance KF, as stored in the JETSET LUDAT3 commonblock,
with the first channel in J = 1, etc.
WIDS(KF,J) : gives a dimensionless suppression factor, which is defined
as the ratio of the total width of channels switched on to the total
width of all possible channels (replace width by width-squared for
a pair of resonances). The on/off status of channels is set by the
MDME switches in JETSET; see section 2.9. The information in WIDS
is used e.g. in cross-section calculations.
KF : standard KF code for resonance considered.
J = 1 : suppression when a pair of resonances of type KF are
produced together. When an antiparticle exists, the
particle-antiparticle pair (like W+ W-) is the relevant
combination, else the particle-particle one (like Z0 Z0).
J = 2 : suppression for a particle of type KF when produced
on its own, or together with a particle of another type.
J = 3 : suppression for an antiparticle of type KF when produced
on its own, or together with a particle of another type.
COMMON/PYINT5/NGEN(0:200,3),XSEC(0:200,3)
Purpose: to store information necessary for cross-section calculation
and differential cross-section maximum violation. These variables
must not be changed by the user.
NGEN(ISUB,1) : gives the number of times that the differential
cross-section (times Jacobian factors) has been evaluates for
subprocess ISUB, with NGEN(0,1) the sum of these.
NGEN(ISUB,2) : gives the number of times that a kinematical setup
for subproces ISUB is accepted in the generation procedure,
with NGEN(0,2) the sum of these.
NGEN(ISUB,3) : gives the number of times an event of subprocess type
ISUB is generated, with NGEN(0,3) the sum of these. Usually
NGEN(ISUB,3) = NGEN(ISUB,2), i.e. an accepted kinematical
configuration normally can be used to produce an event.
XSEC(ISUB,1) : estimated maximum differential cross-section (times the
Jacobian factors used to speed up the generation process) for
the different subprocesses in use, with XSEC(0,1) the sum of these
(except low-p_T, i.e. ISUB = 95).
XSEC(ISUB,2) : gives the sum of differential cross-sections (times
Jacobian factors) for the NGEN(ISUB,1) phase space points
evaluated so far.
XSEC(ISUB,3) : gives the estimated integrated cross-section for
subprocess ISUB, based on the statistics accumulated so far,
with XSEC(0,3) the estimated total cross-section for all
subprocesses included (all in mb). This is exactly the
information obtainable by a PYSTAT(1) call.
COMMON/PYINT6/PROC(0:200)
Purpose: to store character strings for the different possible
subprocesses; used when printing tables.
PROC(ISUB) : name for the different subprocesses, according to ISUB
code. PROC(0) denotes all processes.
_________________
Finally, in addition a number of routines and commonblocks with names
beginning with RK come with the program. These contain the matrix
element evaluation for the process g + g -> Z + q + q~, based on a
program of Ronald Kleiss, with only minor modifications.
______________________________________________________________________
2.9. The JETSET routines
Since the fragmentation and decay of an original parton configuration is
handled by the JETSET routines, using the LUJETS common block described
above, all of the non-e+e- routines and parameters are at the disposal
of the PYTHIA user. This section is intended as a brief reminder of some
of the most useful ones, in particular for hadron physics applications.
Further details may be obtained by consulting the JETSET manual.
The subroutine LULIST can be used to list an event; for output that fits
onto a normal terminal screen CALL LULIST(1) is recommended, whereas
CALL LULIST(2) gives a slightly more comprehensive listing. The event
record nominally contains not only the stable final particles, but also
the original partons and intermediate resonances. With CALL LUEDIT(1)
all partons/particles that have fragmented/decayed are removed from the
event record, and N is updated accordingly. A CALL LUEDIT(2) will remove
neutrinos as well, and CALL LUEDIT(3) will only leave charged, stable
particles. Some of the information indirectly stored in LUJETS can be
accessed more easily by using the PLU and KLU functions, thus e.g.
y = PLU(I,17) gives the true rapidity and eta = PLU(I,19) the
pseudorapidity of particle I.
With MSTP(111) = 0, fragmentation and decays are switched off in
the PYEVNT call. If one desires to fragment the event at some later
opportunity, this can be achieved by a CALL LUEXEC. The action of
the program during this call, which is done automatically from
PYEVNT for MSTP(111) = 1, can be controlled by the switches and
parameters in the LUDAT1 common block:
COMMON/LUDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
Thus, MSTJ(21) = 0 will switch off particle decays but keep jet
fragmentation on, MSTJ(1) can be set to give independent
fragmentation rather than the standard string fragmentation, etc.
Contrary to most PYTHIA variables, the user can freely change the
values during the course of the run, since no initialization is
involved.
The same commonblock also contains a number of conversion factors and
coupling constants. Thus MSTU(101) can be set to used fix or running
alpha_em, PARU(101) to set the alpha_em value, PARU(102) to set
sin^2(theta_W). A host of couplings for particles beyond the
standard model may be found in PARU(121) onwards. These couplings can
not be changed in midrun, since they are used in the initialization to
estimate cross-section maxima.
The common block LUDAT2 contains several pieces of information on
particle and flavour properties:
COMMON/LUDAT2/KCHG(500,3),PMAS(500,4),PARF(2000),VCKM(4,4)
The particle masses may be the most interesting variables here.
Since these masses are directly taken over by the PYTHIA routines,
any mass changes should be made before the PYINIT call. Of
particular interest are the following:
PMAS(6,1) : (D=120. GeV/c^2) mass of top quark.
PMAS(7,1) : (D=120. GeV/c^2) mass of l quark.
PMAS(8,1) : (D=200. GeV/c^2) mass of h quark.
PMAS(17,1) : (D=100. GeV/c^2) mass of chi lepton.
PMAS(23,1) : (D=91.2 GeV/c^2) Z0 nominal mass.
PMAS(24,1) : (D=80.0 GeV/c^2) W+/- nominal mass.
PMAS(25,1) : (D=50. GeV/c^2) H0 nominal mass.
PMAS(32,1) : (D=500. GeV/c^2) Z'0 nominal mass.
PMAS(34,1) : (D=500. GeV/c^2) W'+/- nominal mass.
PMAS(35,1) : (D=300. GeV/c^2) H'0 nominal mass.
PMAS(36,1) : (D=300. GeV/c^2) A0 nominal mass.
PMAS(37,1) : (D=300. GeV/c^2) H+/- nominal mass.
PMAS(39,1) : (D=200. GeV/c^2) LQ nominal mass.
PMAS(40,1) : (D=5000. GeV/c^2) R nominal mass.
For the resonances above (particles 23 - 40), the resonance width is
calculated at PYINIT initialization, taking into account the channels
assumed possible for resonance decays (for H0 this will depend
very much on the mass assumed), so the nominal widths stored in
PMAS(KF,2) are actually overwritten at this point.
Note that the first index in the KCHG and PMAS arrays is not the
standard KF code described in section 2.7, which may be both positive
and negative, and extends to very large values. Instead a compressed
code KC is used, which only runs between 1 and 500, i.e. which contains
no distinction between particles and antiparticles. The KF and KC
codes agree between 1 and 40, i.e. for quarks, leptons and gauge bosons,
while meson and baryon codes have been compressed. Thus, several
hadrons containing heavy quarks (some c and b baryons, most t, l
and h hadrons) have been lumped together, with masses and charges
instead derived from the quark content (and, for the mass, the spin).
The user is not supposed to learn the KC code, but instead use the
LUCOMP function to translate from KF to KC code, i.e. KC = LUCOMP(KF).
PMAS(LUCOMP(663),1) = 110. changes the toponium mass to 110 GeV, to
give one example. LUCOMP will return KC = 0 for any particle code
unknown in the program.
The common block LUDAT3 contains all particle decay data:
COMMON/LUDAT3/MDCY(500,3),MDME(2000,2),BRAT(2000),KFDP(2000,5)
It can be used to switch off the decay of some particles selectively,
and even to switch on and off separate decay modes for a given
particle. These switches are used extensively within PYTHIA, both
to select allowed final state flavours and to calculate the reduction
in cross-section due to these restrictions. MDCY(KC,1) is an on/off
switch (= 0 off, = 1 on) for the decay of compressed particle code KC.
The decay of Z0, both in PYTHIA and in JETSET, is thus switched off
by putting MDCY(23,1) = 0. In order to keep the event record from
becoming uncomfortably large, switching off pi0 decays by putting
MDCY(LUCOMP(111),1) = 0 is probably the single most effective action.
All the decay channels defined in the program are stored sequentially
in MDME, BRAT and KFDP. The entry point into this table and the number
of channels stored for a given particle is given by the last two
entries of the MDCY array. In order to obtain a list of particle data,
including all decay channels, do CALL LULIST(12) (and take out a
hardcopy for future reference). For resonances only, the same
information is available with CALL PYSTAT(2). For a given decay channel
IDC, the switch in MDME(IDC,1) is of special interest. It can be used to
switch on or off that particular channel, or leave it selectively open.
Effective branching ratios are automatically recalculated among the
decay channels left open. If a particle is allowed to decay by the
MDCY(KC,1) value, at least one channel must be left open by the user.
Since the MDME(IDC,1) code is of special interest in PYTHIA for
selecting flavours in single gauge boson or gauge boson pair decays,
the JETSET manual description on this point is reproduced in extenso:
= -1 : this is a non-standard model decay mode, which by default is
assumed not to exist. It is therefore not simulated; neither does
it contribute to the total width of the particle. Normally, this
option is used for decays involving fourth generation or H+-
particles.
= 0 : channel is switched off. It is not simulated, but still
contributes to the total width of the particle.
= 1 : channel is switched on.
= 2 : channel is switched on for a particle but off for an
antiparticle. It is also on for a particle its own antiparticle,
i.e. here it means the same as =1.
= 3 : channel is switched on for an antiparticle but off for a
particle. It is off for a particle its own antiparticle.
= 4 : in the production of a pair of equal or charge conjugate
resonances in PYTHIA, say H0 -> W+ W-, either one of the
resonances is allowed to decay according to this group of
channels, but not both. If the two particles of the pair
are different, the channel is on.
Within JETSET, this option only means that the channel is
switched off.
= 5 : as =4, but an independent group of channels, such that in
a pair of equal or charge conjugate resonances the decay of
either resonance may be specified independently. If the two
particles in the pair are different, the channel is off.
Within JETSET, this option only means that the channel is
switched off.
WARNING: the two values -1 and 0 may look similar, but in fact
are quite different. In neither case the channel so set is
generated, but in the latter case the channel still contributes
to the total width of a resonance, and thus affects both
simulated line shape and the generated cross-section. Thus the
value 0 is appropriate to a channel we assume exists, even if
we are not currently simulating it, while -1 should be used for
channels we believe do not exist. In particular, users are warned
unwittingly to set fourth generation channels 0 (rather than -1),
since by now the support for a fourth generation is small.
All the options above may be freely mixed. The difference,
for those cases where both make sense, between using values
2 and 3 and using 4 and 5 is that the latter automatically
include charge conjugate states, e.g. H0 -> W+ W- ->
e+ nue d ubar or dbar u e- nuebar, but the former only one
of them. In calculations of the joint branching ratio, this
makes a factor 2 difference.
One should note that, whereas the MDME code can be used to set
allowed decay channels both for resonances (Z0, W+/-, H0, Z'0,
W'+/-, H+/-, and R) and for ordinary particles (pi0, K, D, B,
top mesons, etc.), it is only for the resonances that cross-sections
are modified to take into account the allowed channels. The reason
is obvious: resonances are only produced as part of the hard
interaction, or as decay products of other resonances, according
to fairly strict rules, whereas ordinary particles can be produced
anywhere in the event (charm and bottom quark pairs may appear
inside initial and final state showers, e.g.) in a rather haphazard
fashion. It is therefore not possible to generate unbiased events
in a consistent manner for ordinary particles. (One might imagine
doing things properly for the t quark, since it is not likely to be
produced in showers, but this has not been done so far.)
For QCD processes, like the pair production of a new heavy quark
(ISUB = 81 and 82), imagined 'decay modes' are stored for the
gluon among the ordinary ones. Branching ratios here have no meaning,
but the MDME(IDC,1) switch can be used to enable or disable the
production of different new quarks. In the current version of the
massive matrix element treatment, only one heavy flavour can be
generated in a given one; if the user has left several flavours
switched on, only the heaviest is actually used. The same goes for
ISUB = 84, where the photon 'decay modes' can be selected. Similarly,
the 'decay modes' stored for quarks indicate which branchings of the
type q -> q' W are allowed, as used e.g. in W+ W- -> H0 production.
It is also possible to include completely new decay modes, or even a
completely new particle. This is described in the JETSET manual.
The event analysis routines can be called on to evaluate the event
stored in the LUJETS common block. Particularly useful is the cluster
finding routine LUCELL, which assumes a (modifiable) grid of cells in
pseudorapidity and azimuthal angle and carries out a cluster finding
algorithm a la the UA1 one. With a CALL LUCELL(NJET) the number of jets
found is given by NJET and the position and ET of those jets is stored
after the event proper, in lines N+1 through N+NJET.
The routine LUGIVE can be used to feed in commonblock values both to
JETSET and to PYTHIA. Since LUGIVE checks on array boundaries and
also writes on output the variables which have been changes, this may
at times be a better way to set up the values to be used by the program.
A call might look like
CALL LUGIVE('MSEL=0;MSUB(142)=1;PMAS(34,1)=300.;CKIN(1)=200.')
i.e. several commands can be put in one if separated by semicolon
(but it is of course also possible to do one command at a time and
call LUGIVE as often as desired). Note that for people who want data
card driven programs, this facilitates the construction of a driver
routine for PYTHIA.
______________________________________________________________________
2.10. On cross-sections
This section includes two related topics. Thew first is how events are
generated according to the proper matrix elements, the second how
the final cross-sections (available e.g. with PYSTAT) are calculated.
In PYRAND, the variables used in the generation of hard scattering
2 -> n processes are tau (= s-hat/s, i.e. mass-squared of
scattering subsystem scaled by total CM energy-squared), y*
(rapidity of scattering subsystem in the CM frame of the event),
z (= cos(theta-hat), i.e. cosine of scattering angle in CM frame
of the hard interaction; this applies to 2 -> 2 and 2 -> 4
processes only), and tau' (= s'/s, i.e. mass-squared of complete
three- or four-body final state scaled by total CM energy-squared;
this applies to 2 -> 3 and 2 -> 4 processes only). These variables
are generated according to the distributions h1(tau)/tau, h2(y*),
h3(z), and h4(tau')/tau', where
h1(tau) = c0 + I0/I1*c1*1/tau + I0/I2*c2*1/(tau + tau_R) +
+ I0/I3*c3*tau/((s*tau - m_R^2)^2 + (m_R*Gamma_R)^2) +
+ I0/I4*c4*1/(tau + tau_R') + I0/I5*c5 *
* tau/((s*tau - m_R'^2)^2 + (m_R'*Gamma_R')^2)
h2(y*) = I0/I1*c1*(y* - y*_min) + I0/I1*c2*(y*_max - y*) +
+ I0/I3*c3*1/cosh(y*)
h3(z) = c0 + I0/I1*c1*1/(a - z) + I0/I2*c1*1/(a + z) +
+ I0/I3*c3*1/(a - z)^2 + I0/I4*c4*1/(a + z)^2
h4(tau') = c0 + I0/I1*c1*(1 - tau/tau')^3/tau',
where subscripts R and R' in h1 denote values (of tau, mass and
width) pertaining to a resonance R or R', and a in h3 denotes
1 + 2*(m_3*m_4/s-hat)^2 (= 1 for massless products).
In each of the above cases, In denotes the integral over the
quantity multiplying coefficient cn; thus, e.g.,
h1: I0 = Int dtau/tau = ln(tau_max/tau_min)
h3: I0 = Int dz = z-_max - z-_min + z+_max - z+_min
h4: I0 = Int dtau'/tau' = ln(tau'_max/tau'_min)
(for h2, where no c0 occurs, I0 is defined to be Int dy* =
= y*_max - y*_min), and the coefficients cn are normalized, such
that Sum(n) cn = 1. (The symbols Int and Sum are used to denote
integral and sum over, respectively.)
For other types of processes, e.g. elastic and diffractive
scattering, other variables are used in the generation, but the
structure of the generation procedure is similar.
After this primary generation, the variables are then weighted
with the relevant matrix elements, etc., according to standard
Monte Carlo techniques.
The process-dependent weights are coded in the routine PYSIGH, as
follows. The cross-sections for 2 -> 1 and 2 -> 2 processes can be
written, respectively,
sigma(2 -> 1) = Int dx_1 Int dx_2 Sum(ijk) f_i f_j *
* dsigma-hat(ijk)
and
sigma(2 -> 2) = Int dx_1 Int dx_2 Int dt-hat Sum(ijk) f_i f_j *
* dsigma-hat(ijk)/dt-hat,
where Int and Sum denote respectively integral and sum over, f_i
the structure function for parton i, and dsigma-hat the
differential cross-section for the hard scattering. Converting to
the variables tau, y*, and z= cos(theta-hat), the cross-sections
can then be rewritten as
sigma(2 -> 1) = pi/s * Int dtau h1(tau)/tau Int dy* h2(y*) *
* 1/(tau*h1(tau)) * 1/h2(y*) * Sum(ijk) F_i F_j *
* s-hat/pi*dsigma-hat(ijk) =
= pi/s * I(tau) * I(y*) * 1/(tau*h1) * 1/h2 *
* Sum(ijk) F_i F_j * s-hat/pi*dsigma-hat(ijk)
and
sigma(2 -> 2) = pi/s * Int dtau h1(tau)/tau Int dy* h2(y*) *
* Int dz h3(z) * 1/(tau*h1(tau)) * 1/h2(y*) *
* 1/2*beta34/h3(z) * Sum(ijk) F_i F_j *
* s-hat^2/pi*dsigma-hat(ijk)/dt-hat =
= pi/s * I(tau) * I(y*) * I(z) * 1/(tau*h1) *
* 1/h2 * 1/2*beta34/h3 * Sum(ijk) F_1 F_j *
* s-hat^2/pi*dsigma-hat(ijk)/dt-hat,
where F_i denotes x*f_i, beta34 = lambda(1,r_3,r_4) with lambda
the usual lambda-function and r_i = m_i^2/s-hat, and
I(tau) = Int dtau h1(tau)/tau
I(y*) = Int dy* h2(y*)
I(z) = Int dz h3(z),
and the h's are distributions used in the actual generation of
the kinematical variables tau, y* and z (see above), such that
Int dtau h1(tau)/tau = ln(tau_max/tau_min)
Int dy* h2(y*) = y*_max - y*_min
Int dz h3(z) = z-_max - z-_min + z+_max - z+_min.
For each subprocess in PYSIGH, what is coded under the relevant
ISUB is thus the dimensionless quantity
s-hat/pi*dsigma-hat(ijk) (2 -> 1 processes)
s-hat^2/pi*dsigma-hat(ijk)/dt-hat (2 -> 2 processes)
which is then multiplied by the structure functions (or, rather,
x*f) and a common factor, containing pi/s, the conversion factor
from GeV^-2 to mb, I(tau), I(y*) (and, for 2 -> 2 processes,
I(z)) and 1/(tau*h1)), and 1/h2 (and, for 2 -> 2 processes, also
1/2*beta34/h3). It is the latter product which is returned in a
PYSIGH call; we will denote this SIGH(tau,y*,z) in the following.
For 2 -> 2 processes with identical final state particles, the
symmetrization factor of 1/2 is explicitly included at the very
end, when the result is stored in the SIGH array. In the final
cross-section, a factor of 2 is retrieved because of integration
over the full phase space (rather than only half of it).
For 2 -> 1 processes with dsigma-hat given in the zero-width
approximation, the delta function in tau has been replaced by
the (modified) Breit-Wigner
1/pi*s*H_R/((s*tau - m_R^2)^2 - H_R^2),
where
H_R = s*tau/m_R*Gamma_R,
with m_R and Gamma_R the mass and the width of the resonance R,
respectively. The Breit-Wigner thus contains an shat-dependent
width which improves the description of the line shape as compared
to the case with widths calculated at nominal resonance mass. In
the 2 -> 2 processes 71 - 77, where the Higgs occurs in s-, t- and
u-channel exchanges, the width is calculated at nominal Higgs mass,
however.
For effective 2 -> 3, 4 processes (e.g. Z0 Z0 -> H0), the common
factor further includes the integration over the extra variable
tau', i.e. a factor I(tau') * f(tau')/h4(tau'), where
f(tau') = (1 + tau/tau') ln(tau'/tau) - 2 (1 - tau/tau')
and
I(tau') = Int dtau' h4(tau')/tau' = ln(tau'_max/tau'_min).
A user wanting to include new processes (experts only!), will
have to insert then the function s-hat/pi*dsigma-hat(ijk) or
s-hat^2/pi*dsigma-hat(ijk)/dt-hat under the appropriate ISUB
heading (in addition to changes in e.g. SUBROUTINE PYSCAT).
Let us now describe how these differential cross-sections are used
to derive the total cross-sections (within the user-defined cuts)
for the allowed processes. We will first consider high-p_T event
generation only, and later cover low-p_T physics.
At initialization of a 2 -> 2 process, the coefficients cn appearing
in the expressions h1(tau), h2(y*) and h3(z) are first optimized
so as to fit the actual variation of the cross-section as well as can
be achieved, i.e. to see to it that SIGH(tau,y*,z) is a slowly varying
function of its arguments. Thereafter the maximum of SIGH(tau,y*,z)
is searched for. This gives a preliminary cross-section estimate,
which is printed among the initialization information, and which
represents an overestimate of the actual cross-section. It is this
maximum which is used below to find the relative weight for selecting
a processes. For 2 -> 1, 2-> 3 and 2 -> 4 processes, the same
philosophy is used, with trivial changes in the number of variables
that have to be considered.
In particular, note that SIGH is selected so as to allow a direct
comparison between the different classes of events. More precisely,
the SIGH value returned by PYSIGH is the total integrated cross-
section a process would have inside the simulated phase space volume,
if the actual variation of the cross-section were perfectly described
by h1(tau)/tau * h2(y*) * h3(z) (or equivalent), with overall
normalization given by the differential cross-section value in the
point selected.
In order to generate an event, first a process is chosen, thereafter
a set of kinematical variables for this process is selected, as
described above. The probability for accepting this set is given by
the actual quantity SIGH(tau,y*,z), divided by the maximum of this
quantity, as found earlier. In case of rejection, both process and
kinematical variables are to be chosen anew.
For each individual subprocess, there are counters giving the number
of points tried and the sum of the SIGH values in these points. At
the end of each PYEVNT call, the ratio of these two numbers is
calculated, thus providing a Monte Carlo estimate of the actual
cross-section for each process. If only very few events have been
generated, the statistical fluctuations are large. Therefore it is
recommended to call PYSTAT only after all events have been generated,
although a call in principle could be made at any time. In PYSTAT
the cross-sections will be listed together with the number of phase
space points tried and the number of events actually generated.
The cross-sections for each subprocess are calculated taking into
account which initial and final states are allowed, by the KFIN array
in the common block PYSUBS, and the MDME(IDC,1) values in the common
block LUDAT3, respectively. By default, everything up to and including
top is allowed, while it is assumed that no fourth generation exists.
If some initial states are switched off, the cross-sections will be
reduced accordingly. For several processes, like q + g -> q + g,
the flavours in the final state are the ones of the initial state,
and for these MDME has no effect. In processes like
g + g -> q + qb or q + qb -> Z0/gamma* -> mu+ + mu-, on the other hand,
the quark and lepton flavours allowed in the final state can be selected
using MDME, and cross-sections are only based on those channels that
are actually allowed. We emphasize that, for resonances like Z0 or W+/-,
the total widths are assumed given as the sum over three (or,
optionally, four) full generations of quarks and leptons. By only
allowing a few channels, we do thus in no way alter the basic properties
of the resonances, but only restrict the generation of events to those
particular channels we happen to be interested in, with the cross-
sections reduced accordingly. Particles like H0, Z' and W' can either
decay directly into a q-q~ or l-l~ pair, or into a (gauge) boson pair,
which then can decay further. Therefore the production cross-section
given by the program is reduced, not only by switching off some decay
channels of the primary resonance, but also by switching off subsequent
Z, W or H+ decay modes (and the branching ratios of the primary
resonance into the boson pair is reduced accordingly). Even if the
decay of a resonance is switched off completely, by using the MDCY
array in the LUDAT3 commonblock, cross-sections are still based on
the decay channels allowed by the MDME values in the same commonblock.
Note that what is said here about resonances does not hold for
ordinary particles or even the top quark, where a suppression of some
decay channels are never taken into account in PYTHIA cross-sections,
but have to be compensated for by the user himself.
The cross-section values printed by PYSTAT are stored in XSEC(ISUB,3)
in common block PYINT5. In particular, XSEC(0,3) = PARI(1) gives the
total cross-section for all subprocesses studied in current run. Events
are generated with unit weight. At the end of a run, histogram contents
etc. can be converted from number of events to differential cross-
sections by multiplying by XSEC(0,3) and dividing by the total number of
events, the latter stored in MSTI(5). This ratio is given in PARI(2).
For histograms one should also divide by the bin width.
What has been said so far applies when only hard interactions are
studied. The total hadronic cross-section, as well as the elastic,
single diffractive and double diffractive cross-sections are taken from
parametrizations available in the literature. For the total and elastic
cross-sections this means the fits of Block-Cahn [Blo85], for
diffractive events the ansatz of Goulianos [Gou83] is used.
The cross-section for non-diffractive low-pT events is then
obtained by subtraction. When the multiple interaction formalism is used
to generate these latter events, some of them will actually contain a
semihard interaction, and be listed in PYSTAT as such, with only the
events without any semihard interactions being classified as low-p_T
(ISUB = 95). This matter of classification does not remove the
constraint that the total cross-section is fixed, however. Neither can
the KFIN switches be set to suppress some of the flavours of the
initial hadrons (the inclusion of heavy flavours or not can still be
set using MSTP(54), however).
______________________________________________________________________
2.11. Examples
The program is built as a slave system, i.e. the user supplies the main
program, which calls on the PYTHIA and JETSET routines to perform
specific tasks and then resumes control.
A typical program for analysis of collider events at 540 GeV CM energy
with a minimum p_T of 10 GeV/c at the hard scattering (because of
initial state radiation, fragmentation effects, etc., the actual
p_T-cutoff will be smeared around this value) might look like
COMMON/LUJETS/N,K(4000,5),P(4000,5),V(4000,5)
COMMON/PYSUBS/MSEL,MSUB(200),KFIN(2,-40:40),CKIN(200)
COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
... ! set all common block variables that
... ! did not have desired default values
CKIN(3)=10. ! lower p_T cutoff
CALL PYINIT('CMS','p','pbar',540.) ! initialize
... ! initialize analysis statistics
DO 100 IEVENT=1,1000 ! loop over events
CALL PYEVNT ! generate event
IF(IEVENT.EQ.1) CALL LULIST(1) ! list first event
... ! insert desired analysis chain for
... ! each event
100 CONTINUE
CALL PYSTAT(1) ! print cross-sections
... ! user output
END
**********************************************************************
References
Alt89 G. Altarelli, B. Mele, M. Ruiz-Altaba, Z. Physik C45 (1989) 109
Bar90 T.L. Barklow, SLAC-PUB-5364 (1990)
Bau90 U. Baur, M. Spira, P. M. Zerwas, Phys. Rev. D42 (1990) 815
Ben84 H.-U. Bengtsson, Computer Phys. Comm. 31 (1984) 323
Ben84a H.-U. Bengtsson, G. Ingelman, LU TP 84-3, Ref.TH.3820-CERN
Ben85 H.-U. Bengtsson, G. Ingelman, Computer Phys. Comm. 34 (1985) 251
Ben85a H.-U. Bengtsson, W.-S. Hou, A. Soni, D.H. Stork,
Phys. Rev. Lett. 55 (1985) 2762
Ben87 H.-U. Bengtsson, T. Sjostrand, Computer Phys. Comm. 46 (1987) 43
Ber84 E. L. Berger, E. Braaten, R. D. Field, Nucl. Phys.
B239 (1984) 52
Blo85 M. M. Block, R. N. Cahn, Rev. Mod. Phys. 57 (1985) 563
M. M. Block, R. N. Cahn, in Physics Simulations at High Energy,
eds. V. Barger, T. Gottschalk, F. Halzen (World Scientific,
Singapore, 1987), p. 89
BLS87 M. C. Bento, C. H. Llewellyn Smith, Nucl. Phys. B289 (1987) 36
Cah84 R.N. Cahn, S. Dawson, Phys. Lett. 136B (1984) 196
R.N. Cahn, Nucl. Phys. B255 (1985) 341
G. Altarelli, B. Mele, F. Pitolli, Nucl. Phys. B287 (1987) 205
Cha85 M. Chanowitz, M.K. Gaillard, Nucl. Phys. B261 (1985) 379
Cha91 K. Charchula, DESY 91-093 (1991)
Cha91a K. Charchula, ZEUS note 91-73 (1991)
Che75 M.-S. Chen, P. Zerwas, Phys. Rev. D12 (1975) 187;
P. Zerwas, private communication (1991)
Coc90 D. Cocolicchio, F. Feruglio, G.L. Fogli, J. Terron,
CERN-TH.5909/90
F. Feruglio, private communication
Com77 B. L. Combridge, J. Kripfganz, J. Ranft, Phys. Lett. 70B (1977)
234
Con71 V. Constantini, B. de Tollis, G. Pistoni, Nuovo Cim. 2A (1971)
733
Cut78 R. Cutler, D. Sivers, Phys. Rev. D17 (1978) 196
DFL88 M. Diemoz, F. Ferroni, E. Longo, G. Martinelli, Z. Physik C39
(1988) 21
DHT90 A. Dobado, M. J. Herrero, J. Terron, CERN-TH.5670/90 and
FTUAM/08-90, CERN-TH.5813/90
Dic88 D. A. Dicus, S. S. D. Willenbrock, Phys. Rev. D37 (1988) 1801
Dre85 M. Drees, K. Grassie, Z. Physik C28 (1985) 451
Dre89 M. Drees, J. Ellis, D. Zeppenfeld, Phys. Lett. B223 (1989) 454
Duk82 D.W. Duke, J.F. Owens, Phys. Rev. D26 (1982) 1600
Dun86 M. J. Duncan, G. L. Kane, W. W. Repko, Nucl. Phys. B272 (1986)
517
EHL84 E. Eichten, I. Hinchliffe, K. Lane, C. Quigg, Rev. Mod. Phys.
56 (1984) 579; 58 (1985) 1065
Ell86 R. K. Ellis, J. C. Sexton, Nucl. Phys. B269 (1986) 445
Ell88 R. K. Ellis, I. Hinchliffe, M. Soldate, J. J. van der Bij,
Nucl. Phys. B297 (1988) 221
Fad89 V. Fadin, V. Khoze, JETP Lett. 46 (1987) 417, Yad. Fiz 48
(1988) 487, in Proceedings of the 24th Winter School of the
LNPI (Leningrad, 1989), vol. I, p. 3
V. Fadin, V, Khoze, T. Sjostrand, Z. Phys. C48 (1990) 613
Fon81 M. Fontannaz, B. Pire, D. Schiff, Z. Phys. C11 (1981) 211
Gab86 E. Gabrielli, Mod. Phys. Lett. A1 (1986) 465
GHK90 J.F. Gunion, H.E. Haber, G. Kane, S. Dawson, The Higgs Hunter's
Guide (Addison-Wesley, 1990)
A. Djouadi, private communication
Gou83 K. Goulianos, Phys. Rep. 101 (1983) 169
GRV89 M. Gluck, E. Reya, A. Vogt, DO-TH 89/20 (extended version)
Gun86 J.F. Gunion, Z. Kunszt, Phys. Rev. D33 (1986) 665,
errata as private communication from the authors
Hag90 K. Hagiwara, H. Iwasaki, A. Miyamoto, H. Murayama,
D. Zeppenfeld, Durham preprint DTP/90/82 etc. (1990)
Hal78 F. Halzen, D. M. Scott, Phys. Rev. D18 (1978) 3378
Hew88 J.L. Hewett, S. Pakvasa, Phys. Rev. D37 (1988) 3165
and private communication from the authors
Ing87 G. Ingelman et al., in Proceedings of the HERA Workshop,
ed. R.D. Peccei, Vol. 1, p. 3
Kat83 M. Katuya, Phys. Lett. 124B (1983) 421
Kle89 R. Kleiss et al., in Z Physics at LEP 1, eds. G. Altarelli,
R. Kleiss, C. Verzegnassi, CERN 89-08, Vol. 3, p.1
Kle90 R. Kleiss, private communication (1990)
Kun84 Z. Kunszt, Nucl. Phys. B247 (1984) 339
Kun87 Z. Kunszt, CERN 87-07, Vol. 1, p. 123, and private communication
Lan90 K. Lane, private communication (1991)
Mor90 J.G. Morfin, W.-K. Tung, Fermilab-Pub-90/74 and IIT-PHY-90/11
PDG88 Particle Data Group, G. P. Yost et al.,
Phys. Lett. 204B (1988) 1
Plo91 H. Plothow-Besch, in Proceedings of the 3rd Workshop on Detector
and Event Simulation in High Energy Physics, Amsterdam,
8-12 April 1991
Sam90 M.A. Samuel, G. Li, N. Sinha, R. Sinha, M.K. Sundaresan,
Carleton preprint (1990)
Sjo87 T. Sjostrand, M. van Zijl, Phys. Rev. D36 (1987) 2019
Sjo90 T. Sjostrand, JETSET 7.3 program and manual (unpublished,
obtainable on request)
Ter90 J. Terron, private communication (paper in preparation)
Tun87 W.-K. Tung, in Physics Simulations at High Energy,
eds. V. Barger, T. Gottschalk, F. Halzen (World Scientific,
Singapore, 1987), p. 601
W.-K. Tung, private communication
Wud86 J. Wudka, Phys. Lett. 167B (1986) 337
Zer90 P. Zerwas, private communication
**********************************************************************