Four-jet production in single- and double-parton scattering within high-energy factorization

We perform a first study of 4-jet production in a complete high-energy factorization (HEF) framework. We include and discuss contributions from both single-parton scattering (SPS) and double-parton scattering (DPS). The calculations are performed for kinematical situations relevant for two experimental measurements (ATLAS and CMS) at the LHC. We compare our results to those reported by the ATLAS and CMS collaborations for different sets of kinematical cuts. The results of the HEF approach are compared with their counterparts for collinear factorization. For symmetric cuts the DPS HEF result is considerably smaller than the one obtained with collinear factorization. The mechanism leading to this difference is of kinematical nature. We conclude that an analysis of inclusive 4-jet production with asymmetric pT -cuts below 50 GeV would be useful to enhance the DPS contribution relative to the SPS contribution. In contrast to the collinear approach, the HEF approach nicely describes the distribution of the ΔS variable, which involves all four jets and their angular correlations.


Introduction
So far, complete (n ≥ 4)-jet production via single-parton scattering (SPS) was discussed only within collinear factorization. Results up to next-to-leading (NLO) precision can be found in [1,2]. Here we wish to discuss for the first time production of four jets within high-energy (k T -)factorization (HEF) approach with 2 → 4 subprocesses with two off-shell partons. Recently three of us have discussed another reaction with 2 → 4 (gg → cccc) subprocess in the framework of the HEF [3]. For the four-jet production the number of subprocesses is much higher.
Double-parton scattering (DPS) was claimed to have been observed for the first time at the Tevatron [4]. In the LHC era, with much higher collision energies available, the field has received a new impulse and several experimental and theoretical studies address the problem of pinning down DPS effects (for review see [5,6]).
Even just from the purely theoretical point of view, the problem is quite subtle. As for the non perturbative side, it is in principle necessary, when considering a double-parton scattering, to take into account the correlations between the two partons coming from the same protons and involved in the scattering processes. Such an information should be encoded in a set of double parton distribution functions (DPDFs), generalising usual parton distribution functions (PDFs). A benchmarking work on DPDFs was made in ref. [7], where a proper generalisation of the sum rules for DGLAP type evolution equations for DPDSs was provided. Building explicit initial conditions for the evolution equations is challenging. Some successful attempts are becoming to appear only recently [8][9][10]. In JHEP04(2016)175 the meanwhile, phenomenological and experimental studies of double-parton scattering rely on factorized ansatz for the DPDFs, which amount to neglecting longitudinal momentum correlations between partons and treating transversal ones by introducing an effective cross section, σ eff . The latter quantity is usually extracted from experimental data. In the present approach we will use the factorized ansatz and concentrate on the difference between leading-order collinear and high-energy-factorization results. The latter includes effectively higher-order corrections. For most of high-energy reactions the single-parton scattering dominates over the double-parton scattering. The extraordinary example is double production of cc pairs [11,12]. For four-jet production, disentangling the ordinary SPS contributions from the DPS corrections can be quite challenging for several reasons: first of all, it is necessary to define sufficiently sensitive, process-dependent obervables, w.r.t. which the DPS differential cross section manifestly dominates at least in some corners of phase space [13,14]. Nevertheless, even once this is done, one has to be careful about the kinematical regime employed in comparing experimental data to theoretical prediction: in fact, the generally decreasing behaviour of PDFs for large momentum fractions [15] is well known, particularly for gluons, and gluon-initiated processes account for a very large part of the cross section; this implies that, for very energetic final states (characterized by large transverse momenta), it is really unlikely to get contributions from DPS. This is confirmed very well experimentally by the data released by the ATLAS Collaboration for both the 7 and 8 TeV runs [16,17]. This problem is of course slightly tamed by providing high center of mass energy in hadron-hadron scattering, as moderately low values of x should be enough to guarantee observing DPS in a kinematic regime in which perturbative QCD, possibly supplemented by parton showering, still works reasonably well.
In this paper we propose to assess the predictions of HEF for double-parton scattering at the LHC in a leading-order (LO) framework. HEF is an approach introduced in the early 90's in the context of heavy-flavour production, in order to take into account the effect of the colliding parton transverse momentum, which is neglected in the collinear approach [18][19][20]. This implies using off-shell partons, for which the construction of gaugeinvariant scattering amplitudes is not straightforward. However, recent improvements in the understanding of scattering amplitudes have allowed to formulate efficient analytical and numerical algorithms for the computation of such objects [21][22][23][24][25][26][27][28].
With such a machinery, we expand the analysis presented in ref. [14] and assess the differences between the pure collinear approach and the high-energy factorization (HEF) called also k T -factorization framework. We shall focus on the difference between predictions of HEF and standard collinear approach for the DPS contribution.

Single-parton scattering production of four jets
The collinear factorization formula for the calculation of the inclusive partonic 4-jet cross section at the Born level reads Here x 1,2 f i (x 1,2 , µ F ) are the collinear PDFs for the i − th parton, carrying x 1,2 momentum fractions of the proton and evaluated at the factorization scale µ F ; the index l runs over the four partons in the final state, the partonic center of mass energy squared isŝ = 2 x 1 x 2 P i · P j ; the function Θ 4−jet takes into account the kinematic cuts applied and M is the partonic on-shell matrix element, which includes symmetrization effects due to the identity of the particles in the final state. Switching to HEF, the analogous formula to (2.1) looks as follows: Here F i (x k , k T k , µ F ) is a transverse momentum dependent (TMD) distribution function for a given type of parton. Similarly as in the collinear factorization case, x k is the longitudinal momentum fraction, µ F is a factorization scale. The new degrees of freedom are introduced via k T k , which are the parton's transverse momenta, i.e. the momenta perpendicular to the collision axis. The formula is valid when the x's are not too large and not too small when complications from nonlinearities may eventually arise [29]. 1 The TMD parton densities (for a recent review see [32]) can be defined by introducing operators whose expectation values, roughly speaking, count the number of partons [33]. In particular, an evolution equation for TMDs known as CCFM, valid both in the low x and large x domain, [34,35] provides a gluon density depending on x, k T , µ. However, for our purposes this is not enough, since we want to have access to moderate values of x where the CCFM approach needs refinements [36,37]. The alternative is to use the Kimber-Martin-Ryskin (KMR) prescription [38,39] in order to obtain a full set of TMD parton densities. The basic observation is that the k T dependence can be generated at the very last step of the collinear evolution by performing soft gluon resummation between two scales given by k T and µ F , where k T is interpreted as the transverse momentum of the hardest emitted gluon during the partonic evolution, while µ F can be linked to the hard scattering scale. In practical terms, this procedure boils down to applying the Sudakov form factor onto the PDFs (some details can be found in appendix A 2 ). M(i * , j * → 4 part.) is the gauge invariant matrix element for 2 → 4 particle scattering with two initial off-shell legs. In the case of HEF (for recent review see ref. [40]), amplitudes with external off-shell legs in QCD have been computed with different approaches: up to 2 → 2 scattering, they are given for example in [41] and are enough in order to calculate DPS contributions (see section 3). In order to move on to higher multiplicities, which are necessary for the SPS analysis of 2 → 4 parton scattering, it is possible to generate this JHEP04(2016)175 amplitudes analytically applying suitably defined Feynman rules [21,22]. Also recursive methods have been developed for this purpose, like generalised BCFW recursion [24,25] or Wilson lines approaches [23,26,28]. By now also a numerical package implementing numerical BCFW recursion is available [42]. In this case, we rely on a numerical approach implemented in AVHLIB 3 which employs Dyson-Schwinger recursion generalized to treelevel amplitudes with off-shell initial-state particles. Originally proposed in [43,44], this recursive method exists in several explicit implementations with on-shell initial-state particles [45][46][47][48][49], and has even been extended to one-loop amplitudes [50,51]. AVHLIB and the Monte Carlo program therein are also used to perform the phase-space integration. In the collinear case, results were cross-checked by comparing them with the ALPGEN output [46]. We use a running α s provided with the MSTW2008nlo68cl PDF sets and set both the renormalization and factorization scales equal to half the transverse energy, which is defined as the sum of the final state transverse momenta, µ F = µ R =Ĥ T 2 = 1 2 4 l=1 k l T , 4 working in the n F = 5 flavour scheme.
In order to cross-check our numerical tools, we must compare their outputs to results already available in the literature. For this purpose, we compared LO total cross sections for (n ≥ 4)-jet production to those given by the BlackHat collaboration in ref. [1] and crosschecked in ref. [2]. We find excellent agreement, up to phase space integration accuracy. The cuts used in these calculations were those chosen by the ATLAS collaboration in the 2011 analysis of multi-jet events [16], namely p T > 80 GeV for the leading jet and p T > 60 GeV for subleading jets, |η| < 2.8 for the pseudorapidity and jet cone radius parameter ∆R > 0.4. Again we find excellent agreement between the two codes with the LO results reported in the literature, up to phase space integration uncertainties.
To be precise, we reproduce the LO predictions for the total inclusive cross sections where the numbers in brackets stand for the numerical integration uncertainty and the upper and lower errors are obtained by varying the renormalization scale up and down by a factor of two. There are 19 different channels contributing to the cross section at the parton-level: The processes in the first line are the dominant channels, contributing together to ∼ 93% of the total cross section. This stays true in the HEF framework as well.

JHEP04(2016)175
3 Double-parton scattering production of four jets Single-parton scattering contributions are expected to be dominant for high momentum transfer, as it is highly unlikely that two partons from one proton and two from the other one are energetic enough for two hard scatterings to take place, as the behaviour of the PDFs for large x suggests. However, as the cuts on the transverse momenta of the final state are softened, a window opens to possibly observe significant double parton scattering effects, as often stated in the literature on the subject and recently analysed for 4-jet production in the collinear approach in ref. [14]. Our goal here is to perform the same analysis in HEF, in order to assess the difference in the predictions of the two approaches.
First of all, let us recall the formula usually employed for the computation of DPS cross sections, adjusting it to the 4-parton final state, where the σ(ab → cd) cross sections are obtained by restricting formulas (2.1) and (2.2) to a single channel and the symmetry factor m is 1 unless the two hard scatterings are identical, in which case it is 1/2, so as to avoid double counting. Above ξ 1 and ξ 2 stand for generic kinematical variables for the first and second scattering, respectively.
The effective cross section σ eff can be loosely interpreted as a measure of the transverse correlation of the two partons inside the hadrons, whereas the possible longitudinal correlations are usually neglected (for an introduction to this issue, see for example ref. [7]). In this paper we use σ eff provided by the CDF, D0 collaborations and recently confirmed by the LHCb collaboration σ eff = 15 mb, although the latter value may be questioned [52] when all SPS mechanisms of double charm production are included.
As already mentioned in the introduction there are attempts, in the literature, to construct DPDFs which include correlations also between the longitudinal momenta of the two partons and fullfil sum rules. These models are, however, still rather at a preliminary stage. So far they are formulated exlusively in the gluon sector [8] or in the valence quark sector [9]. In addition they are only formulated in a leading order framework which may be not sufficient for many processes. Moreover, as it is expected on physical grounds and confirmed by all the calculations in the various models proposed so far, the longitudinal parton-parton correlations should become far less important as the energy of the collision is increased, due to the increase in the parton multiplicity. For instance, the plots in ref. [8] show that the double gluon distribution obtained with a sum rule approach is essentially equal to the factorized ansatz at the scale Q 2 = 10 GeV 2 down to x = 10 −5 . Looking forward to further improvements in this field, we choose to limit ourselves to a more pragmatic approach for the purpose of this paper, making the following ansatz for DPDF in the collinear-factorization case: is the DPDF and f i (x i , µ) are the ordinary PDFs and the subscripts 1 and 2 simply differentiate between two generic partons in the same proton. This ansatz can be automatically generalised to the case when parton transverse momenta are included.
Coming to DPS contributions, we have to include all the possible 45 channels which can be obtained by coupling in all possible distinct ways the 8 channels for the 2 → 2 SPS process, i.e.

Comparison to the collinear approach and to ATLAS data with hard central kinematic cuts
In the following, we test the HEF calculation against the collinear case and compare it to the 8 TeV data recently reported by the ATLAS collaboration [17]. The kinematic cuts are here slightly different with respect to ref. [16]: p T > 100 GeV for the leading jet and p T > 64 GeV for the first three subleading jets; in addition |η| < 2.8 is the pseudorapidity cut and ∆R > 0.5 is the constraint on the jet cone radius parameter.
As for this framework, we employ, together with the newly obtained TMD PDFs (5 quark flavors and gluon) which we call DLC-2016 (Double-Log-Coherence), the running α s coming with the MSTWnlo200868cl sets. The results of our computation in HEF is shown in figures 1 and 2, where it is apparent that the DPS contribution is completely irrelevant, as expected for final states with high transverse momenta, as it is extremely unlikely that all the four partons in the two couples coming from the colliding protons carry enough energy to produce such a hard final state. A generally good agreement with the ATLAS data can be seen through the transverse momenta spectra of the four jets, thus showing that the HEF approach works reliably well in this region.
First we show the results of the HEF calculation in figures 1 and 2. The prediction is consistent with the ATLAS data for all the p T spectra.
Next we assess the difference between the HEF and collinear predictions at LO as far as SPS is concerned. We see from figures 3 and 4 that collinear factorization performs slightly better for intermediate values and HEF does a better job for the last bins, except for the 4th jet. All in all, both approaches are consistent with the data in this kinematic region.

Comparison to CMS data with softer cuts
As discussed in ref. [14], so far the only experimental analysis of four-jet production relevant for the DPS studies was realized by the CMS collaboration [53]. The cuts used in this analysis are p T > 50 GeV for the first and second jets, p T > 20 GeV for the third and  Figure 1. HEF prediction of the differential cross sections w.r.t. the transverse momenta of the first two leading jets compared to the ATLAS data [17]. The LO calculation describes the data pretty well in this hard regime in which MPIs are irrelevant. In addition we show the ratio of the SPS HEF result to the ATLAS data.  Figure 2. HEF prediction of the differential cross sections w.r.t. the transverse momenta of the 3rd and 4th leading jets compared to the ATLAS data [17]. The LO calculation describes the data pretty well in this hard regime in which MPIs are irrelevant. In addition we show the ratio of the SPS HEF result to the ATLAS data.
fourth jets, |η| < 4.7 and the jet cone radius parameter ∆R > 0.5. In the rest of this section, we present our results for such cuts.
As for the total cross section for the four jet production, the experimental and theoretical LO results are respectively:  the HEF result, compared to collinear factorization, is the dramatic damping of the DPS contribution. The effect of the damping is of kinematical nature and will be explained below.
The effect of the relative damping of the HEF DPS result compared to the leading-order collinear DPS result is of kinematical origin. The main idea can be understood already in a bit simpler case of two-jet production within the HEF approach. For the purpose of this illustration we impose a cut p T > 35 GeV on both jets (leading and subleading). In figure 5 we show the transverse momentum distributions for both leading (long-dashed line) and subleading (long dashed-dotted line) jet. We observe a minimum for the leading jet and maximum for the subleading jet for transverse momenta in the vicintity of the lower cut. The integrated cross section for the leading and subleading jet is of course identical as they are "measured" (identified) in coincidence. For the leading order collinear case both jets have the same distribution and one gets maximum in the vicinity of the transverse momentum threshold in both cases. In this case imposing cuts on both jets does not lead   to "loosing" cross section. In contrast, in the HEF approach, if the leading jet is close to the transverse momentum threshold, then the subleading jet is typically below that threshold, therefore such an event is not counted.
For four-jet DPS production the situation is more complicated and strongly depends on cuts (all identical, two pairs of identical cuts, harder cut for the leading jet and identical for the other, etc.).
In figures 6 and 7 we compare the predictions in HEF to the CMS data. Here both the SPS and DPS contributions are normalized to the total cross section, i.e. the sum of the SPS and DPS contributions. In all cases the renormalized transverse momentum distributions agree with the CMS data. However, the absolute cross sections obtained in this case within the HEF approach are too large.
Not only transverse momentum dependence is interesting. The CMS collaboration extracted also more complicated observables [53]. One of them, which involves all four jets in the final state, is the ∆S variable, defined in ref. [53]  harder and the softer jets, where p T (j i , j k ) stands for the sum of the transverse momenta of the two jets in arguments.
In figure 8 we present our HEF prediction for the normalized to unity distribution in the ∆S variable. Our HEF result approximately agrees with the experimental ∆S distribution. In contrast the LO collinear approach leads to ∆S = 0, i.e. a Kronecker-delta peak at ∆S = 0 for the distribution in ∆S. For the DPS case this is rather trivial. The two hard and two soft jets come in this case from the same scatterings and are back-to-back (LO), so each term in the argument of arccos is zero (jets are balanced in transverse momenta). For the SPS case the transverse momenta of the two jet pairs (with hard jets and soft jets) are identical and have opposite direction (the total transverse momentum of all four jets must be zero from the momentum conservation). Then it is easy to see that the argument of arccos is just -1. This means ∆S = 0. The above relations are not fullfilled in the HEF approach. The SPS contribution clearly dominates and approximately gives the shape of the ∆S distribution. The DPS contribution improves the agreement with the data in the central region, worsening it a little bit for ∆S → 0 while essentially leaving the result unaffected for ∆S → π. It is anyway interesting that we roughly describe the data via pQCD effects within our HEF approach which are in ref. [53] described by parton-showers and soft MPIs.
It would be nice to have more insight into our successful description of the ∆S distribution measured by the CMS experiment. In figure 9 we return to the ∆S spectrum and show also two results with a TMD toy model with the Gaussian smearing of the collinear parton distribution: We take two sets of smearing parameter: σ = 1 GeV (left panel) and 5 GeV (right panel). Taking a bigger value of σ we approach the CMS data. This shows that the transverse momenta bigger than a few GeV are needed to approach the data.

HEF predictions for a possible set of asymmetric cuts
Moving from the previous considerations, in the following subsection we present our results of the DPS employing asymmetric cuts by which we mean here p T > 35 GeV for the leading jet, p T > 20 GeV for the other jets and |η| < 4.7, ∆R > 0.5. Of course it would be interesting to have the results of such an experimental analysis (i.e. with soft enough but asymmetric cuts ) in order to test the predictions of HEF for DPS.
In this case the theoretical total cross sections for four-jet production are:

JHEP04(2016)175 4 Conclusions
In the present paper we have compared the perturbative predictions for four-jet production at the LHC in leading-order collinear and high-energy (k T -)factorization. Both single parton and double parton scattering contributions have been calculated for a first time in the high-energy (k T -)factorization approach. The calculation of the SPS contribution may be considered as a technical achievment. So far only production of the cccc final state (also of the 2 → 4 type) was discussed in the literature in this context. For the four-jet production the number of relevant subproceses is much larger but could be treated in our automatized framework. We find that both the collinear and the (k T -)factorization approaches describe the data for hard central cuts, relevant for the ATLAS experiment, reasonably well. For the harder cuts we get both normalization and shape of the transverse momentum distributions. For the softer cuts used e.g. by the CMS collaboration the tree level result is unreliable. Therefore in this case we have presented results for normalised cross sections. We have presented distributions in transverse momenta for all jets ordered in their transverse momenta. We have found that for symmetric cuts the DPS cross section obtained with more realistic high energy (k T -)factorization approach is smaller than the one obtained in the collinear approach. This is an important result for searches of DPS effects in four-jet production not discussed so far in the literature. We have tried to explain this as kinematical effect due to phase space limitation when simultaneously imposing cuts on all jets but a full explanation is a bit intricate. While we observe, in agreement with ref. [14], that lowering the cut in transverse momenta can significantly enhance the experimental sensitivity to DPS, we also observe that the HEF approach severely tames this effect for symmetric cuts, due to gluon-emission effects which alter the transverse-momentum balance between final state partons. We have found that the damping is not present when cuts are not identical. The discussion how to optimize the cuts will be presented elsewhere.
For other approaches addressing the four-jet production and resummation of BFKL type of singularities [54,55] we refer the Reader to ref. [56]. The authors of ref. [56] define new angular four jet observables to test BFKL approach.
As a side result, we present in appendix B a detailed numerical comparison of results obtained for dijet production with matrix element generated automatically by means of AVHLIB with those obtained analytically within the Parton-Reggeization-Approach (PRA) in ref. [41] and implemented in an independent code. For all (sub)processes we have obtained very good agreement of the corresponding differential distributions. We show corresponding azimuthal angle correlations for different subprocesses which are particularly efficient for such tests, as they sample the situation in a broad range of the phase space. It was shown in the past for some subprocesses that also analytical results coincide [21].

A Construction of the TMDs
The Kimber-Martin-Ryskin (KMR) analysis of coherence effects provides in the DLL limit the following formula for obtaining TMD parton density function (this limit is relevant for us since we consider moderate values of x and rather large scales): where k 2 T min = µ 2 min = 1.7 GeV 2 and T i is appropriate Sudakov form factor. In the gluon channel it reads T g (k T , µ) = exp − µ 2 k 2 T α s (µ 2 ) 2π µ µ+p t 0 dz (z P gg (z ) + n F P qg (z )) (A.2) while for quarks we have: In our calculations we take for xf i the CTEQ10NLO set of pdfs, the P ij are LO splitting functions which to the accuracy we work with are sufficient and in practical calculations we take n F = 5.

B Comparison of 2 → 2 matrix elements
Here we perform a comparison of the cross sections for dijet production obtained both with the QCD amplitudes from the Parton-Reggeization Approach (PRA) [41] and those  Figure 14. Comparison of the AVHLIB and PRA predictions for azimuthal angle correlation between jets. The calculations for g * g * → gg (left panel) and q * q * → gg (right panel) subprocesses are shown as an example. from the AVHLIB, which were in turn cross checked with the results presented in refs. [24,25]. We find perfect agreement modulo phase space integration uncertainty, as shown in figures 13 and 14. The consistency between amplitudes computed in different approaches is, by itself, a non trivial check of the methods employed. Here we have shown only a few examples. In fact we have checked that the agreement holds for all the processes.
Open Access. This article is distributed under the terms of the Creative Commons Attribution License (CC-BY 4.0), which permits any use, distribution and reproduction in any medium, provided the original author(s) and source are credited.