Top Processes

flag  Top:all   (default = off)
Common switch for the group of top production.

flag  Top:gg2ttbar   (default = off)
Scatterings g g → t tbar. Code 601.

flag  Top:qqbar2ttbar   (default = off)
Scatterings q qbar → t tbar by gluon exchange. Code 602.

flag  Top:qq2tq(t:W)   (default = off)
Scatterings q q' → t q'' by t-channel exchange of a W^+- boson. Code 603.

flag  Top:ffbar2ttbar(s:gmZ)   (default = off)
Scatterings f fbar → t tbar by s-channel exchange of a gamma^*/Z^0 boson. Code 604.

flag  Top:ffbar2tqbar(s:W)   (default = off)
Scatterings f fbar' → t q'' by s-channel exchange of a W^+- boson. Code 605.

flag  Top:gmgm2ttbar   (default = off)
Scatterings gamma gamma → t tbar. Code 606.

flag  Top:ggm2ttbar   (default = off)
Scatterings g gamma → t tbar. Code 607 when g gamma → t tbar and 617 when gamma g → t tbar.

By default top always decays to a W and a down-type quark. It is possible to switch on the t → H+ b decay mode. Note that its partial width is calculated using the tan(beta) value stored in HiggsHchg:tanBeta, so that it can be used without having to read in a SUSY parameter file. For the H+ to decay also Higgs:useBSM = on is necessary.

Threshold enhancements

Relevant code is implemented as a new class TopThreshold in SigmaQCD.h/.cc, which is accessed by the internal gg → ttbar and qqbar → ttbar classes in the same files. The implemented scenarios and free parameters within them are as follows.

mode  TopThreshold:model   (default = 0; minimum = 0; maximum = 2)
The choice of threshold behaviour for the g g → t tbar and q qbar → t tbar processes.
option 0 : no modifications to threshold, i.e. pure Born (leading order) matrix elements.
option 1 : the simple Coulomb enhancement, which only works above the threshold, E = mHat(t + tbar) - m(t) - m(tbar) > 0.
option 2 : the Green's function in the whole threshold region, both E > 0 and E < 0, but damped-out below (see thrRegion below). The procedure is to first pick m(t) and m(tbar) by Breit-Wigners, and then an mHat distributed all the way down to m(t) + m(tbar) - 2 * thrRegion. If E = mHat - m(t) - m(tbar) < 0 then a new m'(t) < m(t) and a new m'(tbar) < m(tbar) are picked according to Breit-Wigners under the requirement that E' = mHat - m'(t) - m'(tbar) > 0. Finally the event is accepted or rejected according to the naive cross section reweighted to the Green's function value.
Note: In PYTHIA 8.317 the Green's function expression has been modified in a small but not irrelevant way. The variable E = mHat - 2 m_t (for equal t and tbar masses) is intended to be proportional to the squared velocity beta^2, but that only holds very close to threshold. By using the correct full expression, where beta^2 = 1 - 4m_t^2 / sHat, agreement is obtained with the Coulomb expression some distance above the threshold.
Warning: The first edition of the threshold model framework, in PYTHIA 8.316, contained two intermediate scenarios, then options 2 and 3, that acted as controls and stepping stones towards the recommended option 4. In the current version those intermediate scenarios have been removed, and what used to be option 4 is now option 2. On the other hand, this scenario has been equipped with new variants; read on.

One aspect that was not fully appreciated in the PYTHIA 8.316 implementation is the role of the top quark width. In PYTHIA the top quark has a Breit-Wigner mass distribution, whether the top decay is explicitly considered or not, since the subsequent decays are implicitly going to happen later on. In matrix element calculations a production of a top pair has to be performed with the same mass for the top and the antitop, to preserve gauge invariance, and then the natural choice is the on-shell value. If you want to include the top width you should instead consider the production of its decay products, i.e. b W+ bbar W-, and to also include the W widths their decay products have to be swapped in. The way PYTHIA addresses this issue is to generate separate top and antitop Breit-Wigner masses, but use the event-by-event average in the production matrix elements. More specifically, the momentum three-vector is preserved in the choice of average, which translates into (1 - r_3 - r_4)^2 - 4 * r_3 * r_4 = 1 - 4 * r_avg , where r_3 = m_t^2/sHat, r_4 = m_tbar^2/sHat, and r_avg = m_avg^2/sHat.

The question is now how the top mass issue is handled in the derivation of the threshold cross section expressions. This issue is not addressed in [Fad90], which only does numerical integration of the below-threshold contribution, not Monte Carlo implementation. Going back to more detailed presentation, such as [Fad91], it becomes apparent that the relevant Green's function expressions of option 2 above are intended to include smearing from the top width. So it appears that there is doublecounting, that should be fixed. One could imagine two approaches to this.
(a) Reduce the top quark width to a small number. This is then closer to the general thinking of the theory community, as noted above. The main issue is that the broadening of top quarks above threshold is also inhibited, which is unphysical. Technically, the top width cannot be made arbitrarily small, since a finite width is needed for the handling of the "toponium" below-threshold decay to off-shell top quarks. The top width changed here is the regular one stored in the particle data tables, and can be set as property mWidth of particle 6. Note, however, that PYTHIA will use the chosen top mass to recalculate the top width at initialization, so to inhibit this you must also set 6:doForceWidth = true.
(b) Reduce the top width used in the threshold-behaviour matrix element expressions. This then provides a natural continuation of the above-threshold behaviour traditionally used in PYTHIA, including model options 0 and 1 above, to the below-threshold one. The problem here is a practical one: a shrinking matrix element width gives a spectrum of discrete bound-state peaks in the t+tbar mass spectrum, which becomes messy to handle properly. Notably the initialization step may not fully catch the peaks, which later leads to some events receiving weights above unity to compensate. This kind of top width is regulated by the tWidthGreen parameter below.

parm  TopThreshold:tWidthGreen   (default = 1.34; minimum = 0.1; maximum = 2.)
top quark width as used in the Green's function. See description above for the issue this parameter is intended to address.

parm  TopThreshold:thrRegion   (default = 10.; minimum = 5.; maximum = 20.)
the Green's function, when used, is only assumed valid down to E = -thrRegion. In the region [-2 thrRegion, -thrRegion] it is linearly damped to zero, and below that identically zero.
Note: renamed from width in PYTHIA 8.316, to avoid confusion with the width concept used e.g. in tWidthGreen above. Also, owing to a changed handling of the high-energy behaviour in model 2, damping at high masses is no longer required.

mode  TopThreshold:alphasOrder   (default = 2; minimum = 0; maximum = 2)
the order of the running of the alpha_strong, used for the top threshold factors (and nowhere else).
option 0 : no running.
option 1 : first-order running.
option 2 : second-order running.

parm  TopThreshold:alphasValue   (default = 0.118; minimum = 0.10; maximum = 0.25)
the alpha_strong value at scale M_Z^2, that then runs according to the order defined above.

parm  TopThreshold:ggSingletFrac   (default = 0.28571; minimum = 0.; maximum = 1.)
in the g g → t tbar process, colour factors gives 2/7 singlet and the rest octet, but dynamics might modify this.

parm  TopThreshold:qqSingletFrac   (default = 0.; minimum = 0.; maximum = 1.)
in the q qbar → t tbar process colour arguments gives all octet and no singlet, but again modifications are possible.

flag  TopThreshold:pseudoscalar   (default = off)
distribute angles of the consecutive t, tbar W+ and W- decays as if the t+tbar pair is in a pseudoscalar state, giving specific angular correlations otherwise absent. This only applies in the threshold region, while at higher energies a transition to independent t and tbar decays is assumed, see next two parameters. When on it affects all top threshold models, also the default pure Born level.

parm  TopThreshold:psEBeginDamp   (default = 10.)
when TopThreshold:pseudoscalar = on then the correlated decays are used up until this many GeV above the threshold energy. After that a linear damping is introduced, with the slack being taken over by independent t and tbar decays, until the latter completely takes over at the next scale.

parm  TopThreshold:psEEndDamp   (default = 20.)
the linear damping introduced in the previous parameter ends at the current one, so above this scale (relative to the threshold energy) only independent t and tbar decays remain. Must be chosen larger than the previous parameter.

Note 1: Model 2, for E < 0 redefines the t and tbar mass values in order to achieve a new E' > 0. In order still to have access to the original quantities, before the new masses were selected, the following quantities are saved and retrievable:
pythia.info.toponiumE the original threshold energy E, which may have either sign;
pythia.info.toponiumm3, pythia.info.toponiumm4 the two original t and tbar masses, which are larger than the m' ones found in the event record.

Note 2: Five main programs explore the usage of this code:
main368.cc plots the shape of the pure singlet and octet enhancements;
main369.cc and main370.cc histogram a few quantities for many different scenarios, to allow direct comparisons, with code for parallelization using either OpenMP or PythiaParallel, respectively;
main371.cc histograms a wider set of quantities, but only for one model at a time;
main372.cc studies angular correlations between outgoing fermions.