Event Information
The Info
class collects various one-of-a-kind information,
some relevant for all events and others for the current event.
An object info
is a public member of the Pythia
class, so if you e.g. have declared Pythia pythia
, the
Info
methods can be accessed by
pythia.info.method()
. Most of this is information that
could also be obtained e.g. from the event record, but is here more
directly available. It is primarily intended for processes generated
internally in PYTHIA, but many of the methods would work also for
events fed in via the Les Houches Accord.
Here are the currently available methods related to each event:
method
list()
a listing of most of the information set for the current event.
method
idA(), idB()
the identities of the two beam particles.
method
pzA(), pzB()
the longitudinal momenta of the two beam particles.
method
eA(), eB()
the energies of the two beam particles.
method
mA(), mB()
the masses of the two beam particles.
method
eCM(), s()
the cm energy and its square for the two beams.
method
name(), code()
the name and code of the process that occured.
method
nFinal()
the number of final-state partons in the hard process.
method
isResolved()
are beam particles resolved, i.e. were PDF's used for the process?
method
isDiffractiveA(), isDiffractiveB()
is either beam diffractively excited?
method
isMinBias()
is the process a minimum-bias one?
method
isLHA()
has the process been generated from external Les Houches Accord
information?
method
atEndOfFile()
true if a linked Les Houches class refuses to return any further
events, presumably because it has reached the end of the file from
which events have been read in.
method
hasSub()
does the process have a subprocess classification?
Currently only true for minbias and Les Houches events, where it allows
the hardest collision to be identified.
method
nameSub(), codeSub(), nFinalSub()
the name, code and number of final-state partons in the subprocess
that occured when hasSub()
is true. For a minimum-bias event
the code
would always be 101, while codeSub()
would vary depending on the actual hardest interaction, e.g. 111 for
g g -> g g. For a Les Houches event the code
would
always be 9999, while codeSub()
would be the external
user-defined classification code. The methods below would also provide
information for such particular subcollisions.
method
id1(), id2()
the identities of the two partons coming in to the hard process.
method
x1(), x2()
x fractions of the two partons coming in to the hard process.
method
y(), tau()
rapidity and scaled mass-squared of the hard-process subsystem, as
defined by the above x values.
method
pdf1(), pdf2()
parton densities x*f(x,Q^2 )evaluated for the two incoming
partons; could be used e.g. for reweighting purposes.
method
QFac(), Q2Fac()
the Q^2 or Q^2 factorization scale at which the
densities were evaluated.
method
isValence1(), isValence2()
true
if the two hard incoming partons have been picked
to belong to the valence piece of the parton-density distribution,
else false
. Should be interpreted with caution.
Information is not set if you switch off parton-level processing.
method
alphaS(), alphaEM()
the alpha_strong and alpha_electromagnetic values used
for the hard process.
method
QRen(), Q2Ren()
the Q or Q^2 renormalization scale at which
alpha_strong and alpha_electromagnetic were evaluated.
method
mHat(), sHat()
the invariant mass and its square for the hard process.
method
tHat(), uHat()
the remaining two Mandelstam variables; only defined for 2 -> 2
processes.
method
pTHat(), pT2Hat()
transverse momentum and its square in the rest frame of a 2 -> 2
processes.
method
m3Hat(), m4Hat()
the masses of the two outgoing particles in a 2 -> 2 processes.
method
thetaHat(), phiHat()
the polar and azimuthal scattering angles in the rest frame of
a 2 -> 2 process.
method
weight()
weight assigned to the current event. Is normally 1 and thus uninteresting.
However, for Les Houches events some strategies allow negative weights,
which then after unweighting lead to events with weight -1. There are also
strategies where no unweighting is done, and therefore a nontrivial event
weight must be used e.g. when filling histograms.
method
bMI()
the impact parameter b assumed for the current collision when
multiple interactions are simulated. Is not expressed in any physical
size (like fm), but only rescaled so that the average should be unity
for minimum-bias events (meaning less than that for events with hard
processes).
method
enhanceMI()
The choice of impact parameter implies an enhancement or depletion of
the rate of subsequent interactions, as given by this number. Again
the average is normalized be unity for minimum-bias events (meaning
more than that for events with hard processes).
method
nMI()
the number of hard interactions in the current event. Is 0 for elastic
and diffractive events, and else at least 1, with more possible from
multiple interactions.
method
codeMI(i), pTMI(i)
the process code and transverse momentum of the i
'th
subprocess, with i
in the range from 0 to
nMI() - 1
. The values for subprocess 0 is redundant with
information already provided above.
method
nISR(), nFSRinProc(), nFSRinRes()
the number of emissions in the initial-state showering, in the final-state
showering excluding resonance decys, and in the final-state showering
inside resonance decays, respectively.
Here are the currently available methods related to the event sample
as a whole. While continuously updated during the run, it is recommended
only to study these properties at the end of the event generation,
when the full statistics is available.
method
nTried(), nSelected(), nAccepted()
the total number of tried phase-space points, selected hard processes
and finally accepted events, summed over all allowed subprocesses.
The first number is only intended for a study of the phase-space selection
efficiency. The last two numbers usually only disagree if the user introduces
some veto during the event-generation process; then the former is the number
of acceptable events found by PYTHIA and the latter the number that also
were approved by the user. If you set a
second hard process there may also be a mismatch.
method
sigmaGen(), sigmaErr()
the estimated cross section and its estimated error,
summed over all allowed subprocesses, in units of mb. The numbers refer to
the accepted event sample above, i.e. after any user veto.