In the spirit of section , we want to select simple functions such that the true , and dependence of the cross sections is approximately modelled. In particular, (almost) all conceivable kinematical peaks should be represented by separate terms in the approximate formulae. If this can be achieved, the ratio of the correct to the approximate cross sections will not fluctuate too much, but allow reasonable Monte Carlo efficiency.
Therefore the variables are generated according to the distributions
, and , where normally
The first variable to be selected is . The range of allowed values, , is generally constrained by a number of user-defined requirements. A cut on the allowed mass range is directly reflected in , a cut on the range indirectly. The first two terms of are intended to represent a smooth dependence, as generally obtained in processes which do not receive contributions from -channel resonances. Also -channel exchange of essentially massless particles (, , light quarks and leptons) are accounted for, since these do not produce any separate peaks at non-vanishing . The last four terms of are there to catch the peaks in the cross section from resonance production. These terms are only included when needed. Each resonance is represented by two pieces, a first to cover the interference with graphs which peak at , plus the variation of parton distributions, and a second to approximate the Breit-Wigner shape of the resonance itself. The subscripts and denote values pertaining to the two resonances, with . Currently there is only one process where the full structure with two resonances is used, namely . Otherwise either one or no resonance peak is taken into account.
The kinematically allowed range of values is constrained by , , and you may impose additional cuts. Therefore the allowed range is only constructed after has been selected. The first two terms of give a fairly flat dependence -- for processes which are symmetric in , they will add to give a completely flat spectrum between the allowed limits. In principle, the natural subdivision would have been one term flat in and one forward-backward asymmetric, i.e. proportional to . The latter is disallowed by the requirement of positivity, however. The and terms actually used give the same amount of freedom, but respect positivity. The third term is peaked at around , and represents the bias of parton distributions towards this region.
The allowed range is naïvely . However, most cross sections are divergent for , so some kind of regularization is necessary. Normally one requires , which translates into for massless outgoing particles. Since again the limits depend on , the selection of is done after that of . Additional requirements may constrain the range further. In particular, a constraint may split the allowed range into two, i.e. or . An un-split range is represented by . For massless outgoing particles the parameter in , such that the five terms represent a piece flat in angle and pieces peaked as , , , and , respectively. For non-vanishing masses one has . In this case, the full range is therefore available -- physically, the standard and singularities are regularized by the masses and .
For each of the terms, the coefficients represent the
integral over the quantity multiplying the coefficient ; thus,
The coefficients are normalized to unit sum for , and separately. They have a simple interpretation, as the probability for each of the terms to be used in the preliminary selection of , and , respectively. The variation of the cross section over the allowed phase space is explored in the initialization procedure of a PYTHIA run, and based on this knowledge the are optimized so as to give functions , and that closely follow the general behaviour of the true cross section. For instance, the coefficient in is to be made larger the more the total cross section is dominated by the region around the resonance mass.
The phase-space points tested at initialization are put on a grid, with the number of points in each dimension given by the number of terms in the respective expression, and with the position of each point given by the median value of the distribution of one of the terms. For instance, the distribution gives a median point at , and has the median . Since the allowed and ranges depend on the value selected, then so do the median points defined for these two variables.
With only a limited set of phase-space points studied at the initialization, the `optimal' set of coefficients is not uniquely defined. To be on the safe side, 40% of the total weight is therefore assigned evenly between all allowed , whereas the remaining 60% are assigned according to the relative importance surmised, under the constraint that no coefficient is allowed to receive a negative contribution from this second piece.
After a preliminary choice has been made of , and ,
it is necessary to find the weight of the event, which is to be used
to determine whether to keep it or generate another one.
Using the relation
eq. () may be rewritten as
As can be seen, the cross section for a given Monte Carlo event is given as the product of four factors, as follows:
Of course, the expression in the last line is not strictly necessary to obtain the cross section by Monte Carlo integration. One could also have used eq. () directly, selecting phase-space points evenly in , and , and averaging over those Monte Carlo weights. Clearly this would be much simpler, but the price to be paid is that the weights of individual events could fluctuate wildly. For instance, if the cross section contains a narrow resonance, the few phase-space points that are generated in the resonance region obtain large weights, while the rest do not. With our procedure, a resonance would be included in the factor, so that more events would be generated at around the appropriate value (owing to the numerator in the phase-space expression), but with these events assigned a lower, more normal weight (owing to the factor in the weight expression). Since the weights fluctuate less, fewer phase-space points need be selected to get a reasonable cross-section estimate.
In the program, the cross section is obtained as the average over all phase-space points generated. The events actually handed on to you should have unit weight, however (an option with weighted events exists, but does not represent the mainstream usage). At initialization, after the coefficients have been determined, a search inside the allowed phase-space volume is therefore made to find the maximum of the weight expression in the last line of eq. (). In the subsequent generation of events, a selected phase-space point is then retained with a probability equal to the weight in the point divided by the maximum weight. Only the retained phase-space points are considered further, and generated as complete events.
The search for the maximum is begun by evaluating the weight in the same grid of points as used to determine the coefficients. The point with highest weight is used as starting point for a search towards the maximum. In unfortunate cases, the convergence could be towards a local maximum which is not the global one. To somewhat reduce this risk, also the grid point with second-highest weight is used for another search. After initialization, when events are generated, a warning message will be given by default at any time a phase-space point is selected where the weight is larger than the maximum, and thereafter the maximum weight is adjusted to reflect the new knowledge. This means that events generated before this time have a somewhat erroneous distribution in phase space, but if the maximum violation is rather modest the effects should be negligible. The estimation of the cross section is not affected by any of these considerations, since the maximum weight does not enter into eq. ().
For processes with identical final-state particles, the symmetrization factor of is explicitly included at the end of the calculation. 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). That way, no special provisions are needed in the phase-space integration machinery.