Diagrammatic resummation of leading-logarithmic threshold effects at next-to-leading power

Perturbative cross-sections in QCD are beset by logarithms of kinematic invariants, whose arguments vanish when heavy particles are produced near threshold. Contributions of this type often need to be summed to all orders in the coupling, in order to improve the behaviour of the perturbative expansion, and it has long been known how to do this at leading power in the threshold variable, using a variety of approaches. Recently, the problem of extending this resummation to logarithms suppressed by a single power of the threshold variable has received considerable attention. In this paper, we show that such next-to-leading power (NLP) contributions can indeed be resummed, to leading logarithmic (LL) accuracy, for any QCD process with a colour-singlet final state, using a direct generalisation of the diagrammatic methods available at leading power. We compare our results with other approaches, and comment on the implications for further generalisations beyond leading-logarithmic accuracy.


Introduction
Perturbative calculations of hadronic cross sections in Quantum Chromodynamics (QCD) are the cornerstone of theoretical predictions for all processes of phenomenological interest at particle colliders, such as the Large Hadron Collider (LHC).Furthermore, the everincreasing precision of experimental data demands that theoretical predictions for scattering processes of interest be continually improved.The relevant calculations are carried out using an expansion in powers of the coupling constant α s , and typically proceed on two fronts.First, one may determine the complete behaviour of a given quantity at a fixed order in the coupling constant.The state of the art for most processes is next-to-leading order (NLO) in perturbation theory, with an increasing number of notable exceptions known at NNLO, and even N 3 LO (see e.g.[1] for a review).Whilst successful for many observables, the fixed order approach is only valid provided subleading perturbative corrections are well-behaved.Given that perturbative coefficients depend on the momenta of the scattering particles, this criterion can fail in certain kinematic regimes: a well-known example is the production of heavy particles near threshold.In such processes, one can define a (partonic) threshold variable ξ, which satisfies ξ → 0 when the heavy particles carry all of the energy in the final state.The precise definition of ξ will depend on the process being JHEP11(2019)002 considered: in Drell-Yan production, for example, it takes the form ξ = 1 − Q 2 /ŝ, where Q 2 is the invariant mass of the produced vector boson, and ŝ the partonic centre of mass energy; more generally, ξ will have the form of a dimensionless ratio of kinematic invariants.One may then write a general schematic form for partonic production cross-sections near threshold, as nm log m ξ ξ + + c (δ) n δ(ξ) + c (0) nm log m ξ + O(ξ) . (1.1) Here we denote by σ 0 the Born-level cross section, which may contain additional coupling factors.The first contribution in the square brackets consists of a series of terms, at fixed order in α s , containing powers of the logarithm of the threshold variable, divided by ξ itself.These contributions can be directly traced to soft and collinear singularities of the underlying scattering amplitudes: the cancellation of infrared divergences between virtual corrections and real radiation leaves behind potentially large corrections, which are still singular as ξ → 0, but are regularised by the well-known plus prescription, so that they are integrable; as discussed below, the all-order structure of these terms is well understood.
The second set of terms in eq. ( 1.1) has support localised on the threshold, ξ = 0, and for processes with electroweak final states it is known that such terms can be formally exponentiated (see, for example, ref. [2]).The third set of terms in the square brackets is suppressed by a single power of ξ with respect to the leading-power contribution.These terms are still singular as ξ → 0: while the singularities they generate are integrable, they can still give numerically sizeable contributions in the threshold region.These next-toleading power (NLP) terms are the focus of the present work, while we will neglect all further sub-leading contributions to eq. (1.1), which vanish at threshold.Order by order in perturbation theory, one can distinguish two expansions in eq.(1.1).Firstly, there is an expansion in powers of the threshold variable ξ, in which we can distinguish the plus distributions and delta function terms as being leading power (LP) in ξ, while the remaining logarithms are next-to-leading-power (NLP).Secondly, for each fixed power of ξ, we can consider the expansion in powers of the logarithm, labelling terms proportional to log 2n−1 ξ as leading logarithmic (LL), the next-highest power as next-toleading logarithmic (NLL), and so on. 1 The problematic nature of LP terms was noted already in the early days of QCD (see for example [3]), and it was quickly realised that such terms at LL level could be summed up to all orders in perturbation theory to achieve a well-behaved result as ξ → 0 [4,5].This resummation was subsequently extended to subleading logarithmic accuracy using two equivalent approaches [6][7][8], themselves partially reliant on earlier diagrammatic arguments for the exponentiation of soft behaviour [9][10][11].Since that time, LP threshold resummation has been reinterpreted and clarified using a wide variety of methods, including the use of Wilson lines [12,13], the renormalisation group [14], the connection to factorisation theorems [15], and soft collinear effective theory JHEP11(2019)002 (SCET) [16][17][18][19].The state of the art for resummation at LP is NNLL accuracy in many processes, including some cases of differential distributions.Recent, pedagogical reviews may be found in refs.[20][21][22].
The phenomenological success of LP resummation, together with the increasing precision of contemporary collider data, makes it natural to ponder whether NLP terms in the threshold expansion can also be classified and resummed, particularly since they have been shown to be numerically significant, for example in the case of Higgs boson production [23,24].Indeed, the study of such contributions has a long history.Subleading corrections involving soft momenta were first investigated in the classic works of refs.[25,26], which dealt exclusively with massive particles in QED.The analysis of ref. [27] updated this to include massless particles.Some years later, the topic was investigated using pathintegral methods in ref. [28], which derived a set of effective Feynman rules for the emission of gauge bosons at next-to-soft level, and argued that a large class of NLP contributions exponentiates.The results were subsequently confirmed by an all-order analysis of Feynman diagrams [29], but concerned massive partons only.In a different approach, NLP effects in certain processes were argued to be resummable, based on well-motivated physical assumptions [30][31][32][33][34] (see also refs.[35][36][37][38][39] for other work related to elucidating all-order properties).
More recently, there has been a revival of interest in studying NLP effects at amplitude level, partly motivated by more formal work relating soft radiation to asymptotic symmetries of the S-matrix in gauge theories and in gravity [40,41].Thus, in addition to the phenomenological applications mentioned above, the study of subleading threshold effects in quantum field theory can have a role to play in finding new representations of, and relations between, gauge and gravity theories [42][43][44][45][46][47], whilst also finding applications in transplanckian scattering [48][49][50][51].In the latter context (as potentially in gauge theories), resummation plays a key role.
In QCD (and related gauge theories), threshold resummation at leading power is known to be a consequence of the universal factorisation of soft and collinear divergences in scattering amplitudes (see for example ref. [15] for a dedicated discussion of this point).This has motivated attempts to construct a factorisation formula for NLP effects.References [52,53] use a diagrammatic approach, building on the earlier work of ref. [27], to describe the effect of dressing a general non-radiative amplitude with an additional gluon emission up to NLO, and NLP in the threshold expansion.This formula contains universal functions similar to those appearing at LP level, but including extra contributions that describe, for example, the emission of wide-angle soft gluons from within jets.A more complete analysis for scalar theories coupled to electromagnetism was undertaken in refs.[54][55][56], which again stress the importance of new quantities (both universal and non-universal) that appear beyond LP order in emitted gluon momentum.Related analyses of NLP effects, envisaging a wide range of applications, both within and beyond a threshold resummation context, have been carried out in SCET [57][58][59][60][61][62][63] (see ref. [64] for earlier work in the context of flavour physics), and results using either diagrammatic or effective theory methods have been shown to be potentially useful for improving the accuracy of fixed-order calculations [59,[65][66][67][68][69][70][71][72][73][74].Recently, the SCET framework has been used to demonstrate that the leading-logarithmic JHEP11(2019)002 (LL) NLP contributions can be resummed, first for event shapes [75], and then for Drell-Yan production [76], where the results agree with the predictions of the physical evolution kernel approach of refs.[30][31][32][33][34].
Our aim in this paper is to show how a similar resummation of LL NLP effects can be achieved using the diagrammatic approach developed in refs.[28,29], and itself analogous to the original LP resummations of refs.[6,7,[9][10][11].As in the SCET approach of ref. [76] (and as observed in refs.[52][53][54][55][56]), we will see that, while it is true that a number of new functions appear at NLP level in the threshold expansion, many of them are irrelevant for discussing the highest power of the NLP logarithm at any given order in perturbation theory.Thus, the resummation of LL NLP contributions is remarkably straightforward.Importantly, this method is sufficiently simple and universal that it can be directly applied to any hadronic cross section with colour-singlet final states: indeed, we explicitly discuss applications to Higgs boson production in the gluon fusion channel, and the formalism can readily be generalised to multi-boson final states.There are a number of motivations for the present analysis.First, they are a natural application of the programme of work commenced in refs.[28,29], where it is was shown that a broad subclass of NLP effects indeed exponentiates.Second, the history of LP resummation suggests that it is highly useful to have more than one formalism for describing equivalent physics: one may note, for example, that comparing logarithmic accuracies between different approaches is non trivial (see for example refs.[77][78][79]), approaches to the choice of factorisation scales may differ, and resummations can be performed in different kinematic spaces [16].Comparing different approaches may therefore provide interesting estimates of the theoretical uncertainty of physical predictions.Third, our diagrammatic approach will provide an alternative starting point for generalising the NLP resummation formalism beyond leading-logarithmic accuracy.
The structure of our paper is as follows.In section 2, we review the resummation of LP threshold contributions, introducing notation that will be useful for what follows.In particular, we will relate our calculation to the path-integral methods of ref. [28], which provide a particularly elegant proof of exponentiation.In section 3, we show how the picture can be naturally extended to NLP level, using existing results.We will argue in detail that potential additional contributions to NLP behaviour, including hard collinear effects, nonuniversal behaviour and phase-space correlations between gluons, can be ignored at LL. Armed with this knowledge, we will then perform an explicit calculation that resums the LL NLP terms in Drell-Yan, comparing our results with others in the literature [30][31][32][33][34].We will then comment on the general applicability of our framework to the production of an arbitrary number of colour singlet particles, before examining Higgs production in the large top mass limit as a further example.In section 3.4, we briefly compare our framework with the recent analysis of ref. [76], in the framework of Soft Collinear Effective Theory (SCET).Finally, we discuss our results in section 4 before concluding.Technical details are collected in three appendices.

Threshold resummation at leading power
In this section, we review the resummation of terms at leading power in the threshold variable, using factorisation methods.Given that our aim in what follows is to sum leadinglogarithmic terms only at NLP, we will mostly concern ourselves here with LL terms also at LP. Furthermore, we will phrase our discussion in terms of methods and notation that allow a straightforward generalisation to subleading power in the threshold expansion.While our discussion applies to general colour-singlet final states, we will first explicitly consider the Drell-Yan production of a massive (or off-shell) vector boson, which at LO corresponds to the partonic scattering process depicted in figure 1.We will not explicitly consider here the quark-gluon production channel, where NLP logarithms are present, but constitute in fact the leading power, since LP logarithms are absent; in the gluon-gluon channel for the Drell-Yan process, only NNLP logarithms can arise.We write the invariant mass distribution in the q q channel as where we restrict ourselves to a single quark flavour for simplicity.Here σ 0 (Q2 ) is the LO total partonic cross section, whose precise value will depend on the nature of the vector boson.Furthermore, α s (µ 2 R ) is the strong coupling at the renormalisation scale µ R , q(x, µ 2 F ) is a quark distribution function with longitudinal momentum fraction x and factorisation scale µ F , while q is the equivalent for an antiquark.Given that scale choice effects contribute to only subleading logarithms (see for example [7]), we will simply choose µ F = µ R = Q from now on, and simplify notation accordingly.In eq. ( 2.2) we defined the invariants where ŝ = (p 1 + p 2 ) 2 is the squared partonic centre of mass energy, and s = ŝ/x 1 x 2 is the hadronic centre of mass energy.The ratio z represents the fraction of ŝ carried by the final state vector boson.At LO this must be unity, so that one has 2 (2.4)

JHEP11(2019)002
The invariant mass distribution in eq.(2.2) is a convolution in z, and can be diagonalised by taking Mellin moments with respect to τ , with the result where is the transformed quark distribution (and similarly for the antiquark), and we have defined Note that, here and in the following, we use a notation in which functions and their Mellin transforms are distinguished only by their arguments.Beyond LO, eq.(2.7) receives potentially large threshold corrections.In Mellin space, these appear as contributions of the form which in momentum space are associated with plus distributions of the form defined such that (2.10) When computed in perturbation theory from quark scattering, ∆ N, Q 2 is affected by collinear divergences, which must be reabsorbed in the quark distributions: below, we will mostly work with the 'bare' ∆, before renormalisation of the coupling α s , and before the factorisation of collinear divergences, which will be regulated using dimensional regularisation in d = 4 − 2 .For clarity, we will denote this bare partonic cross section with ∆(z, Q 2 , ) in momentum space, and with ∆(N, Q 2 , ) in Mellin space.Collinear factorisation is understood to be performed in the MS scheme.
For any QCD process with a colour-singlet final state produced near threshold, ∆ has a factorised structure, and can be written as [6,80] where H Q 2 is an amplitude-level finite hard function containing off-shell virtual contributions, S(N, Q 2 , ) is a soft function collecting all soft enhancements associated with (real or virtual) soft radiation, and ψ i (N, Q 2 , ) is a perturbative (anti-)quark distribution function, collecting collinear singularities associated with initial line i; finally, given that infrared enhancements of both soft and collinear origin are included twice (both in the soft JHEP11(2019)002 and quark distribution functions), one may remove the double counting by dividing each quark distribution by its own eikonal approximation ψ eik,i (N, Q 2 , ).Formal definitions of the (eikonal) quark distributions and of the soft function are given, for example, in ref. [80]: sometimes, eikonal quark distributions are absorbed into the soft function to build the socalled reduced soft function, organising wide-angle soft radiation.On the other hand, one may consider the factor for each initial parton line: this has the effect of removing the soft behaviour from each quark distribution, leaving hard collinear behaviour only.This arrangement is particularly convenient if one wishes to focus only on leading logarithms, as we do in this paper: indeed, at any fixed order in α s , leading logarithms at leading power arise only when the maximum number of singular integrations is performed, yielding the highest inverse power of .Thus, the factor ψ h,i (N, Q 2 , ) for each external line contributes only at subleading logarithmic accuracy, and can be put equal to unity at LL.We are then left with the simple result implying that leading logarithms in the DY cross-section at arbitrary orders in perturbation theory are governed purely by the soft function [6,7,13], on which we now focus.
For any QCD process with a colour-singlet final state, the soft function has a formal definition as a vacuum expectation value of Wilson line operators associated with the colliding partons.Defining the dimensionless four-vectors β i via one may write the soft function (in momentum space) as Here the trace is over colour indices, and the Wilson line operators are defined by where T a is a colour generator in the fundamental representation; furthermore, eq.(2.15) includes a sum over final states containing n partons generated by the Wilson lines, including the appropriate phase space integration, and subject to the constraint that the total energy radiated in the final state equals (1 − z)ŝ; finally, the division by the number of colours N c corrects for the fact that this factor has already been included in the LO cross-section σ 0 in eq.(2.13).Introducing the momentum space gauge field Ãµ (k), one may write the Wilson line exponent as

JHEP11(2019)002
where the square-bracketed factor on the right constitutes the momentum-space factor associated to the emission of a gluon from the Wilson line.We recognise this as the well-known eikonal Feynman rule for soft gluon emission, so that finding the soft function amounts to calculating the cross section for the incoming partons in the eikonal approximation.This cross-section is known to exponentiate, which relies on two properties: first, vacuum expectation values of Wilson lines exponentiate before any phase space integrations are carried out, which may be shown diagrammatically [9][10][11], or using renormalisation group arguments, themselves relying on the multiplicative renormalisability of Wilson line operators [81][82][83][84][85][86]; second, the phase space for the emission of n soft partons factorises into n decoupled one-parton phase space integrals, given that momentum conservation can be ignored at leading power in the threshold expansion.Combining these two properties, one finds that the complete soft function, at crosssection level, has an exponential form, and the exponent can be directly computed in terms of a special class of Feynman diagrams known as webs [9][10][11].These results have been reinterpreted more recently using a path integral approach [28], which incorporated statistical physics methods (the replica trick ) to provide a particularly streamlined proof of diagrammatic exponentiation.These methods have in turn allowed the web language to be generalised to multiparton scattering [87][88][89][90][91][92][93][94][95] (see also [96,97], or ref. [98] for a pedagogical review).We review the replica trick here in appendix A, given that it can also be used to demonstrate directly the exponentiation of a large class of contributions at NLP in the threshold expansion.
Concentrating on leading logarithms, it is important to note that the pattern of exponentiation of soft and collinear singularities is non-trivial, in that the exponent is singlelogarithmic (containing terms of the form α n s log m N with m ≤ n + 1), while the cross section is double-logarithmic, as noted in eq.(2.8).The leading logarithms for the cross sections are therefore completely determined by a one-loop evaluation, which we briefly review below.The eikonal cross-section, up to NLO and in momentum space, can be written as: The real radiation contribution can be obtained from the graphs of figure 2 using eikonal Feynman rules, and one finds The momentum integral is easily evaluated introducing the Sudakov decomposition where k µ T is a four-vector transverse to β µ 1 and β µ 2 , Real emission diagrams for the eikonal cross-section of eq.(2.19), where all emission vertices are assumed to be eikonal.

Contracting eq. (2.20) with β µ
1 and β µ 2 , it is straightforward to verify that the components k ± are given by furthermore, the integration measure in eq.(2.19) becomes where dΩ m is the element of solid angle in m spatial dimensions.Eq. (2.19) then becomes The remaining integrals can be easily carried out using the reparameterisation which expresses eq.(2.24) as Taking into account also that where μ is the MS renormalisation scale, μ2 = 4π e −γ E µ 2 , and the Q dependence is easily restored by recalling that ŝ = Q 2 /z.In experiments one measures the Drell-Yan cross section at fixed Q, thus for z → 1 we expand the cross section around Q, and eq.(2.28) becomes4

JHEP11(2019)002
Given that we focus on leading logarithmic corrections, in the following we will drop the dependent term in the square brackets.
The virtual contribution at O(α s ) can be obtained by direct calculation, or by imposing the soft gluon unitarity requirement reflecting the requirement that soft divergences from the virtual and real contributions must cancel, and the fact that Wilson line correlators are pure counterterms in dimensional regularisation.This requirement implies which applied to eq. (2.29) immediately gives . (2.32) Summing real and virtual correction as in eq. ( 2.18) we obtain the eikonal cross-section at O(α s ): At this point, we can safely take the Mellin transform which gives Expanding in one finds where ψ (n−1) denotes the n-th derivative of the logarithm of the Γ function.Keeping the dominant logarithmic behaviour as N → ∞ one finds the simple result -10 -

JHEP11(2019)002
where we kept the leading power of the logarithm separately for the divergent and for the finite contributions. 5As discussed above, we may exponentiate this result to obtain the leading logarithmic behaviour at all orders.Upon doing so, we may absorb the resulting collinear poles into the parton distributions, using the MS scheme.This amounts to defining renormalised and resummed quark distributions via and similarly for the antiquark, so that eq.(2.5) becomes This formula explicitly sums up leading logarithms in N to all orders.It can easily be verified that eq. ( 2.39) reproduces the well-known results of earlier studies, see for example [6,7,31], both in Mellin space and in momentum space.We note in passing that in our analysis that the dimensional regularisation scale µ appears only through the factor µ 2 , as must be the case on dimensional grounds.Given that µ is then identified with the renormalisation and factorisation scales, it follows that logarithms of these scales (which may be chosen to depend on z) must be suppressed by a single power of , and thus do not contribute to the leading logarithmic behaviour in the threshold variable (1 − z), as could be expected.The same argument will hold at NLP level.We also note that, in going beyond LP level, we will have to keep track of subleading terms in eq.(2.36).Expanding this to NLP order one finds which will be useful later on.
3 Threshold resummation at next-to-leading power In the previous section, we reviewed the exponentiation of leading logarithmic threshold contributions to the Drell-Yan cross-section at leading power.We now discuss how to extend this procedure to next-to-leading power, and we will keep our remarks general enough to apply to both quark and gluon-initiated processes, and for general colour-singlet final states.Recall that LP resummation at LL accuracy relied on two facts: the exponentiation of the soft function before integration over phase space (at squared matrix element level), and the factorisation of phase space for m parton emissions into m decoupled singleparton phase space integrals.This motivates the following schematic decomposition of the partonic cross-section up to NLP order, which was already shown to be useful in ref. [29]: JHEP11( 2019)002 The first term on the right-hand side of eq.(3.1) gives the leading-power result of section 2, integrating the leading-power squared matrix element with leading-power phase space, i.e. neglecting correlations between radiated partons.The second term consists of the NLP contribution to the squared matrix element, integrated with LP phase space.The third term consists of the LP matrix element, but where the phase space includes the effect of parton correlations at NLP.Finally, the ellipsis denotes terms which are NNLP and beyond in the threshold expansion.Based on this classification, the task of determining whether LL NLP terms can be resummed amounts to elucidating the relevant structure of the NLP matrix element, as well as considering whether NLP corrections to the LP phase space are significant.Let us consider each of these issues in turn.

Structure of the NLP squared matrix element at LL
For simplicity, let us first describe the structure of squared matrix elements at NLP level when the hard emitters are massive, following refs.[28,29] (themselves building on refs.[25,26]).In figure 3(a), we draw a non-radiative amplitude with an incoming quark and antiquark, which interact via a hard interaction H. Radiation can then be divided into two types of contribution: 1. External emissions.In this case, radiation couples directly to the incoming hard lines, as exemplified in figure 3(b).Notice that this case includes all radiation that does not resolve the structure of the hard interaction, and incorporates intricate diagrammatic cancellations that lead to a factorised form of the amplitude: indeed, all radiation at leading power falls in this category.
2. Internal emissions.At next-to-leading power, non-factorisable contributions from next-to-soft partons arise, which can be depicted as originating from inside the hard interaction, as in figure 3(c).This corresponds to the insertion of sub-leading power operators in an effective field theory language, and it is the first level of interaction where soft radiation begins to unravel the structure of the hard scattering.
As shown for the first time in ref. [28], external emissions can be described by generalised Wilson lines, which extend the definition given in eq.(2.16) to next-to-leading power in the soft expansion.Along the lines of eq. ( 2.17), we may write this operator in momentum JHEP11(2019)002 space as [53] F for a generalised semi-infinite straight Wilson in the direction of four-momentum p.Here T a is a colour generator in the appropriate representation, and S µν is the generator of Lorentz transformations for the parton under consideration, vanishing for scalar fields, while it is given by for spin-1/2 fields, and for vector fields.The first term in eq. ( 3.2) corresponds to the eikonal Feynman rule of eq.(2.17), and the remaining terms (suppressed by one power of the gluon momentum k) correspond to effective next-to-eikonal Feynman rules, describing the emission of next-tosoft gauge bosons [28].The ellipsis in eq.(3.2) refers to terms involving two gluons being emitted from the same point along the Wilson line, through seagull-type vertices.These vertices start contributing to the cross section at NNLO (either through double radiation, or through one-loop corrections to single radiation), therefore they cannot contribute at leading logarithmic accuracy at NLP (as was the case at LP).Indeed, our proposed resummation rests upon an amplitude-level factorisation theorem [27,52,53], which implies the existence of evolution equations, which in turn can be understood in terms of renormalisation of suitable operator matrix elements: the solution of evolution equations of this type always leads to a non-trivial pattern of exponentiation, with single logarithms in the exponent generating double logarithms in the cross section.Such a pattern of exponentiation implies that all leading logarithms are generated by the one-loop exponent.A test of this argument is provided by ref. [36], where the exponentiation of leading logarithms at NLP was explictly tested at NNLO; finally, as a further check, we verify in appendix C that next-to-soft Feynman rules for double radiation in eq.(3.2) do no contribute to leading logarithms in the case of Drell-Yan production at two loops.
Let us now consider the contribution of internal emissions.When massive external particles are being considered, the hard interaction is analytic in the total momentum K of the emitted radiation, and can safely be expanded about the soft limit K µ → 0. One may then show, using Ward identities, that the effect of a single internal emission is given by derivatives of the non-radiative amplitude with respect to its external momenta.As has been noted in the context of the so-called next-to-soft theorems of refs.[40,41], 6 these derivatives can be organised in terms of the orbital angular momentum operator associated with each external leg, which, in momentum space, has the form

JHEP11(2019)002
On each hard leg, this combines with the spin angular momentum contribution to construct the total angular momentum operator In ref. [28], the orbital angular momentum contribution was not included in the generalised Wilson line operator of eq.(3.2), despite the fact that it might make sense to do so, given that the internal and external emission contributions are not separately gauge-invariant, but instead combine into a gauge-invariant object, the total angular momentum.For practical purposes, however, it remains convenient to keep the orbital angular momentum separate, given that it involves derivatives which have yet to act on the hard interaction.How to keep track of such contributions will be discussed explicitly below.
Armed with the operator defined in eq. ( 3.2), we may construct a next-to-soft function by analogy with the LP soft function of eq.(2.15), as Here we have replaced the Wilson line operators in the LP soft function by their NLP counterparts, and the subscript in the sum over final states indicates that all phase space integrals are to be carried out with LP phase space only (i.e. with a measure of integration consisting of a product of single-gluon phase space integrals).Corrections to this will be considered in section 3.2.Note also that all generalised Wilson lines are semi-infinite straight lines proceeding from the origin in position space.At NLP accuracy the cross section is sensitive to a potential non-zero initial position, but this is related to the derivative contributions above [28], which are to be dealt with separately.As was the case at LP, the next-to-soft function in eq.(3.7) can be shown to exponentiate using replica trick arguments (see appendix A).At NLP, however, we must carefully disentangle what this means, given that the generalised Wilson line of eq.(3.2) is matrix-valued in the spin space of the external hard particles.As an example, consider the spin-1/2 case, and let us write the non-radiative amplitude with an incoming fermion and antifermion with explicit spin indices {α, β}, as so that the spin matrix M is defined by stripping off the initial state wave functions from the full amplitude.The next-to-soft function of eq.(3.7) can then be explicitly written as a spin operator (3.9) where the ordering of spinor indices is depicted in figure 4. The spin matrix M αβ in eq.(3.8) factorises into a product of hard and next-to-soft factors, so that the integrated squared matrix element, dressed by arbitrary amounts of radiation from the next-to-soft function, -14 -JHEP11(2019)002 can be written as ) where dΦ (m) denotes the m-particle Lorentz-invariant phase space measure, and the integration over the phase space of the heavy vector boson has been singled out, relying upon the factorisation of the n-body phase space at LP.Note that eq.(3.10) also applies to the case of incoming particles of spin one, if the spinor wave functions are replaced with polarisation vectors, and spinor indices by vector indices.
The discussion so far applies strictly only to the case of massive external particles.When massless particles are involved, it is no longer true that the hard function H is analytic in the momentum carried by soft radiation: it develops logarithmic singularities due to the presence of collinear divergences.As discussed in the Introduction, this was first explored in a QED context in ref. [27], which presented a factorisation formula at amplitude level, extending the Low-Burnett-Kroll theorem to include collinear effects.Similar ideas have recently been extended to full QCD [52,53], and analysed using SCET [57-62, 68-72, 75, 76, 101-103], while an alternative, first-principles, diagrammatic approach has been explored in ref. [54].Common to all these approaches is the factorisation of collinear contributions into universal functions, sensitive to the spin of the colliding particles but otherwise independent of the details of the hard scattering.More precisely, one may recall that, at leading power, collinear radiation is accounted for by means of jet-type functions, such as the parton distributions introduced in eq.(2.11).At NLP, one must generalise this analysis, expressing radiative amplitudes in terms of new types of jet functions, describing soft emissions from collinearly enhanced configurations.The first such radiative jet function was proposed in ref. [27], and was recently calculated at one loop JHEP11(2019)002 in QCD for external quarks in refs.[52,53], where it was used to reproduce known NLP threshold logarithms in Drell-Yan production at NNLO.The amplitude-level factorisation proposed in [52] is expected to apply only to annihilation processes involving two colliding hard partons producing colourless final states: on the other hand, the analysis of ref. [54] (which focuses on scalar theories, but considers more general scattering processes), and the results of refs.[57,76,101,102] suggest that further types of jet emission function are necessary in QCD, which have yet to be calculated.
Importantly, for the purposes of the present paper, radiative jet functions can be ignored, since it can be shown that leading logarithms at NLP can only arise from momentum regions of integration that are already fully accounted for by the next-to-soft function introduced in eq.(3.2).To illustrate this point, consider first, for comparison, the wellunderstood situation at leading power. 7In that case, threshold singularities, inducing nonanalytic behaviour at z → 1, are directly related to infrared singularities of the amplitude; these, in turn, arise from integrations of the relevant momentum components ('normal variables') near singular surfaces in momentum space, which can be completely characterised to all orders in perturbation theory by means of Landau equations and power counting techniques [104].For massless theories, it can be shown in general that infrared singularities arise only from soft and collinear momentum configurations.At leading power, therefore, one finds that at n loops there are precisely 2n normal variables that must be integrated with a logarithmic measure: in a suitable frame, these can be taken to be n parton energies E i , with a leading-power integration measure dE i /E i , and n transverse momenta with respect to the directions defined by external particles, k iT , with a leading-power integration measure dk iT /k iT .Threshold logarithms in general arise when different combinations of normal variables become small at different rates, but leading logarithms arise only with a very specific scaling, when all energies and transverse momenta are strongly ordered,8 say E 1 . . .E n and k 1T . . .k nT .In that limit, at LP, the 2n logarithmic integrations yield contributions of the form ln 2n−1 (1 − z)/(1 − z), since the last logarithmic integration must not performed when computing dσ/dz.At NLP, either the phase-space measure or the squared matrix element provide a single power of one of the normal variables, so that only 2n − 1 momentum components need be integrated with a logarithmic measure.Once again, leading logarithms will arise from the configuration where the remaining normal variables are strongly ordered, with 2n − 1 integrations leading to contributions of the form ln 2n−1 (1 − z), while the z integration will not introduce further singularities at NLP.Now, two possibilities arise.On the one hand, the normal variable whose integration has become non-singular can be a transverse momentum, in which case the corresponding parton is soft, but not strictly collinear: this configuration is accounted for by the leading-power JHEP11(2019)002 soft function, which contains 'wide-angle' soft gluons.On the other hand, the suppressed variable can be an energy: in this case, all transverse momenta must be strongly ordered; such next-to-soft, collinear configurations are accounted for by the next-to-soft function defined in eq.(3.7).Notice that radiative jet functions such as the one computed in [52] also contain the next-to-soft, collinear configuration: this, however, contributes to a double counting that must be explicitly subtracted, either by introducing eikonal jets, as done in eq.(2.11), or by defining an appropriate counterterm, as done for example in ref. [53].The subtracted radiative jet function then contains only hard collinear configurations for all radiated partons, and cannot contribute at leading logarithmic accuracy.
An explicit example and test of the above discussion is provided in ref. [53], where the non-abelian radiative jet function for quarks was computed at one-loop order, and the overlap between (next-to-)soft and collinear emissions was explicitly identified.Furthermore, a large class of (N)LP threshold effects has been calculated in Drell-Yan production at NNLO [66] and N 3 LO [67] using the method of regions [109][110][111], which allows for a precise identification of the (next-to)soft and/or collinear origin of all contributions to the cross section.The role of hard collinear effects is indeed found in these studies to be associated with NLL terms and beyond,9 while all LL contributions can be traced to the (next-to)soft function, if the results are recast in the present framework.Notice that, as discussed above, upon exponentiation leading logarithms at NLP must be generated by one-loop contributions: the results of [53,66,67] therefore provide a complete test of our argument for the Drell-Yan process.
To summarise, NLP contributions to squared matrix elements can be categorised into two main types, as follows.
(i) (Next-to-)soft emissions.These are captured by the next-to-soft function, defined in terms of generalised Wilson lines in eq.(3.7), together with the orbital angular momentum contributions associated with internal emissions in figure 3.As at LP, the next-to-soft function exponentiates (see appendix A).
(ii) Collinear contributions.These are described by radiative jet functions, which overlap with the next-to-soft function.Upon removing the double counting, the remaining collinear effects do not contribute at LL accuracy.
In this section, we have discussed the second term on the right-hand side of eq.(3.1), and argued that the next-to-soft function underlies all contributions to the NLP matrix element that can result in LL terms in the cross-section.This is only part of the story: we must also check whether or not LL terms can arise from the LP matrix element, once correlations between radiated gluons (a NLP effect) are included.This is the subject of the following section.

Corrections to the LP phase space
The third term in eq.(3.1) consists of the LP matrix element integrated over the NLP phase space.To see whether or not this term can give LL contributions at NLP, it is JHEP11(2019)002 sufficient to take the LL contribution to the LP matrix element at each order, and then to evaluate the phase space integral up to NLP order.The LL contributions to the matrix element have already been discussed in section 2, and involve exponentiating the NLO eikonal squared matrix element.This generates terms with n ≥ 1 gluon emissions, and, according to eq. ( 2.15), one must then integrate each such term over the n-gluon phase space.Considering all possible contributions to an n-gluon final state yields a squared matrix element of the form where the prefactor f (α s , , µ 2 ) collects coupling dependence, possible poles in due to the integration over loop momenta, and combinatorial factors from the exponentiation of the squared matrix element.The explicit form of this function is irrelevant for what follows.
We must now integrate eq.(3.11) over the (n + 1)-body phase space, consisting of n gluons, as well as the electroweak vector boson that defines the final state at LO.The integration measure is given by where the integration of the vector boson momentum has already been carried out, using the overall momentum conservation δ function.In order to compute the integral, it is particularly convenient to use the Sudakov decomposition of eq.(2.20) for each momentum k i .One finds where we absorbed the θ functions from the factors δ + (k 2 i ) into the integration limits for k i± .In order to proceed, we can represent the δ function in the second line of eq.(3.13) using We may then rewrite eq.(3.13) as where in the second line we Taylor-expanded the term in the exponent that is quadratic in soft momentum, anticipating that higher order contributions in T in the last line will -18 -

JHEP11(2019)002
correspond to subleading powers of (1 − z) in the final result.We will verify this fact later, but, for the moment, note that the term at O(T ) corresponds to a phase space correlation between pairs of gluons that is absent at LP. Thus, this term constitutes the "NLP phase space" correction referred to in eq.(3.1).In the Sudakov decomposition, the dot product of gluon momenta reads The term involving the transverse momenta leads to an odd integrand in each k iT in eq.(3.15), and will therefore give a vanishing contribution to the final result.We can then carry out the remaining transverse momentum integrals to obtain Next, the k i± integrals can be straightforwardly carried out to give The integral in T is recognisable as an inverse Laplace transform, which yields the result Note that the terms O(T 2 ) we have neglected in expanding the exponential factor in eq.(3.15) give subleading power corrections in (1 − z), justifying the approximation made above.Eq. (3.19) is the final result of integrating the LP contribution to the matrix element responsible for LL terms, with the multigluon phase-space measure expanded to NLP order.The second term in the last line of eq.(3.19) is the desired NLP correction, as can be seen by the fact that it is suppressed by a single power of (1 − z).Furthermore, it contains an explicit factor of d − 4 = −2 , which directly implies that the phase space correction does not affect LL terms, which are associated with the most singular poles in .
In summary, we have shown that the third term in eq.(3.1), consisting of the LP matrix element dressed with NLP phase space corrections, does not contribute to LL terms at NLP order.It can thus be neglected for the purposes of this paper.Combining this observation with the results of the previous section, we now have everything we need to perform an explicit resummation of LL NLP threshold logarithms in Drell-Yan production.We turn to this task in the next section.

JHEP11(2019)002
H H (b) Figure 5. Diagrams contributing to a squared amplitude with a q q initial state, arising from the next-to-soft function S acting on the LO hard interaction H, as in eq.(3.10).Further diagrams are obtained by reflection about the final state cut, or by interchanging p 1 ↔ p 2 .

Resummation of leading NLP logs in Drell-Yan production
In the previous sections, we have seen that LL contributions at NLP level are governed by next-to-soft radiation.This in turn is captured by the next-to-soft function defined in eq.(3.7), possibly complemented by contributions involving the orbital angular momentum of each incoming parton.In this section, we apply these ideas to resum LL NLP terms in Drell-Yan production.While clearly very interesting for its own sake, this example is also a useful warm-up case: first, it will allow us to make contact with the LP treatment of section 2; furthermore, in this case the hard interaction is especially simple, so that its derivatives with respect to the external momenta vanish at leading order.Thus, we do not have to worry about orbital angular momentum contributions at LL accuracy, and it is sufficient to calculate the next-to-soft function.Once this has been calculated for single radiation, it may be exponentiated (as at LP), yielding the resummation formula that we are seeking.
For the first steps of our derivation, we do not need to specify the final-state particle content of the process we are studying.Rather, we will consider a general hard interaction H connecting to an incoming q q pair, such that the LO amplitude is given by eq.(3.8).Representative diagrams contributing to the squared amplitude arising from the next-tosoft function at NLO, eq.(3.10), are shown in figure 5. We may directly evaluate them using the Feynman rules arising from eq. (3.2).First, we may note that contributions involving k 2 vanish, since the radiated gluon is on shell.Next, it is convenient to combine the scalar-like and spin-dependent emission vertices as Then, the diagrams of figure 5 yield a NLP contribution where a factor of two for complex conjugate diagrams has already been included.In order to extract the LO squared amplitude, we may use an argument similar to one presented

JHEP11(2019)002
recently in ref. [65].Writing the decomposition (cf.eq.(2.20)) and substituting into eq.(3.21) reveals that the term involving transverse momentum occurs linearly in the squared matrix element, leading to an odd integrand which vanishes upon integrating over k T .This contribution can thus be ignored, leading effectively to the expression Notably, in eq.(3.23) the LO squared matrix element is factored out.Combining this with diagrams obtained from those of figure 5 by interchanging p 1 with p 2 , summing over spins and colours, and dividing out the LO cross section one easily obtains an expression for the real emission contribution to the one-loop next-to-soft function.Notice that at NLP singularities as z → 1 are integrable: thus, there is no need to combine real emission with virtual corrections in order to generate LL contributions, and the NLP soft function reads The integration over the real gluon phase space can be carried out straightforwardly using the Sudakov decomposition of eq.(2.20), and the subsequent change of variables in eq.(2.25), with the result Taking the Mellin transform we find The leading behaviour as N → ∞ is where the ellipsis denotes terms which are non-singular in and non-logarithmic in N , as well as terms suppressed by further powers of N .We see that the NLP soft function generates contributions which are suppressed by (at least) a single power of N compared to LP, as expected.Eq. (3.27) must now be combined with the LP soft function given eq.(2.40), -21 -

JHEP11(2019)002
which itself includes subleading terms in N space arising from the Mellin transformation from z space.The result is As explained above, we may directly exponentiate eq.(3.28), and combine it with the LO cross-section.Furthermore, the collinear pole can be absorbed in the quark distributions, for which we again use the MS factorisation scheme.To this end, we generalise eq.(2.38) to and similarly for the antiquark.Note that it is important in eq.(3.29) that we correctly kept track of subleading terms in the Mellin transform of the LP soft function.The crosssection at NLP in the threshold expansion and at LL accuracy then becomes Upon expanding the exponential factor in powers of α s , we may now perform the inverse Mellin transform of the partonic cross-section order by order, using the results of appendix B, to get This is in complete agreement with (and indeed provides an independent proof of) the result of ref. [31], which argued (consistently with previous observations [23,36]) that the LL NLP terms at any order have a coefficient which is always the negative of that of the corresponding leading logarithmic plus distribution.The origin of this phenomenon can be traced to the coefficient of the pole in eq.(3.27).Given that this pole represents a collinear singularity that must be absorbed in the parton distributions, it must emerge from the NLP contribution to the LO DGLAP splitting kernel that governs such terms.More specifically, the collinear poles in the NLO Drell-Yan cross-section have the form (see for example ref. [112]) where the factor of 2 arises from having collinear singularities associated with either of the incoming partons.The splitting function can be expanded near threshold as .33)where the second term gives the NLP contribution in z-space, whose Mellin transform is

JHEP11(2019)002
We thus expect the collinear pole of the NLP contribution to the next-to-soft function in Mellin space to be given by 2α which is indeed observed in eq.(3.27), and, in momentum space, in eq.(3.25).We see that the next-to-soft function generates the correct NLP correction to the splitting kernel as expected.This in turn dictates the LL behaviour in the finite part: indeed, in z-space, this contribution arises completely from an overall -dependent power of (1 − z), dressing the pole term.Thus, ensuring that the NLP behaviour of the pole term is correct is sufficient to describe also the finite part. 10Note that the fact that all NLP information in the DGLAP splitting function is correctly generated by the next-to-soft expansion provides a test of the statement made earlier, that all LL threshold effects at (N)LP arise from radiation that is (next-to) soft, in addition to being collinear.This in turn confirms that, at LL accuracy, one may neglect radiative jet functions [27,[52][53][54][55][56].
Eq. (3.30) resums the leading-logarithmic behaviour of the Drell-Yan cross section at LP and NLP in the threshold expansion: it completely agrees with expectations from the literature [31,36], and thus with the recent SCET analysis of ref. [76], which cross-checked against the same references.We emphasise that, of course, at leading power there is no need to limit the resummation to leading logarithms: this was done in eq.(3.30)only for simplicity, and to underline the close connection between leading logarithms at LP and NLP.Indeed, because of the link discussed above between NLP leading logarithms and the DGLAP kernels, it is straightforward to incorporate our results in the standard LP resummation formalism: it is sufficient to include NLP terms in the quark splitting function.This was argued to be appropriate in refs.[23,31,36], and, with the mild assumptions discussed in section 3.1, it is now proven.For completeness, we include here the general resummation ansatz introduced in ref. [36], which implements this change in the classic threshold resummation formula of [6][7][8], together with other proposed modifications that have effects on subleading NLP logarithms.In Mellin space, the result of ref. [36] for the Drell-Yan process can be written as In eq.(3.36), D(α s ) is the well-known LP wide-angle soft function for the Drell-Yan process, which has been computed up to three loops in [113][114][115][116]; F DY (α s ) resums N -independent contributions following ref.[2]; P LP+NLP qq (z, α s ) is the soft expansion of the DGLAP splitting function up to NLP, order by order in perturbation theory, which was derived in [36] starting from the results of ref. [117]; furthermore, the 'plus' prescription is defined to JHEP11(2019)002 apply only to LP contributions, that are singular as z → 1.Leading NLP logarithms in eq.(3.36) are generated by the one-loop NLP contribution to P LP+NLP qq , as discussed in this section.Higher-order terms in the NLP splitting function will contribute to, but not exhaust, subleading NLP logarithms; indeed, the shifts in the phase space boundary and in the argument of the coupling, proposed in eq.(3.36), and corresponding to a NLP-accurate definition of the soft scale of the process, also contribute to subleading logarithms at NLP.In ref. [36], the accuracy of eq.(3.36) was tested by comparing its expansion to NNLO with existing exact results: as expected from our current discussion, leading NLP logarithms are exactly predicted; furthermore, one observes that next-to-leading NLP logarithms are predicted very accurately, and they mostly arise from the NLO contribution to the NLP splitting function.The small discrepancy arising at this level of accuracy (NLL at NLP) between the resummation and the finite order result is the first footprint of the need to include radiative jet functions at NLP.

A brief comparison with the SCET approach
In section 3.3, we have achieved the resummation of leading NLP logarithms (jointly with all LP logarithms) by applying essentially diagrammatic arguments, based on the previous analysis of ref. [28], summarised here in appendix A. Clearly, these diagrammatic arguments are in turn based on an underlying factorisation [52,53], but the argument for resummation is greatly simplified by the diagrammatic exponentiation properties of the (next-to-)soft function.Recently, the resummation of these same contributions has been achieved for the Drell-Yan process also within an effective field theory approach based on SCET [76] (see also ref. [75]).In this section, we briefly compare our methods with the SCET analysis, whose physics must ultimately be equivalent.
The SCET approach relies upon a factorisation of the partonic cross section ∆(z), obtained by expanding the Drell-Yan QCD current into operators defined in terms of effective soft and collinear fields, with a different collinear sector associated with each external parton.Hard modes of the field contribute through short-distance coefficients of SCET operators, which can be obtained by matching to full QCD.Under the assumption that Glauber-type modes of the gluon field do not contribute to the relevant observable, 11soft and collinear modes can be factorised into (universal) matrix elements, which define collinear and soft functions, in direct correspondence with the jet and the soft functions emerging from the diagrammatic approach considered in this work.
Restricting to the terms relevant for the resummation of leading logarithms up to NLP, the momentum-space SCET factorisation for the partonic cross-section ∆ introduced in eq.(2.5) takes the form where H(Q 2 , µ) is the hard function, S DY represents the leading-power soft function, equivalent to the one defined in eq.(2.15), and S 2ξ represents that part of the NLP soft function

JHEP11(2019)002
which contributes at leading logarithmic accuracy, to be compared with LL form of S NLP in section 3.3.In eq.(3.37) no collinear functions appear explicitly, since hard collinear modes contribute only starting at NLL accuracy.This conclusion is obtained within SCET by an analysis of all possible operators contributing at NLP, and confirmed a posteriori, as we will see below.At a given factorisation scale µ, as written in eq.(3.37), either the hard or the soft function (or both) develop large logarithms.The effective field theory formulation allows one to evaluate each function at their characteristic momentum scale, thus setting µ = µ h ∼ Q in the hard function, and µ = µ s ∼ Q(1 − z) in the soft function, so that no large logarithms appear.The independence of physical observables on factorisation scales yields then a renormalisation group equation, whose solution allows one to evolve the hard and soft scale to the common scale µ, thus resumming the large logarithms.To this end, in ref. [76] the common scale was chosen to be µ To LL accuracy, the evolved hard and soft functions can be written as [76] Here the evolution factor E LL resums the logarithms, and can be expressed in terms of the strong coupling evaluated at different scales, as where we have introduced the QCD cusp anomalous dimension Γ cusp (α s ) and the beta function β(α s ), with the normalisation Using the fact that, for Drell-Yan production, H(Q 2 , µ h ) = 1 + O(α s ), one can compute ∆(Q 2 , z, µ c ) in eq.(3.37) with all factors evaluated at the scale µ c .In order to avail oneself of the standard collinear factorisation machinery, one must then evolve both ∆(Q 2 , z, µ c ) and the parton distributions to a generic factorisation scale µ, exploiting the RG invariance of the physical cross-section.One then finds The fact that the evolution of parton distributions (dictated by DGLAP splitting functions) is consistent at LL accuracy with the evolution of the partonic cross section to a generic scale provides an independent check of the fact that collinear function (absent in eq.(3.37)) cannot contain leading NLP logarithms.This is directly analogous to the observation made here in section 3.3, where we noted that the effect of including next-to-soft radiation led to -25 -

JHEP11(2019)002
reproducing the NLP contribution to the DGLAP kernels, testing our arguments for not including radiative jet functions in the derivation leading to eq. (3.30).
In order to compare the result in eq.(3.41) with eq.(3.31), one needs to expand the ratios of running couplings in eq.(3.41) and eq.(3.39) in powers of α s (µ).When this is done, the NLP term in eq.(3.41) reduces to (3.42) Upon setting the hard and soft scales to their natural values, µ h = Q and µ s = Q(1 − z), and choosing (as above) a factorisation scale of µ = Q, we find which is equivalent to eq. (3.31).To see this, one must Mellin transform eq.(3.43), expand in α s , and then perform the inverse transform back to momentum space.

Resummation for general quark-initiated colour-singlet production
In section 3.3 we have seen how to resum the highest power of NLP logs, for the specific case of Drell-Yan production.In fact, the result can easily be generalised to the production of N colour singlet particles (which may be loop-induced at leading order), with a q q initial state.Crucial to our arguments will be exponentiation of the next-to-soft function in terms of webs [28,29], which implies that the next-to-soft function has the schematic form where the first sum is over leading power webs composed with eikonal Feynman rules, and the second sum is over next-to-leading power webs, containing eikonal Feynman rules with at most one next-to-eikonal vertex.Next, we may note, as was remarked in refs.[28,29], that if we are only interested in NLP terms in the final result for the cross-section, we do not in fact have to exponentiate the NLP webs: upon expanding eq.(3.44) in powers of the coupling, quadratic and higher powers of the NLP term will give NNLP contributions and beyond.Thus, we may formally replace eq.(3.44) with the equivalent expression (up to NLP level) This expression shows us that, in order to generate a contribution to the highest power of the NLP logarithm at any given order, we must take the leading logarithmic behaviour from the NLP web term, namely the contribution proportional to JHEP11( 2019)002 in Mellin space, and dress this with the leading logarithms coming from the leading-power soft function.Note in particular that the webs W NLP do not contain terms of the type such terms will arise in the cross section only through the expansion of the exponential in eq.(3.44), precisely through the interference between leading-power and next-to-leading power webs.We can see this directly in eq.(3.30) for Drell-Yan production: upon Taylorexpanding in α s , the leading logarithm at NLP comes from a single instance of the leading NLP log at O(α s ), dressed by arbitrary powers of the leading logarithm at leading power.
For arbitrary processes, we must broaden the discussion presented for the Drell-Yan case to include an additional next-to-soft contribution, associated with the orbital angular momentum of incoming particles, which combines with the spin angular momentum present in the next-to-soft function to build a gauge-invariant result.To this end, let us consider the effect of a single emission from the non-radiative amplitude; this has been examined recently in ref. [65], and we will now present a short summary of that discussion, before drawing consequences for the present paper.We label momenta as shown in figure 6(a), and we write the LO non-radiative amplitude for a q q-initiated process as where {p i } are the incoming parton momenta, and H LO is the LO hard function, which coincides with the LO stripped matrix element M (q q) LO .Let us now consider the radiative amplitude with external wave functions removed, which we denote by M (q qg) σ .As shown in ref. [65], this amplitude, up to NLP order, can be decomposed as where the first (second) term on the right-hand side originates from the spin-independent (spin-dependent) part of the next-to-soft function, while the third term corresponds to the orbital angular momentum contribution discussed above.The squared real emission -27 -

JHEP11(2019)002
amplitude, summed over colours and spins, is then given by12 where we have used the gluon polarisation sum since the contribution of unphysical polarisations vanishes when just a single gluon is radiated.The various contributions to eq. (3.49) have been calculated explicitly in ref. [65], and one obtains for squared matrix element, summed over colours and spins, the expression where initial state momenta in the LO matrix element have been shifted according to13 (3.53)In words, the NLP squared amplitude for single real emission (summed over colours and spins) consists of an overall eikonal factor dressing the LO squared amplitude, whose incoming momenta are shifted according to eq. (3.53).These shifts have the effect of rescaling the partonic Mandelstam invariant ŝ according to where the threshold variable z is defined by satisfying the momentum conservation condition Crucially for what follows, all NLP effects in the matrix element are absorbed in the momentum shift, so that the prefactor in eq.(3.52) simply dresses the shifted matrix element with a leading-power soft emission.We may obtain the partonic cross-section for the single real emission contribution by integrating over phase space and including flux and spin/colour averaging factors.The phase space for the (N + 1)-body final state, with momenta labelled as in figure 6(a), may be written in factorised form as

JHEP11(2019)002
namely as the convolution of a two-body phase space for the gluon momentum k and the total momentum P carried by colour singlet particles, with the subsequent decay of the latter into the individual colour singlet momenta {p i }.Parametrising momenta according to (1, 0, . . ., sin χ, cos χ) , Eq. (3.57) becomes [65] dΦ N +1 (P + k; and where dΦ (z) N denotes the phase space for the N colour singlet particles, but with kinematics shifted according to eq. (3.54).Then, the partonic cross-section including a single additional emission (up to NLP level) takes the form14 where while the LO cross-section with shifted kinematics is given by (3.63) To calculate this in practice, one must implement the momentum shifts in the squared amplitude, and then integrate over the real-emission phase space, whose momentum conservation condition includes NLP corrections: see ref. [65] for a full discussion.Note that we have changed the flux factor on the right-hand side so as to give the full shifted crosssection on the left-hand side.As discussed above, the generalisation of eq.(3.61) to all orders is obtained by dressing the single-emission cross-section with a further arbitrary number of leading-power soft gluon emissions.In eq.(3.52), this has the effect of replacing the prefactor -whose form is obtained from the soft function at O(α s ) -with that obtained from the all-order leading-power soft function.Furthermore, the (N + m)-body phase space for the emission of N colour singlet particles and m additional gluons, with momenta {p i } and {k j } respectively, factorises as in eq.(3.59), and one may write JHEP11(2019)002 so that eq.(3.61) can be straightforwardly replaced with ∆ (q q) NLP (z, ) = z S LP (z, ) ∆ (q q) LO zŝ , (3.65) where the factor of z on the right-hand side originates from having shifted the flux factor in eq.(3.63).In eq.(3.65), S LP (z, ) is the leading-power soft function, defined to include integration over the soft gluon phase space, as in eq.(3.7): it will contain residual collinear poles in that must be absorbed into the quark distribution functions, as was done in eq.(3.29).We may then resum leading-logarithmic LP and NLP terms in the partonic crosssection as follows.First, we notice that the leading-order partonic cross section with shifted kinematics becomes a function of zŝ = Q 2 , which is the physically measured invariant mass that must be kept fixed: we can therefore treat the factor ∆(qq) LO as independent of z, writing ∆ (q q) NLP (z, ) = z S LP (z, ) ∆ Taking the Mellin transform we find Since the leading-power soft function is insensitive to the details of the hard process, we can directly use eq.(2.40) for the soft factor.Removing collinear poles, using exponentiation, and keeping track of NLP terms that arise from the Mellin transform of the LP soft function, we find This simple result resums leading logarithmic terms in Mellin space at both LP and NLP, in the partonic cross-section, for a general quark-induced colour singlet production process.
In the Drell-Yan case, it agrees with eq.(3.30), 15 thus providing an important cross-check of eq.(3.68).As was the case for the Drell-Yan process, eq.(3.68) can be generalised to include the complete known result for the resummation of leading-power subleading logarithms, yielding an expression identical to eq. (3.36) for the resummed partonic cross section ∆(N, Q 2 ).Indeed, the orbital angular momentum contribution that was trivial for the Drell-Yan cross section, due to the point-like nature of the Born process, will result in a shift of the center-of-mass energy ŝ, which must be applied to the Born cross section, with consequences that will depend on the particular process and observable being considered (it would for example be non-trivial for loop-induced processes).The Sudakov exponent, on the other hand, will be unaffected, so that eq.(3.36) will still apply.
In this section, we have seen that resummation of LL terms is possible at both LP and NLP for the general production of N colour-singlet particles, in the q q channel.Similar arguments may be made for gluon-initiated processes, as we discuss in the following section.

Resummation for general gluon-initiated colour-singlet production
In section 3.5 we considered the production of a generic colour singlet final state in quarkantiquark scattering.A similar analysis can be made for gluon-initiated processes: one may obtain leading logarithmic NLP contributions by combining the next-to-soft function with orbital angular momentum contributions.As for the quark case of section 3.5, we can then dress the effect of a single gluon emission at NLP with an arbitrary number of leading-power soft gluon emissions.The case of single emission has been studied alongside the quark case in ref. [65], leading to a result identical in form to eq. (3.52) for the squared amplitude.Indeed one finds As in the quark case, this takes the form of the LO non-radiative transition probability, with kinematics shifted according to eq. (3.53), dressed by a single leading-power soft emission, whose colour factor in this case reflects the emission from an initial-state gluon rather than an initial-state (anti)-quark.The factorisation of phase space will be identical to the previous section, given that this is independent of the particle species.One then obtains the resummed result where the soft function on the right-hand side is defined in terms of Wilson lines in the adjoint representation.One may then follow similar arguments to those leading to eq. (3.68), yielding A check of this result is that it reproduces known LP, and conjectured NLP results for Higgs boson production, in the large top mass limit.As is well-known, the LO process consists of an effective coupling between the Higgs boson and a pair of gluons, as shown in figure 7. Higher-order contributions near threshold have been discussed for example in -31 -JHEP11(2019)002 ref. [33], which expressed the hadronic cross section for the gg channel as Here g(x i , µ 2 ) is the gluon distribution, we have set the factorisation and the renormalisation scales to the common value µ, and c gg a perturbative coefficient function.Furthermore, we have introduced the quantities which normalises the LO cross-section, and where v is the Higgs field vacuum expectation value.With the normalisation adopted in eq. ( 2.2) we have σ 0 = m 2 H σ0 , and where collinear poles in S LP (z, ) have already been factorised into the gluon distributions, leaving a finite remainder S LP,fin.(z).Eq. (3.74) is valid in general, i.e. also away from the infinite top mass approximation for the ggH vertex.In this case one has where F m 2 H /(4m 2 t ) = F m 2 H /(4m 2 t ), | =0 is defined e.g. in eq. ( 5.3) of [65].In the infinite top mass approximation F m 2 H /(4m 2 t ) → 1, and we can easily identify This in turn implies that the coefficient of the LL NLP term in c gg at a given order in α s is related to the LL LP term by a minus sign.Furthermore, both sets of terms are related to their counterparts in Drell-Yan production by the simple replacement C F → C A , given that the LP soft functions in both cases obey 'Casimir scaling' to the relevant order.We thus reproduce the results of ref. [33] for the resummation of LL NLP logarithms in single Higgs production in the large top mass limit.We stress, however, that the main result of this section (eq.(3.71)) is more general: it applies also away from the large top mass limit, and for other gluon-induced processes.Its application to specific collider processes of interest will be the subject of future study.Also for gluon-initiated processes, we note that subleading LP logarithms can be included, and the result will take the general form of eq.(3.36): in this case, the gluon DGLAP splitting functions will be involved, while the soft function for gluon annihilation can be obtained from the quark case by Casimir scaling, at least up to three loops.

JHEP11(2019)002 4 Conclusion
In this paper, we have developed a formalism for resumming leading-logarithmic (LL) threshold contributions to perturbative hadronic cross-sections, at next-to-leading power (NLP) in the threshold variable.This generalises previous approaches at leading power (see for example refs.[6][7][8][12][13][14][15][16][17][18][19]), and applies to the production of an arbitrary coloursinglet final state at LO.Our method builds upon the previous work of refs.[28,29] (and subsequent studies [52,53]), which describes leading NLP effects in terms of a next-to-soft function, which can be shown to exponentiate at the diagram level, so that the logarithm of the next-to-soft function can be directly expressed in terms of Feynman diagrams dubbed next-to-soft webs.In general processes, the next-to-soft function must then be supplemented by terms involving derivatives acting on the non-radiative amplitude, which can be interpreted in terms of the orbital angular momentum of the colliding partons.Leadinglogarithmic accuracy can then be achieved by dressing the effect of a single emission, computed up to NLP level, with the LP soft function.In this sense, our results provide a non-trivial generalisation of the so-called next-to-soft theorems [40,41], which have recently been intensively studied in a more formal context, for both gauge theories and gravity.
We have explicitly reproduced previously conjectured results for both Drell-Yan production [31] and Higgs boson production in the large top mass limit [33].In particular, we have verified the observation that the LL NLP contribution at a given order in perturbation theory is generated by including a subleading term in the DGLAP kernels that accompany the leading pole in in the unsubtracted cross-section.Our reasoning provides a proof of one of the ingredients building up the resummation ansatz proposed in ref. [36], which was partly based on the idea of exponentiating NLP contributions to DGLAP splitting functions.We note again that it is natural, in this context, to exponentiate NLP contributions to the splitting functions also beyond leading order in perturbation theory: this step is strongly suggested by the arguments in ref. [117], which were, in turn, based on the idea of reciprocity between time-like and space-like splitting kernels.Ref. [36] verified that the inclusion in the Sudakov exponent of NLP terms in the NLO DGLAP kernel is responsible for the bulk of next-to-leading logarithms at NLP in the Drell-Yan and DIS cross sections.On the other hand, it is clear that, beyond leading NLP logarithms, hard collinear effects and phase space corrections become relevant, and a full resummation can only be achieved by including in the initial factorisation the contributions of radiative jet functions, as done for example in refs.[52,53].For the Drell-Yan cross section, we have also compared our results with a recent analysis based on Soft-Collinear Effective Theory techniques [76] (see also ref. [75]), finding complete agreement.
There are many directions for further work.First, of course, is the extension of the present results to subleading logarithmic accuracy at NLP.This will require a proper treatment of non-factorising phase-space effects for real emission contributions, and a thorough study of the radiative jet functions introduced in [27,[52][53][54][55][56].The latter have yet to be fully classified in QCD, while considerable progress was recently achieved in SCET [121].We note that the quark radiative jet function needed for quark annihilation processes into electroweak final states is currently known to one-loop order [52,53], which will consti-JHEP11(2019)002 tute a key ingredient to extend the present work to subleading NLP logarithms.A second direction for further studies is the inclusion of processes with final state partons at Born level: in these cases, additional threshold contributions associated with hard collinear real radiation are expected, as happens at leading power in the threshold variable.An analysis of processes of this kind was performed very recently in refs.[73,74].When more than one parton is present in the final state, further complications due to non-trivial colour flow will have to be handled, as was the case at leading power.
In order to move towards phenomenological applications of this formalism, another required step will be the inclusion of threshold contributions arising beyond leading order from different partonic channels, that are not available at Born level.For example, for the Drell-Yan process, the quark-gluon channel enters at NLO, and it generates Sudakov logarithms suppressed by an overall power of the threshold variable, because of the required radiation of a final state fermion.The inclusion of such contributions is necessary for consistent treatment of (resummed) NLP threshold effects.Another important issue that will need to be studied in detail, in order to gauge the impact of NLP resummation on phenomenology, is related to exponentiation: as pointed out in this paper, when including NLP corrections, exponentiation has to be understood in a limited sense, since NLP terms in the Sudakov exponent will generate a large set of potentially spurious contributions at NNLP and beyond upon expanding the exponential to any finite order.A precise way to limit the resummation to relevant and well-understood contributions must therefore be devised, for example by expanding the NLP part of the Sudakov exponent to fixed order, as was done in this paper.This issue is closely related to that of matching the resummation to finite order results, which is likely to be particularly relevant at NLP.Once these issues are understood, NLP resummation will provide a new versatile tool to gauge the impact of high-order corrections for a range of highly topical Standard Model and BSM processes at the Large Hadron Collider and beyond, significantly enhancing our mastery of precision high-energy phenomenology. of n Wilson line operators, as would be appropriate for contributions to the soft function involving virtual radiation.We take a number of semi-infinite Wilson lines emanating from a common vertex, and write

JHEP11(2019)002
where In path-integral language, this matrix element may be written as Note that the sum in the Wilson line term in eq.(A.4) is over both the replica numbers and the external lines, since all replicated gauge fields may interact with any given Wilson line.Furthermore, the fact that the action for the replicated theory is just the sum of the actions of individual replicas follows from the fact that replicas are non-interacting.Carrying out

JHEP11(2019)002
Considering now the application of these results to eq. (3.30), we note that the partonic factor for the (q q)-channel of the resummed Drell-Yan cross-section at LL accuracy, at O(α m s ), and in Mellin space, takes the form 2α s C F π

C Two gluon emission from the generalised Wilson line
In section 3.1, we defined a next-to-soft function in terms of generalised Wilson lines, which have been introduced and discussed extensively in refs.[28,29].These operators generate effective Feynman rules for the emission of (next-to-)soft gluons from a given hard particle, and the one-gluon emission terms required for describing radiation at O(α s ) are shown in eq.(3.2).However, as refs.[28,29] make clear, the required Feynman rules also involve effective vertices describing the emission of two gluons from the same point.These are neglected in the analysis of this paper, for reasons discussed in section 3.1.It is therefore appropriate to check explicitly in a simple example that such vertices cannot contribute to leading-logarithmic NLP terms at higher orders in perturbation theory.
Ignoring coupling and colour factors, the form of the two-gluon emission vertex from a hard scalar particle, in momentum space, is given by [28]  where p is the hard momentum of the emitting particle, and (k, l) are the soft momenta of the emitted gluons.The latter may also be sums of individual gluon momenta, which will not affect the following.Throughout the paper, we have considered processes with two incoming massless hard partons carrying four-momenta p 1 and p 2 .Without loss of generality, let us consider the two-gluon emission vertex as occuring on leg p 1 .Then, as we have argued in section 3.1, leading logarithmic effects can only come from radiation that is maximally (next-to) soft, as well as collinear.This in turn means that either k or l must be proportional to p µ 1 or p µ 2 .From eq. (C.1), it is straightforward to show that R(p; k, l) vanishes if k ∝ p or if l ∝ p.Thus, for a non-zero contribution, both k and l must be proportional to p 2 , yielding

JHEP11(2019)002
In a squared matrix element summed over final state gluon polarisations, the Lorentz indices µ and ν must ultimately be contracted with one of the external momenta p 1 or p 2 , or with a further soft momentum.However, the combination in the square brackets in eq.(C.2) acts as a projection tensor, that removes the component of any four-momentum that is collinear with p 1 or p 2 .We have already seen that leading log behaviour can only arise from soft gluon emissions that are maximally (next-to) soft and collinear.We thus find that the two-gluon emission vertex is irrelevant at LL accuracy.
To be more precise, the above discussion relates only to emissions from a scalar particle.In the case of non-zero spin, an extra contribution to the two-gluon emission vertex appears, that involves the spin generator of the emitting particle.There is however an independent line of argument that allows us to discard this contribution, and that indeed could be applied to the first term of eq.(C.1): in position space, a four-point vertex for double gluon emission necessarily involves both gluons being emitted from the same point on the emitting Wilson line, and thus involves one less propagator than contributions involving two separate gluon emissions.As a result, such contributions will not contribute a leading logarithm which, as discussed in section 3.1, requires a maximal number of integrations over normal variables.We therefore conclude, also in the case of spinning hard particles, that the two-gluon next-to-soft emission vertex can be neglected at LL accuracy.
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.

Figure 3 .
Figure 3. (a) Non-radiative amplitude with two incoming particles and a hard interaction H; (b) external emission contribution; (c) internal emission contribution.

Figure 4 .
Figure 4. Labelling of spin indices for the squared amplitude, where H is the hard function, and F i a generalised Wilson line.

Figure 6 .
Figure 6.Production of N colour singlet particles with (a) a q q initial state; (b) a gg initial state.

Figure 7 .
Figure 7. Higgs boson production via gluon-gluon fusion, where • denotes the effective coupling resulting from the integration of the top quark loop.

Figure 8 .
Figure 8.(a) Example diagram generated by the path integral in eq.(A.3), with straight semiinfinite Wilson lines for illustration, containing three subdiagrams.(b) Example diagram in the replicated theory, with different colours denoting different replicas.

Φ
i e iS(Aµ, ψ,ψ)= DA µ exp n i=1 ie dx µ i A µ (x i ) + iS A µ , ψ, ψ .(A.3)where S A µ , ψ, ψ is the QED action.Carrying out the path integral generates Feynman diagrams in which multiple Wilson lines are connected by subdiagrams consisting of photons and fermion loops, as shown for example in figure8(a).Now let us generate N independent copies or replicas of the gauge and fermion fields, labelled by {A (j) µ } and {ψ (j) }, such that particle species with different replica number j never interact.The soft function in such a theory, involving the same n Wilson lines (which are not replicated) is given by S n,R = DA(1)  µ . . .DA(N )