The Main Subroutines

There are two routines that you 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 are frequently useful. The other two routines described here, PYFRAM and PYKCUT, are of more specialized interest.

Purpose:

to initialize the generation procedure. Normally it is foreseen that this call will be followed by many PYEVNT ones, to generate a sample of the event kind specified by the PYINIT call. (For problems with cross section estimates in runs of very few events per PYINIT call, see the description for PYSTAT(1) below in this subsection.)

FRAME :

a character variable used to specify the frame of the experiment. Upper-case and lower-case letters may be freely mixed.

= 'CMS' :

colliding beam experiment in c.m. 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.

= '3MOM' :

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 PYJETS. Particles are assumed on the mass shell, and energies are calculated accordingly.

= '4MOM' :

as '3MOM', except also energies should be specified, in P(1,4) and P(2,4), respectively. The particles need not be on the mass shell; effective masses are calculated from energy and momentum. (But note that numerical precision may suffer; if you know the masses the option '5MOM' below is preferable.)

= '5MOM' :

as '3MOM', except also energies and masses should be specified, i.e the full momentum information in P(1,1) - P(1,5) and P(2,1) - P(2,5) should be given for beam and target, respectively. Particles need not be on the mass shell. Space-like virtualities should be stored as $-\sqrt{-m^2}$. Especially useful for physics with virtual photons. (The virtuality could be varied from one event to the next, but then it is convenient to initialize for the lowest virtuality likely to be encountered.) Four-momentum and mass information must match.

= 'USER' :

a run primarily intended to involve external, user-defined processes, see subsection . Information on incoming beam particles and energies is read from the HEPRUP common block. In this option, the BEAM, TARGET and WIN arguments are dummy.

= '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 cannot 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. Upper-case and lower-case letters may be freely mixed. An antiparticle can be denoted by `bar' at the end of the name (`$\sim$' is a valid alternative for reasons of backwards compatibility). It is also possible to leave out the underscore (`_') directly after `nu' in neutrino names, and the charge for proton and neutron. The arguments are dummy when the FRAME argument above is either 'USER' or 'NONE'.

= 'e-' :

electron.

= 'e+' :

positron.

= 'nu_e' :

$\nu_{\mathrm{e}}$.

= 'nu_ebar' :

$\overline{\nu}_{\mathrm{e}}$.

= 'mu-' :

$\mu^-$.

= 'mu+' :

$\mu^+$.

= 'nu_mu' :

$\nu_{\mu}$.

= 'nu_mubar' :

$\overline{\nu}_{\mu}$.

= 'tau-' :

$\tau^-$.

= 'tau+' :

$\tau^+$.

= 'nu_tau' :

$\nu_{\tau}$.

= 'nu_taubar' :

$\overline{\nu}_{\tau}$.

= 'gamma' :

photon (real, i.e. on the mass shell).

= 'gamma/e-' :

photon generated by the virtual-photon flux in an electron beam; WIN below refers to electron, while photon energy and virtuality varies between events according to what is allowed by CKIN(61) - CKIN(78).

= 'gamma/e+' :

as above for a positron beam.

= 'gamma/mu-' :

as above for a $\mu^-$ beam.

= 'gamma/mu+' :

as above for a $\mu^+$ beam.

= 'gamma/tau-' :

as above for a $\tau^-$ beam.

= 'gamma/tau+' :

as above for a $\tau^+$ beam.

= 'pi0' :

$\pi^0$.

= 'pi+' :

$\pi^+$.

= 'pi-' :

$\pi^-$.

= 'n0' :

neutron.

= 'nbar0' :

antineutron.

= 'p+' :

proton.

= 'pbar-' :

antiproton.

= 'K+' :

$\mathrm{K}^+$ meson; since parton distributions for strange hadrons are not available, very simple and untrustworthy recipes are used for this and subsequent hadrons, see subsection .

= 'K-' :

$\mathrm{K}^-$ meson.

= 'KS0' :

$\mathrm{K}_{\mathrm{S}}^0$ meson.

= 'KL0' :

$\mathrm{K}_{\mathrm{L}}^0$ meson.

= 'Lambda0' :

$\Lambda$ baryon.

= 'Sigma-' :

$\Sigma^-$ baryon.

= 'Sigma0' :

$\Sigma^0$ baryon.

= 'Sigma+' :

$\Sigma^+$ baryon.

= 'Xi-' :

$\Xi^-$ baryon.

= 'Xi0' :

$\Xi^0$ baryon.

= 'Omega-' :

$\Omega^-$ baryon.

= 'pomeron' :

the pomeron $\mathrm{I}\!\mathrm{P}$; since pomeron parton distribution functions have not been defined this option can not be used currently.

= 'reggeon' :

the reggeon $\mathrm{I}\!\mathrm{R}$, with comments as for the pomeron above.

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='3MOM', '4MOM', '5MOM' :

dummy (information is taken from the P vectors, see above).

FRAME='USER' :

dummy (information is taken from the HEPRUP common block, see above).

FRAME='NONE' :

dummy (no information required).

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.)

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 desired information.

= 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 you in PYKCUT.
Note that no errors are given on the cross sections. In most cases a cross section is obtained by Monte Carlo integration during the course of the run. (Exceptions include e.g. total and elastic hadron-hadron cross sections, which are parameterized and thus known from the very onset.) A rule of thumb would then be that the statistical error of a given subprocess scales like $\delta \sigma / \sigma \approx 1/\sqrt{n}$, where $n$ is the number of events generated of this kind. In principle, the numerator of this relation could be decreased by making use of the full information accumulated during the run, i.e. also on the cross section in those phase space points that are eventually rejected. This is actually the way the cross section itself is calculated. However, once you introduce further cuts so that only some fraction of the generated events survive to the final analysis, you would be back to the simple $1/\sqrt{n}$ scaling rule for that number of surviving events. Statistical errors are therefore usually better evaluated within the context of a specific analysis. Furthermore, systematic errors often dominate over the statistical ones.
Also note that runs with very few events, in addition to having large errors, tend to have a bias towards overestimating the cross sections. In a typical case, the average cross section obtained with many runs of only one event each may be twice that of the correct answer of a single run with many events. The reason is a `quit while you are ahead' phenomenon, that an upwards fluctuation in the differential cross section in an early try gives an acceptable event and thus terminates the run, while a downwards one leads to rejection and a continuation of the run.

= 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 PYTHIA decay tables, the decay products (usually two but sometimes three), the partial decay width, branching ratio and effective branching ratio (in the event some channels have been excluded by you).

= 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 you 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.

= 6 :

prints a table of all subprocesses implemented in the program.

= 7 :

prints two tables related to $R$-violating supersymmetry, where lepton and/or baryon number is not conserved. The first is a collection of semi-inclusive branching ratios where the entries have a form like ~chi_10 -> nu + q + q, where a sum has been performed over all lepton and quark flavours. In the rightmost coloumn of the table, the number of modes that went into the sum is given. The purpose of this table is to give a quick overview of the branching fractions, since there are currently more than 1500 individual $R$-violating processes included in the generator. Note that only the pure $1\to 3$ parts of the 3-body modes are included in this sum. If a process can also proceed via two successive $1 \to 2$ branchings (i.e. the intermediate resonance is on shell) the product of these branchings should be added to the number given in this table. A small list at the bottom of the table shows the total number of $R$-violating processes in the generator, the number with non-zero branching ratios in the current run, and the number with branching ratios larger than $10^{-3}$. The second table which is printed by this call merely lists the $R$-violating $\lambda$, $\lambda'$, and $\lambda''$ couplings.

Purpose:

to transform an event listing between different reference frames, if so desired. The use of this routine assumes you do not do any boosts yourself.

IFRAME :

specification of frame the event is to be boosted to.

= 1 :

frame specified by you in the PYINIT call.

= 2 :

c.m. frame of incoming particles.

= 3 :

hadronic c.m. frame of lepton-hadron interaction events. Mainly intended for Deeply Inelastic Scattering, but can also be used in photoproduction. Is not guaranteed to work with the 'gamma/lepton' options, however, and so of limited use. Note that both the lepton and any photons radiated off the lepton remain in the event listing, and have to be removed separately if you only want to study the hadronic subsystem.

Purpose:

to enable you to reject a given set of kinematic variables at an early stage of the generation procedure (before evaluation of 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 you require is by setting MSTP(141)=1, and never if `minimum-bias'-type events (including elastic and diffractive scattering) are to be generated as well. Furthermore it is never called for user-defined external processes. 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 common block 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 information available includes subprocess type ISUB, $E_{\mathrm{cm}}$, $\hat{s}$, $\hat{t}$, $\hat{u}$, $\hat{p}_{\perp}$, $x_1$, $x_2$, $x_{\mathrm{F}}$, $\tau$, $y$, $\tau'$, $\cos\hat{\theta}$, and a few more. Some of these may not be relevant for the process under study, and are then set to zero.