Universality of next-to-leading power threshold effects for colourless final states in hadronic collisions

TL;DR

This paper demonstrates a universal form for next-to-leading power threshold effects in the production of colourless particles in hadronic collisions, applicable at NLO and beyond, with new analytic results for Higgs pair production.

Contribution

It introduces a universal expression for NLP threshold effects in colourless final states, applicable to both tree-level and loop-induced processes, and provides new analytic results for Higgs pair production.

Findings

Universal NLP threshold effect formula at NLO

Re-derivation of Higgs and gauge boson production results

New analytic results for Higgs pair production with top-mass dependence

Abstract

We consider the production of an arbitrary number of colour-singlet particles near partonic threshold, and show that next-to-leading order cross sections for this class of processes have a simple universal form at next-to-leading power (NLP) in the energy of the emitted gluon radiation. Our analysis relies on a recently derived factorisation formula for NLP threshold effects at amplitude level, and therefore applies both if the leading-order process is tree-level and if it is loop-induced. It holds for differential distributions as well. The results can furthermore be seen as applications of recently derived next-to-soft theorems for gauge theory amplitudes. We use our universal expression to re-derive known results for the production of up to three Higgs bosons at NLO in the large top mass limit, and for the hadro-production of a pair of electroweak gauge bosons. Finally, we present new analytic results for Higgs boson pair production at NLO and NLP, with exact top-mass dependence.

Equations

$$\frac{d\sigma }{d\xi }=K_{\mathrm{ew}}{\left(4\pi {\alpha }_s\right)}^{n_0}\sum _{n=0}^{\infty }{\left(\frac{{\alpha }_s}{\pi }\right)}^n\sum _{m=0}^{2n-1}\left[c_{nm}^{(-1)}{\left(\frac{{\log }^m\xi }{\xi }\right)}_{+}+c_n^{(\delta )}\delta (\xi )+c_{nm}^{(0)}{\log }^m\xi +\dots \right],$$

$$\mathcal{A}(\{p_i\})=\tilde{H}\left(\{p_i\},n_1,n_2\right)\tilde{\mathcal{S}}\left(p_1,p_2\right)\frac{{\prod }_{k=1}^{2}J(p_k,n_k)}{{\prod }_{k=1}^2\tilde{\mathcal{J}}(p_k,n_k)}.$$

$${\mathcal{A}}_{\mu }^a\left(\{p_i\},k\right)$$

$$G_l^{\mu \nu }={\eta }^{\mu \nu }-\frac{(2p_l-k{)}^{\nu }}{2p_l\cdot k-k^2}{k}^{\mu }.$$

$$\tilde{\mathcal{S}}\left(p_1,p_2\right)=J(p_i,n_i)=\tilde{\mathcal{J}}\left(p_i,n_i\right)=1,n_i^2=0,$$

$${\mathcal{A}}_{\mu }^a\left(\{p_i\},k\right)$$

$${\cal A}^{(1)}_{\mu,a}\left(\{p_{i}\},k\right)\,=\,\sum_{l=1}^{2}\Bigg{[}g_{s}\,{\bf T}_{l,\,a}\,G^{\nu}_{l,\mu}\,\frac{\partial}{\partial p_{l}^{\nu}}+J^{(1)}_{\mu,\,a}\left(p_{l},n_{l},k\right)\Bigg{]}{\cal A}^{(0)}\left(\{p_{i}\}\right)\,,$$

$$J_{\mu }^a(p,n,k)=g_s{\mathbf{T}}^a\left[\frac{(2p-k{)}_{\mu }}{2p\cdot k}+\frac{\mathrm{i}{k}^{\beta }}{p\cdot k}{S}_{\beta \mu }\right],S_{\beta \mu }=\frac{\mathrm{i}}{4}\left[{\gamma }_{\beta },{\gamma }_{\mu }\right],$$

$$J_{\mu ,\rho \sigma }^a(p,n,k)=g_s{\mathbf{T}}^a\left[\frac{(2p-k{)}_{\mu }}{2p\cdot k}{\eta }_{\rho \sigma }-\frac{\mathrm{i}{k}^{\beta }}{p\cdot k}{M}_{\beta \mu ,\rho \sigma }\right],$$

$$M_{\beta\mu,\,\rho\sigma}\,=\,{\rm i}\Big{(}\eta_{\beta\rho}\eta_{\mu\sigma}-\eta_{\beta\sigma}\eta_{\mu\rho}\Big{)}\,.$$

$${\left[{\mathbf{T}}_a\right]}_{bc}=\mathrm{i}{f}_{abc}.$$

$$G_{l,\mu }^{\nu }\frac{\partial }{\partial p_l^{\nu }}=\frac{k^{\nu }}{p_l\cdot k}\left[p_{l,\nu }\frac{\partial }{\partial p_l^{\mu }}-p_{l,\mu }\frac{\partial }{\partial p_l^{\nu }}\right]=-\frac{\mathrm{i}{k}^{\nu }{\mathrm{L}}_{\nu \mu }^{(l)}}{p_i\cdot k},$$

$${\mathcal{A}}_{\mu ,a}^{(1)}\left(\{p_i\},k\right)$$

$$\sum _{i=1}^2{p}_i^{\mu }=\sum _{i=3}^{N+2}{p}_i^{\mu }\equiv P^{\mu },$$

$$z=\frac{P^2}{s},$$

$${\mathcal{A}}_{\mathrm{NLP}}^{\sigma ,\mu \nu }={\mathcal{A}}_{scal.}^{\sigma ,\mu \nu }+{\mathcal{A}}_{\mathrm{spin}}^{\sigma ,\mu \nu }+{\mathcal{A}}_{orb.}^{\sigma ,\mu \nu },$$

$${\mathcal{A}}_{scal.}^{\sigma ,\mu \nu }=\mathrm{i}{g}_s{f}_{abc}\left[\frac{{\left(2p_1-k\right)}^{\sigma }}{2p_1\cdot k}-\frac{{\left(2p_2-k\right)}^{\sigma }}{2p_2\cdot k}\right]{\mathcal{A}}^{\mu \nu }$$

$${\mathcal{A}}_{\mathrm{spin}}^{\sigma ,\mu \nu }=\mathrm{i}{g}_s{f}_{abc}\left[\frac{{\mathcal{A}}^{\alpha \nu }}{p_1\cdot k}\left(k^{\mu }{\delta }_{\alpha }^{\sigma }-k_{\alpha }{\eta }^{\mu \sigma }\right)-\frac{{\mathcal{A}}^{\mu \alpha }}{p_2\cdot k}\left(k^{\nu }{\delta }_{\sigma }^{\alpha }-k_{\alpha }{\eta }^{\nu \sigma }\right)\right].$$

$${\mathcal{A}}_{orb.}^{\sigma ,\mu \nu }=\mathrm{i}{g}_s{f}_{abc}\left[G_1^{\alpha \sigma }\frac{\partial {\mathcal{A}}^{\mu \nu }}{\partial p_1^{\alpha }}-G_2^{\alpha \sigma }\frac{\partial {\mathcal{A}}^{\mu \nu }}{\partial p_2^{\alpha }}\right].$$

$${\left|{\mathcal{A}}_{\mathrm{NLP}}\right|}^2$$

$${\mathcal{P}}_{\alpha \beta }(p,l)\equiv \sum _{\lambda }{ϵ}_{\alpha }^{(\lambda )}(p){ϵ}_{\beta }^{(\lambda )\ast }(p)=-{\eta }_{\alpha \beta }+\frac{p_{\alpha }{l}_{\beta }+p_{\beta }{l}_{\alpha }}{p\cdot l},$$

$${\left|{\mathcal{A}}_{\mathrm{NLP}}\right|}^2=\sum _{\mathrm{colours}}\left\{{\left|{\mathcal{A}}_{scal.}^{\sigma ,\mu \nu }\right|}^2+2\mathrm{Re}\left[{\left({\mathcal{A}}_{\mathrm{spin}}^{\sigma ,\mu \nu }+{\mathcal{A}}_{orb.}^{\sigma ,\mu \nu }\right)}^{\ast }{\mathcal{A}}_{scal.\sigma ,\mu \nu }\right]\right\},$$

$$\sum _{\mathrm{colours}}{\left|{\mathcal{A}}_{scal.}^{\sigma ,\mu \nu }\right|}^2=2g_s^2{N}_c\left(N_c^2-1\right)\frac{p_1\cdot p_2}{p_1\cdot k{p}_2\cdot k}|{\mathcal{A}}_{\mu \nu }{|}^2,$$

$$\sum _{\mathrm{colours}}2\mathrm{Re}\left[{\mathcal{A}}_{\mathrm{spin}}^{\sigma ,\mu \nu }{\mathcal{A}}_{scal.\sigma ,\mu \nu }^{\ast }\right]$$

$$\sum _{\mathrm{colours}}2\mathrm{Re}\left[{\mathcal{A}}_{orb.}^{\sigma ,\mu \nu }{\mathcal{A}}_{scal.\sigma ,\mu \nu }\right]$$

$$\delta p_1^{\alpha }=-\frac{1}{2}\left(\frac{p_2\cdot k}{p_1\cdot p_2}{p}_1^{\alpha }-\frac{p_1\cdot k}{p_1\cdot p_2}{p}_2^{\alpha }+k^{\alpha }\right),\delta p_2^{\alpha }=-\frac{1}{2}\left(\frac{p_1\cdot k}{p_1\cdot p_2}{p}_2^{\alpha }-\frac{p_2\cdot k}{p_1\cdot p_2}{p}_1^{\alpha }+k^{\alpha }\right).$$

$${\left|{\mathcal{A}}_{\mathrm{NLP}}\right|}^2=\frac{2g_s^2{N}_c\left(N_c^2-1\right)p_1\cdot p_2}{p_1\cdot k{p}_2\cdot k}{\left|{\mathcal{A}}_{\mu \nu }\left(p_1+\delta p_1,p_2+\delta p_2\right)\right|}^2.$$

$$s\to {\left(p_1+p_2+\delta p_1+\delta p_2\right)}^2=s+2\left(\delta p_1+\delta p_2\right)\cdot (p_1+p_2).$$

$$p_1+p_2=k+P,$$

$$s\to zs.$$

Full text

NIKHEF/2017-25 Edinburgh 2017/10 QMUL-PH-17-07 ARC-17-03

Universality of next-to-leading power threshold effects for colourless final states in hadronic collisions

V. Del Duca1,2, E. Laenen3,4,5, L. Magnea6, L. Vernazza7, C. D. White8

1*ETH Zurich, Institut fur theoretische Physik, Wolfgang-Paulistr. 27, 8093, Zurich, Switzerland

2*INFN Laboratori Nazionali di Frascati, 00044 Frascati (Roma), Italy

3Nikhef, Science Park 105, NL–1098 XG Amsterdam, The Netherlands

4ITFA, University of Amsterdam, Science Park 904, Amsterdam, The Netherlands

5ITF, Utrecht University, Leuvenlaan 4, Utrecht, The Netherlands

6*Dipartimento di Fisica and Arnold-Regge Center, Università di Torino,

and INFN, Sezione di Torino, Via P. Giuria 1, I-10125 Torino, Italy*

7*Higgs Centre for Theoretical Physics, School of Physics and Astronomy, The University of Edinburgh, Edinburgh EH9 3JZ, Scotland, UK

8Centre for Research in String Theory, School of Physics and Astronomy, Queen Mary University of London, 327 Mile End Road, London E1 4NS, UK

We consider the production of an arbitrary number of colour-singlet particles near partonic threshold, and show that next-to-leading order cross sections for this class of processes have a simple universal form at next-to-leading power (NLP) in the energy of the emitted gluon radiation. Our analysis relies on a recently derived factorisation formula for NLP threshold effects at amplitude level, and therefore applies both if the leading-order process is tree-level and if it is loop-induced. It holds for differential distributions as well. The results can furthermore be seen as applications of recently derived next-to-soft theorems for gauge theory amplitudes. We use our universal expression to re-derive known results for the production of up to three Higgs bosons at NLO in the large top mass limit, and for the hadro-production of a pair of electroweak gauge bosons. Finally, we present new analytic results for Higgs boson pair production at NLO and NLP, with exact top-mass dependence.

1 Introduction

The remarkable experimental precision and high statistics offered by current and future colliders necessitates the calculation of high-order effects in QCD perturbation theory for processes of increasing complexity. The calculation of next-to-leading order (NLO) QCD corrections has reached the stage where automated tools are being used and multi-jet production rates are evaluated (see, for example [1]). Attention then focuses on QCD corrections at NNLO or above, and on the inclusion of higher-order electroweak effects. In this context, processes which are induced by loop effects present special difficulties: typically, the leading order contribution involves a loop containing heavy particles, such as top quarks or electroweak vector bosons, and is often computed in the context of an effective field theory. NLO QCD corrections with exact dependence on the heavy particle masses then involve intricate two-loop calculations: for multi-particle final states, these corrections are not known, and even for two-particle final states they are typically only known approximately, or in some cases numerically. Clearly, in all these cases it is desirable, where possible, to improve upon existing results by providing analytic information. Even partial information can be useful, since it can be used, for example, to speed up numerical codes, and also to provide additional consistency checks.

In this paper, we consider the production of a generic colour-singlet final state in hadronic collisions, and we study the effects of additional gluon radiation near partonic threshold, where emitted gluons have a low energy or transverse momentum with respect to the incoming partons. In such cases, one commonly defines a dimensionless threshold variable $\xi$, vanishing at the threshold, and it is well known that the corresponding differential cross section has the generic form

[TABLE]

where the overall factor is associated with the leading-order cross section: for loop-induced processes, this may be proportional to a power of the strong coupling, as indicated, while $K_{\rm ew}$ contains electroweak couplings. The first two sets of terms on the right-hand side of Eq. (1.1) originate from soft and collinear radiation (real or virtual), and correspond to the leading power in the threshold variable, and to corrections localised at the threshold, respectively. These contributions are known to have a universal (process-independent) form, that permits their resummation to all orders in perturbation theory. This resummation is well understood and widely applied, and can be performed within a variety of approaches (see for example [2, 3, 4, 5, 6, 7, 8, 9, 10, 11]). Even without a full resummation, a fixed order evaluation of these terms can be useful in estimating higher-order corrections to the cross-section, when these are not known (see for example [12] for a review).

The third set of terms on the right-hand side of Eq. (1.1) defines next-to-leading power (NLP) contributions in the threshold variable, roughly corresponding to gluon radiation that can be next-to-soft or collinear. Although power-suppressed, these terms are still singular as $\xi\rightarrow 0$, and can be numerically significant: for example, they contribute significantly to the theoretical uncertainty for Higgs boson production in the gluon fusion channel [13], and their numerical impact has been explicitly confirmed by the recent calculation of this process at N3LO [14, 15, 16, 17, 18]. A full, generally applicable resummation prescription for NLP contributions is not presently known, even in the relatively simple case of parton annihilation into electroweak final states: the problem has been intensively studied in recent years, and partial progress has been made using a variety of methods [19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39], buiding upon the earlier work of [40, 41, 42]. Even without a full resummation, however, knowledge of NLP contributions at fixed order can provide useful analytic information where this is missing, as well as furnishing improved approximations for unknown higher-order cross sections. This is especially true if one can derive universal properties of NLP contributions, applicable at a given order for a broad class of processes.

In this paper, we examine the universal properties of NLP radiation in the hadro-production of an arbitrary number of colour-singlet particles at NLO. We will prove that the NLO cross-section for this class of processes can be written in terms of the leading-order cross-section with shifted kinematics, convolved with a simple universal $K$-factor. More precisely, we find that the $K$-factor does not depend on the spin of the emitting parton, and depends on the color representation only through a trivial replacement of colour factors. Our starting point will be a factorisation formula for NLP effects at amplitude level recently derived in Refs. [23, 24], which expresses the effect of adding an additional gluon to an arbitrary hard process with an electroweak final state in terms of universal functions. Our results provide a useful testing ground for this formula, and an examination of its simpler consequences, aside from its role as a basis for resummation. Furthermore, our results have more practical consequences, providing an analytic approximation to the NLO cross section for a number of interesting loop-induced processes for which only limited information is available. Interestingly, the universality of the result extends to differential distributions, provided the shift in LO kinematics is properly understood. From a theoretical point of view, the universality and simplicity of our results at NLO can be seen as a consequence of recently derived next-to-soft theorems [43, 44, 45] for radiative tree-level gauge theory amplitudes: we note however that our factorisation formula is an all-order result, and therefore will yield more general results, once the appropriate ingredients are computed at the relevant orders.

From the point of view of phenomenology, the most interesting applications of our results will concern the production of Higgs bosons in the gluon fusion channel, possibly in association with electroweak gauge bosons. We will however show explicitly that the formalism can be used also with (anti)quark initial states, and compare our results with existing calculations. In the gluon fusion channel, we will begin by showing how known properties of the single Higgs boson cross section emerge as a special case of our result. We next move on to multiple Higgs boson production, which has been the focus of much recent research. Beyond the leading-order result [46], we note that analytic expressions for the cross section are known only in the large top mass limit for Higgs pair production [47, 48] and for triple Higgs boson production [49]. In the case of Higgs pair production, leading order results with full top mass dependence were obtained in refs. [50, 51], and numerical results have recently been presented at NLO [52] (see also [53]); leading power threshold corrections have been considered in Ref. [54], and corrections to the large top mass approximation in Ref. [55]. In the case of triple Higgs boson production, numerical results with full top mass dependence at leading order and for real radiation at NLO were obtained in Ref. [56] and in Ref. [57], respectively. Associated production of electroweak bosons and Higgs bosons in the gluon fusion channel, discussed in Ref. [58], also falls within the scope of our method, although we will not discuss it in detail here.

In this paper, we go beyond previous analytic results for Higgs boson pair production, by providing NLO corrections, up to NLP in the threshold variable, with full top mass dependence. As a further illustration and check, we demonstrate consistency with known results for triple Higgs boson production in the large top mass limit [49], diphoton production [59], and the production of $W^{+}W^{-}$ pairs [60]. These results serve as an illustration of the method: we postpone a detailed phenomenological analysis, as well as applications to triple Higgs production with full top mass dependencce, and to associated production of Higgs bosons with Z bosons, to future work. We note in passing that the present results also provide a strong consistency check on all future NLO analytic computations of loop-induced processes of colour-singlet particles, which of course must agree with the simple factorised expressions we derive at NLP.

The structure of our paper is as follows. In Section 2, we briefly review the next-to-soft factorisation formula of Refs. [23, 24], before describing how it can be extended for use in gluon-induced processes. In Section 3, we derive an explicit expression for the NLO cross-section of a colour-singlet final state in gluon fusion, valid up to NLP level. A similar result for the quark channel is derived in Section 4. In Section 5, we show how known results in single Higgs boson production are reproduced, before examining multiple Higgs boson production in Section 6. Vector boson pair production is considered in Section 7. Finally, we discuss our results and future prospects in Section 8.

2 NLP amplitude factorisation

In this section, we briefly review the results of Refs. [23, 24], which derive a factorisation formula for QCD radiation up to NLP level, and we provide a generalisation of these results to the case of external incoming gluons. Consider an amplitude with two incoming partons of momenta $p_{1}$ and $p_{2}$, and any number $N$ of final state colour singlet particles, with momenta $p_{i}$, $3\leq i\leq N+2$. As described in detail in Ref. [24], any such amplitude may be written in the factorised form

[TABLE]

Here $\tilde{\cal S}(p_{1},p_{2})$ is a next-to-soft function, dressing the hard interaction with virtual exchanges of (next-to-)soft gluons, which couple to the external partons through the next-to-eikonal Feynman rules described in Refs. [20, 21]. Associated with each hard parton is a jet function $J(p_{i},n_{i})$, which collects collinear singularities, and which depends on an auxiliary vector $n_{i}$. The next-to-eikonal jet $\,\tilde{\cal J}(p_{i},n_{i})$ corrects for the double counting of contributions from gluons which are both (next-to-)soft and collinear, and finally the hard function $\tilde{H}(\{p_{i}\},n_{1},n_{2})$ is defined by matching to the amplitude on the left-hand side of Eq. (2.1), so that all dependence on the auxiliary vectors $\{n_{i}\}$ cancels out. If one ignores the presence of next-to-eikonal Feynman rules, Eq. (2.1) reduces to the well-known soft-collinear factorisation formula, describing the dressing of a given hard interaction process with leading-power soft and collinear radiation (see, for example, [61]). The form of Eq. (2.1), however, is a crucial intermediate step in considering the emission of an additional gluon, of momentum $k^{\mu}$ and colour $a$. Up to NLP in this momentum, the resulting amplitude is given by [24]

[TABLE]

where $g_{s}$ is the QCD coupling111We absorb a factor $\mu^{\epsilon}$, where $\mu$ is the dimensional regularisation scale, into the coupling for simplicity., ${\bf T}_{i}^{\,a}$ a colour generator on line $i$ with adjoint index $a$, and we have introduced the tensor [62]

[TABLE]

In addition to the functions already appearing in Eq. (2.1), Eq. (2) contains two more universal functions. First, the radiative next-to-soft function $\tilde{\cal S}_{\mu}^{\,a}$ is a matrix element of next-to-eikonal Wilson lines directed along the directions of the incoming partons, like the virtual next-to-soft function $\tilde{\cal S}$, but with a single gluon present in the final state. Furthermore, Eq. (2) includes a radiative jet function $J_{\mu}^{\,a}$ collecting all contributions associated with the emission of a gluon from the $i^{\rm th}$ parton, and enhanced by virtual collinear poles. This function was first introduced in the context of abelian gauge theory in Ref. [42], and its definition was recently generalised to non-abelian theories in Ref. [24]. The radiative functions can be defined in terms of operator matrix elements, but for our present NLO analysis, where radiative functions enter only at tree level, a diagrammatic definition is sufficient.

The final term on the right-hand side of Eq. (2) is a subtraction term that removes any double counting of contributions occuring in both the radiative next-to-soft emission function, and in the radiative jet emission functions: it can be obtained simply by taking the next-to-soft limit of the radiative jet function. As was done in Ref. [24], Eq. (2) can be considerably simplified by using renormalisation group arguments and computing the right-hand side in the bare theory and with light-like reference vectors $n_{i}^{2}=0$ for the jets. With these choices, one can use the bare quantities

[TABLE]

and the amplitude can be written as

[TABLE]

If we now focus on the NLO contributions to the cross section, there is a further significant simplification: indeed, the leading-order term in the next-to-soft emission function, $\tilde{\cal S}_{\mu,\,a}^{(1)}$ consists of single gluon emissions from the hard incoming partons, and these contributions are completely cancelled [24] by the leading-order subtraction term ${\cal A}^{(1),\tilde{\cal J}_{l}}_{\mu,\,a}$, leaving

[TABLE]

which expresses the complete one-gluon radiative amplitude at NLO and NLP in terms of the Born amplitude. Note that in Eq. (2.6) the non-radiative amplitude and jet emission functions are understood as carrying implicit spin indices, depending on the identity of the particle species in each jet.

The quark radiative jet function at leading order is simply given by the emission of a single gluon from the incoming (anti)quark [42, 24], as shown in Fig. 1(a). Evaluating the diagram gives222Note that, by definition, the radiative jet function does not include the spinor wave function for the external quark line.

[TABLE]

where we have decomposed the result into spin-dependent and spin-independent parts, introducing the generator $S_{\beta\mu}$ of Lorentz transformations on spinors. Note that at leading order the quark jet function is independent of the auxiliary vector $n^{\mu}$, consistently with Eq. (2.6), which represents a physical amplitude and cannot depend on $n$. For the gluon radiative jet function, at leading order, we can simply use a diagrammatic definition, analogous to the radiative quark jet, and shown in Fig. 1(b).

Restoring explicit spin indices for the external gluon, we can write the result of this diagram as

[TABLE]

where we have introduced the generator of Lorentz transformation acting on vector fields,

[TABLE]

Once again, we have decomposed the kinematic part into its spin-dependent and spin-independent parts (see for example [63]). The colour operator for the gluon case can be explicitly interpreted as

[TABLE]

Turning now to the derivative contribution to Eq. (2.6), one may note that the action of the projector $G^{\mu\nu}_{l}$ defined in Eq. (2.3), up to NLP order, can be recast in terms of the orbital angular momentum of parton $l$. Indeed, to this accuracy

[TABLE]

where ${\rm L}^{(l)}_{\nu\mu}$ is the orbital angular momentum operator associated with the $l^{\rm th}$ parton. Using Eqs. (2.8, 2.9), we can now rewrite Eq. (2.6) in a unified notation for quarks and gluons, as

[TABLE]

where in the first line ${\rm\Sigma}^{(l)}_{\nu\mu}$ is the spin angular momentum operator for parton $l$, in the relevant representation of the Lorentz group, acting as $-S^{(l)}_{\nu\mu}$ for spin one half, and as $M^{(l)}_{\nu\mu}$ for spin one, while ${\rm J}^{(l)}_{\nu\mu}$ is the total angular momentum operator. Furthermore, in the second line, we have omitted the term proportional to $k_{\mu}$, which gives a vanishing contribution when contracted with a physical polarisation vector for the emitted gluon.

Equation (2.12) is recognisable as the recently derived next-to-soft theorem [43], which mirrors a similar result derived in gravity [44, 45]. As noted, this formula encompasses both the quark and gluon cases, provided the spin operator is interpreted appropriately, validating our diagrammatic definition for the leading order gluon radiative jet function. For the NLO analysis performed in this paper, we could in fact have simply adopted Eq. (2.12) as the starting point for our following analysis; note, however, that Eq. (2) and Eq. (2.5) are much more general results, applicable in principle to any order in perturbation theory.

3 Colour-singlet particle production in the gluon channel

In this section, we apply the result of Eq. (2.12) to obtain a general expression for the NLO cross-section for the production of $N$ colour-singlet particles near threshold. We begin by considering the gluon-induced process shown in Fig. 2, while we will turn to the quark-induced process in Section 4. At Born level, the momenta introduced in Fig. 2 satisfy the leading-order momentum conservation condition

[TABLE]

with the Born-level centre-of-mass energy squared given by $s\,=\,P^{2}$. Beyond Born level, we may define the dimensionless variable

[TABLE]

which represents the fraction of the partonic centre-of-mass energy carried into the final state by all colour singlet particles. At leading order obviously $z=1$; beyond leading order, additional real radiation may be emitted, in which case $0\leq z\leq 1$, and $\xi\equiv 1-z$ is a dimensionless threshold variable of the kind introduced in Eq. (1.1). In particular, at NLO only a single gluon can be emitted, and all contributions up to NLP in the emitted momentum $k$ are captured by Eq. (2.12). We can then use this to obtain a cross-section formula that is correct up to the first sub-leading order in $\xi$. To this end, it is useful to write the complete radiative amplitude (before contraction with external gluon polarisation vectors) as

[TABLE]

where $\sigma$ is the Lorentz index of the emitted gluon, while $\mu$ and $\nu$ are the Lorentz indices of the incoming gluons, and, for simplicity, we have suppressed momentum dependence, colour indices and the superscript denoting the perturbative order. The three terms on the right-hand side correspond to the scalar, spin-dependent and orbital angular momentum terms in Eq. (2.12) respectively. The colour indices of the incoming gluons are displayed in Fig. 2, and we note that by colour conservation the leading order amplitude must be proportional to $\delta^{bc}$. Following Eqs. (2.10, 2.12), one may write the scalar-like contribution to the amplitude as

[TABLE]

where, as above, we omitted the superscript denoting the perturbative order for the Born amplitude ${\cal A}^{\,\mu\nu}$. Using Eqs. (2.9, 2.12), the spin-dependent contribution to the amplitude is given by

[TABLE]

Finally, the orbital angular momentum contribution is

[TABLE]

After including polarisation vectors for the two incoming gluons, the squared matrix element, accurate to NLP level and summed over polarisations and colours, is

[TABLE]

where we defined the polarisation sum

[TABLE]

with $l$ is an arbitrary light-like reference vector used to define physical polarisation states, whose dependence must cancel in the final result. Alternatively, one could sum over all polarisations, using ${\cal P}_{\alpha\beta}=-\eta_{\alpha\beta}$, and correct for this by including external ghost contributions. Following this second approach, it is fairly easy to conclude that ghost contributions vanish at NLP: indeed, final state ghost emission is suppressed by a power of the energy at amplitude level, and thus contributes at NNLP at cross section level; furthermore, diagrams with a ghost-antighost pair in the initial state do not couple directly to fermions or to the Higgs boson, and are strongly suppressed. These expectations are borne out by a direct calculation, showing that all terms proportional to the vector $l^{\mu}$ in Eq. (3.7) are beyond the required accuracy. We conclude that we can perform a sum over all polarisation, so that Eq. (3.7) simplifies to

[TABLE]

where in the second term we need to keep only those terms which are leading power in the scalar part of the amplitude. It is straightforward to show that the first term on the right-hand side yields

[TABLE]

so that only the leading power term survives. For the gluon-initiated process we are considering in this section, the spin term in Eq. (3.9) can be shown to vanish upon summing over all polarisations. For example, the spin contribution from the first leg is

[TABLE]

The prefactor in the second line is antisymmetric under the interchange of $\alpha$ and $\mu$, and thus vanishes when contracted with the squared Born amplitude, which is symmetric; the same argument applies to the second incoming gluon. Note that the argument applies also when the Born amplitude is loop induced, and thus may acquire an imaginary part (as is the case here). The orbital angular momentum contributions give

[TABLE]

where we defined

[TABLE]

Note that these shifts are proportional to the soft momentum $k$ and transverse to their respective momenta, $p_{i}\cdot\delta p_{i}=0$. This second property follows from the fact that the $i^{\rm th}$ momentum shift is derived from the orbital angular momentum operator of Eq. (2.11), which generates an infinitesimal Lorentz transformation transverse to the momentum $p_{i}$. Combining Eq. (3) with Eq. (3.10), we can write

[TABLE]

Eq. (3.14) is a focal point of this paper: it shows that all NLP contributions to the NLO squared matrix element for the production of an arbitrary colour-singlet final state can be absorbed into a shift in the kinematics of the Born contribution. Corrections to this shifting procedure involve terms at least quadratic in $\delta p_{i}$, and thus beyond the NLP approximation. In Section 4, we will show that the same property is shared by quark-initiated processes. Note that Eq. (3.14) is fully differential in final state momenta, and can be applied to generate distributions valid at NLO and to NLP accuracy. On the other hand, using simple properties of phase space, one can also derive a similarly simple expression for the inclusive cross section. In order to do so, note that the effect of the required momentum shifts on the partonic centre-of-mass energy is given by a simple rescaling. Indeed,

[TABLE]

Substituting the definitions of Eq. (3.13) in Eq. (3.15), and using Eq. (3.2), together with the NLO momentum conservation condition

[TABLE]

it is easy to show that Eq. (3.15) can be written simply as

[TABLE]

To construct the partonic cross-section, we must now introduce the appropriate factors to average over initial state colours and spins, integrate over the ($N$ + 1)-body final state phase space, and include the flux factor. We find

[TABLE]

where

[TABLE]

denotes the $n$-body Lorentz-invariant phase space for a process with total final state momentum $Q=\sum_{i}q_{i}$, and q_{i}^{\mu}=(E_{i},\mbox{\boldmathq}_{i}) in a suitable frame. For the phase space, we may use the well-known result

[TABLE]

factorising the phase space of the $N$ colour-singlet particles from a two-body phase space involving the total momentum of the colourless system, and the additional gluon momentum $k$. The latter can be written more explicitly by parametrising

[TABLE]

Introducing the variable

[TABLE]

one then finds (see for example Ref. [21] for a recent derivation)

[TABLE]

Using Eq. (3.23) together with Eqs. (3.2, 3.14) in Eq. (3.18), one then finds

[TABLE]

where we reinstated the explicit dependence on the dimensional regularisation scale $\mu$, and we denoted by $d\Phi_{N}^{(z)}$ the phase space for $N$ (colour-singlet) particles with a partonic centre-of-mass energy shifted according to Eq. (3.17). We may easily rewrite this result in terms of the leading-order cross section with shifted kinematics, which is given by

[TABLE]

This leads us to our second central result: a simple factorised expression for the inclusive cross section, valid at NLO and NLP for the production of a generic colour-singlet system, which can be written as

[TABLE]

where the next-to-leading power $K$ factor is easily computed, with the result

[TABLE]

In the second line of Eq. (3.27) we expanded the result to NLP in $(1-z)$ and to finite order in $\epsilon$, and we introduced the $\overline{\rm MS}$ scale $\overline{\mu}^{2}=\mu^{2}\,{\rm e}^{\,\ln(4\pi)-\gamma_{E}}$ and the plus distributions

[TABLE]

Eqs. (3.26, 3.27) show explicitly that the NLO K-factor for the production of $N$ colour-singlet particles in the gluon channel is simple and universal, up to next-to-leading power in the threshold variable. This is a powerful constraint, and we will discuss some specific examples in the following sections. First, however, we consider an analogous formula in the quark channel.

4 Colour-singlet particle production in the quark channel

In the previous section, we have derived an explicit universal K-factor for multiple colour-singlet particle production in the gluon-gluon channel. In this section, we consider the cross section for quark-induced production of colour-singlet particles, and show that an identical result holds, up to a trivial replacement of colour factors. The universality of the result is not obvious from the outset, and it comes about through an interesting reshuffling of the contributions of spin and angular momentum operators, as compared to the gluon-induced process. We take the leading order process shown in Fig. 3, and consider the radiation of an additional gluon from the incoming quark and antiquark lines.

One may write the LO amplitude as

[TABLE]

where $i$, $j$ are the colour indices of the incoming quark and antiquark, and the factor ${\cal A}(\{p_{i}\})$, matrix-valued in spinor space, is the quantity entering Eq. (2.6), namely the leading-order amplitude with external wave functions removed. Following the procedure adopted in the gluon case, we may decompose the NLO amplitude, before contraction with external spinors, according to

[TABLE]

where the three terms on the right-hand side denote the scalar-like, spin, and orbital angular momentum contributions, and we have suppressed spinor indices (as well as color labels) for brevity. For the scalar and orbital angular momentum contributions, which do not depend explicitly on the spin (apart from replacing the vector indices on the leading-order amplitude with spinor indices), the arguments of the previous section may be repeated, and one may write333Following the convention of Eq. (4.1), we do not include the colour generator $t^{a}_{ji}$ in the fundamental representation in the definition of the stripped amplitude ${\cal A}$.

[TABLE]

where the momentum shifts are defined in Eq. (3.13). Including the (anti)quark wave functions and performing color and spin sums, we then find

[TABLE]

where in the second line we have introduced the dimensionless vectors

[TABLE]

By comparing Eq. (4.4) with its LO counterpart,

[TABLE]

and using Eq. (3.17) we may promote the momentum shift in Eq. (4.4) to apply to the entire squared amplitude. This leads to

[TABLE]

Note the close resemblance of Eq. (4.7) and Eq. (3.14): they differ only by the color factor and a rescaling by a factor of $z$. We must still, however, add to Eq. (4.7) the interference between the spin-dependent part of the NLO amplitude, and the eikonal amplitude. In the gluon case, this turned out to vanish in Feynman gauge, upon summing over all gluon polarisations, which was allowed at NLP accuracy. For an incoming fermion, we find that the spin contribution does not vanish, and indeed it precisely compensates for the $z$ rescaling observed in Eq. (4.7), recovering the universality of the result.

The spin-dependent part of the NLO amplitude is given by the diagrams of Fig. 4, which evaluate to

[TABLE]

In the second line of Eq. (4.8), we anticommuted Dirac matrices and used the physical polarisation condition $k^{\sigma}\epsilon_{\sigma}(k)=0$ to write the result in a form which will be more convenient in what follows. Up to NLP accuracy in the squared amplitude, we only need to consider the interference of Eq. (4.8) with the (leading power) scalar part of the NLO amplitude. Furthermore, since we are considering the emission of a single gluon, we can sum over all polarisations, rather than restricting to physical polarisations only. The relevant contribution to the squared matrix element is then

[TABLE]

To simplify this further, we may expand the emitted gluon momentum in the Sudakov decomposition

[TABLE]

We then observe that, to linear order in $k^{\mu}$, the Dirac trace in Eq. (4.10) cannot depend on $k_{T}$. Indeed, one easily finds

[TABLE]

By comparing with the squared scalar part of the amplitude

[TABLE]

we see that the spin-dependent contribution to the squared amplitude can be obtained from the part which is leading power in the gluon momentum, simply through rescaling by the factor

[TABLE]

where we have used the momentum parametrisation of Eq. (3.21).

Combining Eq. (4.7) with Eq. (4.12), we see that the rescaling factors cancel at NLP in $(1-z)$. Indeed one may write

[TABLE]

Expanding now in powers of $(1-z)$, one gets to first order

[TABLE]

and one observes that the second line is effectively ${\cal O}(1-z)^{2}$. We find then

[TABLE]

which is precisely analogous to Eq. (3.14), except for the replacement of the colour factor, which here is associated with the fundamental rather than adjoint representation of the gauge group. Once again, at NLP, Eq. (4.16) can be used in a fully differential implementation for the final state kinematics.

Having obtained Eq. (4.16), one may form the cross-section by integrating with the $(N+1)$-body phase space, exactly as was done in the gluon case. One finds then

[TABLE]

with the same factor $K_{\rm NLP}\left(z,\epsilon\right)$, given in Eq. (3.27).

A first check on this result is that it reproduces the NLO K-factor for Drell-Yan production of a vector boson of invariant mass $Q^{2}$, where one has

[TABLE]

In this case, as for any $2\to 1$ process, the LO partonic cross section has support only on the partonic threshold: for Drell-Yan production,

[TABLE]

so that the LO cross section with shifted kinematics is

[TABLE]

The delta-function imposes the correct definition of the threshold variable at NLO, while the rest of the cross section is unaffected by the shift in kinematics. To compare with standard results, we must note that the $\overline{\rm MS}$ scale $\overline{\mu}^{2}$ is usually set equal to the final state invariant mass $Q^{2}$. To this end, one may write

[TABLE]

so that Eq. (3.26) becomes

[TABLE]

which precisely agrees with the well-known results quoted for example in Refs. [64, 21].

5 Single Higgs boson production via gluon fusion

Having presented our results for both quark- and gluon-induced colour-singlet particle production, we now examine a first significant application of the gluon result, Eq. (3.27): the single production of Higgs bosons in the gluon fusion channel. As is well known, this is is the principal production mode for Higgs bosons at the LHC, and higher-order QCD corrections have been studied in great detail and with great efforts in recent years. In the effective field theory with the top quark integrated out, they have been calculated up to N3LO in perturbation theory [65, 66, 67, 68, 15, 18, 17, 16, 14]. Top-mass effects are know exactly at NLO [69], and have been studied at NNLO as a power expansion in $m_{h}^{2}/m_{t}^{2}$ [70, 71]. Here we will see how the intricate top mass dependence at NLO simplifies considerably in the threshold region, including NLP corrections.

At leading order, the incoming gluons couple to the Higgs boson via a top-quark loop, as shown in Fig. 5(a).

The leading order cross-section for this process (see for example [65]) can be written as

[TABLE]

where $m_{h}$ and $v$ are the Higgs mass and vacuum expectation value respectively444In Eq. (5.1) we have omitted scale factors relating to the $d$-dimensional coupling $\alpha_{s}$, which amounts to the choice $\overline{\mu}=m_{h}$.. The form factor $F(\tau,\epsilon)$ depends on the dimensionless variable

[TABLE]

and it is given by [65]

[TABLE]

with a normalisation chosen so that $F(\tau,\epsilon)\to 1$ as $\tau\to 0$. The cross section with kinematics shifted according to Eq. (3.15) can then be written as

[TABLE]

Substituting this result into Eq. (3.26) and expanding in powers of $(1-z)$ and $\epsilon$ one finds

[TABLE]

It is easy to check that Eq. (5.5) agrees with the known analytic NLO result of Ref. [65] in the $m_{t}\to\infty$ limit. We note, however, that the result of Eq. (5.5) is much more informative: it includes the full dependence on the top quark mass up to NLP order, and can thus be applied for arbitrary $m_{t}$. This is a remarkable simplification of the intricate result of Ref. [69] for the full $m_{t}$ dependence: after shifting the kinematics of the leading order result, the resulting K-factor is entirely independent of the top quark mass, which makes the formula especially simple to apply in practical applications 555Indeed, we have checked that eq. (5.5) reproduces the K-factor reported in ref. [70], which features a double expansion in threshold parameter and top mass..

It is interesting to examine the anatomy of the result in Eq. (5.5) in slightly more detail. If one were to calculate the NLO cross section by starting manifestly in the large top mass limit (i.e. by using an effective field theory), the leading order graph would contain an effective point-like interaction coupling the two incoming gluons to a Higgs, as shown in Fig. 5(b). At NLO, one can radiate the extra gluon from either of the incoming gluons, and one must also include the additional effective coupling shown in Fig. 5(c), namely a point-like interaction between three gluons and a Higgs boson. If one resolves the top quark loop as in Fig. 5(a), this extra interaction corresponds to emissions from inside the top quark loop. In the NLP calculation, there is no need to include any additional diagrams to capture these contributions: they are generated precisely by the orbital angular momentum contributions in Eq. (3.6): therefore, as the above analysis reveals, we can choose to associate these terms with a shift in the kinematics of the leading order result, up to corrections subleading in soft momentum. Seen from the point of view of the effective field theory at large $m_{t}$, it is highly non-trivial that such a shift captures the contribution of higher-order operators in the effective Lagrangian.

6 Multiple Higgs boson production

In the previous section we have tested our main result, given by Eq. (3.26) for gluon scattering, by reproducing known results in the cross section for single Higgs boson production via gluon fusion. We now consider the case of Higgs boson pair production, a process of ongoing interest at the LHC, due to its potential role in extracting the Higgs boson self-coupling. Analytic results for this process are known up to NNLO in the large top mass limit [46, 47, 48], but only at leading order with full top mass dependence [50, 51]. Further studies have looked at systematically improving the effective field theory results by including leading-power threshold effects [54], or contributions suppressed by powers of the top mass [55]. Recently, numerical results at NLO accuracy with full top mass dependence have become available [52] (see also [53]). This, however, does not preclude the desire for analytic results, which can serve to improve the efficiency of numerical computations, whilst also providing clues regarding higher-order structures in perturbation theory. This is especially true in Higgs boson pair production, given that the large top mass limit is not a good approximation, unlike the case of single Higgs production. The leading order diagrams for Higgs pair production are shown in Fig. 6, and the leading order amplitude may be written as

[TABLE]

where (in the Standard Model)

[TABLE]

and $F_{\triangle}$, $F_{\Box}$, $G_{\Box}$ are form factors arising from the triangle and box graphs, as indicated by the subscripts. They depend on the Higgs boson and top masses, as well as the partonic centre of mass energy $s$ and the other Mandelstam invariants

[TABLE]

The basis tensors $T_{s,\,\mu\nu}$, with $s=0,2$, in Eq. (6.1) are associated with the exchange of spin 0 and spin 2 in the $s$ channel, respectively. Denoting the gluon momenta by $p_{1}$ and $p_{2}$ and the Higgs boson momenta by $p_{3}$ and $p_{4}$, their explicit forms are

[TABLE]

where

[TABLE]

With these notations, the leading-order distribution in the Mandelstam invariant $t$ can be written as [46]

[TABLE]

This expression simplifies considerably in the large top mass limit, where

[TABLE]

so that Eq. (6.6) becomes

[TABLE]

We observe that, in the large top mass limit, the leading-order cross section vanishes at threshold, as $s\to 4m_{h}^{2}$, due to the cancellation between the box and triangle contributions. This property is one of the reasons that make the large top mass limit a poor approximation in Higgs boson pair production, necessitating the calculation of higher order corrections with full $m_{t}$ dependence. It also causes problems when trying to define a K-factor as a function of the variable $z$. Ordinarily, one would divide the NLO cross section by the LO one, however this becomes ill-defined in the threshold region $z\to 1$. In Ref. [46] this problem is circumvented by dividing by the LO cross section with kinematics shifted according to Eq. (3.17). The resulting K-factor thus matches precisely the quantity defined in Eq. (3.26).

With this convention, the NLO cross section for Higgs boson pair production, up to NLP accuracy, can be written as

[TABLE]

where we have extracted an explicit factor of $z$ on the left-hand side, to match the conventions adopted in Ref. [46]. In the large top mass limit, Eq. (6.9) reproduces the results of Ref. [46]. As in the case of single Higgs production, however, the result is much more powerful, in that it applies to the full top mass dependence. Eq. (6.9) thus provides an explicit analytic form of the cross section, at NLO in perturbation theory, and up to NLP in the threshold variable $(1-z)$. This extends the results of Ref. [54], which considered supplementing the fixed-order cross section with threshold effects at leading power only. Ref. [55], on the other hand, studied the systematic improvement of the large top mass limit, by including corrections expressed as a power series in $m_{h}^{2}/m_{t}^{2}$. The authors found that the convergence of this expansion could be dramatically improved by factorising the leading order cross-section with exact top mass dependence. The results of this paper explain why this is so: indeed, we find that, to first subleading order in the threshold expansion, the NLO cross section can be completely expressed in terms of the leading-order cross section with shifted kinematics, and with full top mass dependence.

Returning to the large top mass limit, we can go further and consider triple Higgs production, a process for which analytic NLO corrections were presented recently in Ref. [49]. If we denote the invariant mass of the triple Higgs system by $M_{3h}$, and we set the $\overline{\rm MS}$ scale according to $\mu^{2}=M_{3h}^{2}$, the NLO cross section in the gluon fusion channel can be written as [49]

[TABLE]

with

[TABLE]

where $C^{(\delta)}_{3h}$ can be read off from Ref. [49], and does not affect our arguments. Expanding to NLP in $(1-z)$, we may rewrite Eq. (6.10) as

[TABLE]

This result is indeed reproduced from Eq. (3.26): implementing the scale choice as in Eq. (4.21), one can rewrite Eq. (3.26) in the present case as

[TABLE]

extracting the finite part in the $\overline{\rm MS}$ scheme, one finds precise agreement with Eq. (6.12).

7 Vector boson pair production

In the preceding two sections, we have illustrated the application of our general expression for the NLO K-factor in the gluon channel. In order to verify and illustrate the quark result, Eq. (4.17), in a non-trivial case, it is instructive to see how known matrix elements in vector boson pair production can be reproduced. This calculation is also an important illustration of the fact that our prediction, while based on power counting in the threshold variable $z$, applies to the fully differential squared amplitude, and not only to the integrated cross section. As an example, we consider di-boson production,

[TABLE]

where $V$ is an electroweak gauge boson. Let us start with di-photon production, the amplitudes for which can be found in Ref. [59] up to NLO, in four spacetime dimensions. The squared LO matrix element for this process, summed and averaged over colours and spins, is given by

[TABLE]

where $e_{q}$ is the electromagnetic charge of the quark. One may now consider the radiation of an additional gluon, leading to the process

[TABLE]

Computing the di-photon cross section at NLO, one must include squared Born matrix element for the process in Eq. (7.3), again summed and averaged over colours and spins. It is given by [59]

[TABLE]

Our aim is to show that, up to NLP accuracy, this matrix element can be obtained by shifting kinematics in the LO squared amplitude, as dictated by Eq. (4.16). To this end, one may rescale the gluon momentum as $k\to\lambda k$ in Eq. (7.4), and expand to next-to-leading power in $\lambda$, before setting $\lambda\to 1$. Next, one can parametrise the momenta in the centre of mass frame of the $VV$ system, as (see for example Ref. [60])

[TABLE]

where

[TABLE]

and we have introduced the invariants

[TABLE]

Note that in these notations the threshold variable is $\xi=1-z=1-s_{2}/s$. With these definitions, we can now expand to first subleading power in the gluon energy, and set $m_{V}=0$ for the di-photon case. The result is

[TABLE]

By performing the same procedure, one may easily show that the LO amplitude of Eq. (7.2), evaluated with the kinematic shifts defined in Eq. (3.13), yields

[TABLE]

We therefore see that Eq. (7.8) explicitly confirms the expectations raised by Eq. (4.16).

A similar exercise may be carried out for $W^{+}W^{-}$ production: the NLO squared amplitudes for this process (again in four spacetime dimensions) may be found in Ref. [60]. The Born-level squared matrix element, summed and averaged over colours and spins, is given by

[TABLE]

where $c^{tt}_{i}$, $c^{ss}_{i}$ and $c^{ts}_{i}$ are coefficients associated with $t$-channel, $s$-channel and interference graphs, respectively, for a quark of flavour $i$, and can be found in Ref. [60]. The remaining functions of Mandelstam invariants are given by

[TABLE]

Similarly, the NLO squared matrix element including the radiation of a gluon, for the $q\bar{q}$ initial state, also summed and averaged over colours and spins, is given by

[TABLE]

where the functions $\widehat{X}_{i}$, $\widehat{Y}_{i}$ and $\widehat{Z}_{i}$ have been obtained by rescaling the corresponding functions $X_{i}$, $Y_{i}$ and $Z_{i}$, given in Appendix D of Ref. [60], extracting the singular prefactor $-4s/(t_{k}u_{k})$. Parametrising momenta as in Eq. (7.5), one may expand each function up to NLP in $\xi=(1-z)$, using the same procedure as outlined above for the di-photon matrix element. Introducing the notation $\rho=m_{W}^{2}/s$, the results can be written as

[TABLE]

[TABLE]

[TABLE]

We have explicitly checked that the same results are obtained from Eq. (4.16), where for the right-hand side one must use the LO squared amplitude given in Eq. (7.10), with momenta shifted according to Eq. (3.13). This is a highly non-trivial cross-check: comparing Eqs. (7 - 7.15) with the di-photon case given in Eq. (7.8), one sees that the $WW$ case involves a much more complicated dependence on the opening angle $\theta_{1}$, and the partonic centre of mass energy $s$. A similar analysis could be carried out for $ZZ$ production [72], and also to provide analytic information for triple vector boson production, numerical results for which have been presented in refs. [73, 74].

8 Conclusion

In this paper, we have considered the hadro-production of an arbitrary heavy colourless system, in both the gluon-fusion and quark-antiquark-annihilation channels, near partonic threshold for the production of the selected final state. Our starting point is an all-order factorisation formula for the relevant scattering amplitudes, introduced in Refs. [23, 24], given here in Eq. (2), and valid to next-to-leading power in the threshold expansion. Specialising this formula to NLO in the strong coupling, we have observed how the general expression simplifies, and takes the form of a next-to-soft theorem, as derived for example in [43, 44, 45]. This simple expression, in turn, leads to a universal form for the transition probability, completely differential in the final state variables, and proportional to the Born-level transition probability, computed with a specific shift for the initial parton momenta. The result is the same for quarks and gluons (up to a trivial substitution of color factors), and is reported in Eqs. (3.14, 4.16). When the transition probability is integrated over final state variables, one finds that the inclusive cross section for the selected process can also be written in a simple and universal factorised form, given here in Eqs. (3.26, 4.17). More precisely, at NLO, and up to next-to-leading power in the threshold expansion, the cross-sections can be written as a universal $K$-factor, multiplying the leading order cross-section with a shifted partonic centre-of-mass energy. All these results apply regardless of whether the leading order process is tree-level or loop-induced, and the resulting K-factors are independent of hard scales such as heavy quark masses.

We have checked our results by reproducing known expressions for the production of up to three Higgs bosons at NLO, in the large top mass limit. Away from this limit, our formula provides new analytic information in the case of Higgs pair production at NLO, where only numerical results are presently known. Furthermore, we explain the observation, made previously in Ref. [55], that the convergence of the large top mass expansion can be improved in this process by factoring out the LO cross-section with exact top mass dependence. In the quark channel, we have shown how our formula is consistent with previous results for the production of photon and $W$ boson pairs, again at the level of differential distributions.

The results we have presented show that, to NLP accuracy, differential and inclusive NLO cross sections for colour-singlet final states are dramatically simpler than exact results, and we expect that they will be very easy to implement in numerical codes, providing checks of existing calculations, and improved approximations for differential distributions when complete results are not available, as is the case for loop induced processes with multi-particle electroweak final states. A detailed phenomenological analysis, including a study of the accuracy of the NLP approximation in different processes and kinematic domains, has been left to future work.

We emphasise that the simple universal expressions that we find at NLO can be systematically improved upon by relying on Eq. (2): in particular, a NLO calculation of the radiative soft function and of the radiative jet functions for quarks and gluons, which is under way, will lead to NLP approximations for cross sections of the type studied in this paper at NNLO level, exploring uncharted territory, in particular for processes which are loop induced and have so far been studied predominantly in the context of effective field theory approximations. We also note that the universal expressions we have derived for inclusive cross sections are also applicable to colourless final states arising beyond the Standard Model, in which case our results can provide a controlled approximation to estimate the size of higher-order corrections in selected models without the need to perform calculations that would be expensive for loop-induced processes.

Acknowledgments

We thank Fabrizio Caola, Jort Sinninghe-Damsté, Falko Dulat, Bernhard Mistlberger and Adrian Signer for discussions. We also thank Michael Spira for kindly providing us with personal notes relevant to our results in Section 3. This research was supported by the Research Executive Agency (REA) of the European Union through the contract PITN-GA-2012-316704 (HIGGSTOOLS). EL was supported by the Netherlands Foundation for Fundamental Research of Matter (FOM) programme 156, “Higgs as Probe and Portal”, and by the National Organisation for Scientific Research (NWO). CDW is supported by the UK Science and Technology Facilities Council (STFC), and thanks ETH, Zurich and the Higgs Centre for Theoretical Physics, University of Edinburgh, for hospitality. LV is supported by the People Programme (Marie Curie Actions) of the European Union Horizon 2020 Framework H2020-MSCA-IF-2014, under REA grant N. 656463 (Soft Gluons).