Two-loop helicity amplitudes for the production of two off-shell electroweak bosons in quark-antiquark collisions

Knowledge of two-loop QCD amplitudes for processes qq′→V1V2→l1l¯1′l2l¯2′\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ qq^{\prime}\to {V}_1{V}_2\to \left({l}_1{\overline{l}}_1^{\prime}\right)\left({l}_2{\overline{l}}_2^{\prime}\right) $$\end{document} is important for improving the theoretical description of four-lepton production in hadron collisions. In this paper we compute these helicity amplitudes for all intermediate vector bosons, V1V2 = γ*γ*, W+W−, ZZ, W±Z, W±γ*, including off-shell effects and decays to leptons.


I. INTRODUCTION
Production of vector boson pairs at the LHC is an interesting process for a variety of reasons.For definiteness, let us consider the W + W − final state.In that case, ATLAS and CMS collaborations have recently observed that the W + W − production cross section is about twenty percent higher than existing theoretical predictions [1][2][3].This observation prompted speculations about the possibility to explain this excess by physics beyond the Standard Model [4] and, at the same time, strongly emphasized the need to improve predictions for W + W − production within the Standard Model itself [5,6].In addition, production of W + W − pairs is an important process for studying anomalous couplings of electroweak gauge boson.Although current limits are already quite impressive [7], it is clear that studies of anomalous gauge boson couplings will intensify once the LHC Run II is underway.Making use of higher experimental precision will require improved modeling of W + W − production in the Standard Model.Finally, the pp → W W * process with one W -boson on the mass shell and the other off the mass-shell, is an important background to Higgs boson production in pp → H → W W channel. Better understanding of this background should allow improved measurements of Higgs boson couplings to W -bosons -including the anomalous ones -in the next run of the LHC.Similar arguments can be given for processes with other vector bosons in the final state.
It follows from these examples, that higher theoretical accuracy for vector-boson pair production in hadron collisions is essential.It can be achieved by extending existing computations of cross sections and kinematic distributions of pp → V 1 V 2 processes [5,6,11,12,[14][15][16] to next-to-next-to-leading order (NNLO) in perturbative QCD. 1 In particular, predictions for fiducial volume cross sections, where kinematic restrictions on final state particles are taken into account exactly, are crucial.A NNLO QCD prediction for pp → V 1 V 2 needs three ingredients: i) real-emission matrix elements for and q q → V 1 V 2 q q; ii) one-loop matrix elements for q q → V 1 V 2 g and qg → V 1 V 2 q and, finally, iii) two-loop amplitudes for tree-level process q q′ → V 1 V 2 .Once these three ingredients become available, they need to be put together in a self-consistent manner using existing methods for fully differential NNLO computations [17].For the V 1 V 2 production, the major unknown is the two-loop amplitude for q q′ → V 1 V 2 ; the goal of this paper is to provide it.
The remainder of the paper is organized as follows.In Section II we explain the general setup of the calculation.Since two-loop scalar master integrals required for this computation have been recently computed for equal electroweak boson masses in Refs.[18,19] and for unequal masses in Refs.[20,21], we mainly focus on the procedure that allows us to express the various contributions to scattering amplitudes in terms of these integrals.In Section III we describe checks on our computation, evaluate the scattering amplitudes numerically and discuss numerical stability of the results.We conclude in Section IV.

II. THE SETUP OF THE COMPUTATION
We are interested in the production of a four-lepton final states in proton collisions pp → (l 1 l′ 1 )(l 2 l′ 2 ).The four-lepton final states can be produced in three distinct ways: • through the production of a pair of off-shell vector bosons, pp • through the production of a single vector boson that decays to two vector bosons that, in turn, decay to lepton pairs, pp • through the production of a single vector boson that decays into a pair of leptons (one of them of-shell), followed by the emission of a vector boson by an of-shell lepton, e.g.
Each of these three processes depends on different combinations of electroweak couplings; therefore, they are separately invariant under QCD gauge transformations and we can compute QCD corrections to each of them separately.We note in this respect that QCD corrections to processes mediated by single gauge-boson production are simple since they are directly related to QCD corrections to the quark form factor of the vector current F V (s) whose perturbative expansion through NNLO is well-known [22].
However, computation of NNLO QCD corrections to process pp , where vector bosons couple only to fermions and not to other vector bosons or primary leptons, is non-trivial.
Calculation of this contribution to two-loop QCD amplitudes for four-lepton production processes is the main focus of this paper.
We consider quark-antiquark annihilation in to vector bosons and work in the approximation where quarks of the first two generations are massless and quarks of the third generation are consistently neglected.We also set the CKM matrix to an identity matrix.
Since we work with massless quarks, helicity is a conserved quantum number; therefore, once the helicity of the incoming quark is specified, the helicity of the incoming anti-quark is completely fixed.We will use this observation when writing the amplitude for quark-antiquark annihilation process in Eq.( 1).
The partonic process in Eq.( 1) can proceed in two different ways since the vector bosons can either couple directly to external fermions q q′ or to closed loops of virtual fermions, see Fig. 1.
Consequently, we write the scattering amplitude for process in Eq.(1) as where g W = e/ sin θ W is the SU (2) weak coupling, is the V i -boson propagator, λ q , λ 5 , λ 7 are helicities of the incoming quark and outgoing leptons, respectively, 1 are helicity-dependent couplings of vector bosons to quarks and leptons, and ǫ 3,4 are matrix elements for leptonic decays of V 1 and V 2 that we will specify shortly.The amplitudes A (d)   describe direct coupling of the vector bosons to external fermions and the amplitude A ng describes contributions of diagrams where vector bosons couple to loops of virtual fermions, see Fig. 1.The factor C ng V 1 V 2 involves sums over couplings of virtual fermions to gauge bosons.We note that full amplitude for q q′ → V 1 V 2 in Eq.( 2) is written as the sum of two A (d) amplitudes where two vector bosons appear in different order.These amplitudes are, therefore, quantities similar to ordered or primitive amplitudes often used in perturbative QCD computations.For us they are useful because we only need to compute one of them and then obtain the other one from p 3 ↔ p 4 permutation.We stress that amplitudes A (d) and A (ng) in Eq.(2) are computed assuming that vector bosons couple to a single quark generation and that the V qq coupling is vector-like with unit coefficient, i.e. γ µ .As we explain below, such amplitudes are sufficient to obtain physical amplitudes for pair production of all electroweak gauge bosons.
We begin by discussing the parametrization of A (d) µν ǫ µ 3 ǫ ν 4 amplitude in terms of Lorentz-invariant form factors.A representative diagram that contributes to A (d) at two loops in perturbative QCD is shown in Fig. 1.Our goal is to re-write these diagrams in such a way that all Feynman integrals can be dealt with using the integration-by-parts technique [23].To achieve this goal, we need to express the amplitude in terms of invariant form factors.To this end, we note that the most general µν is where γ [ν p⊥ γ µ] = γ ν p⊥ γ µ − γ µ p⊥ γ ν and the vector p ⊥ is defined by the Sudakov decomposition of the momenta p 3,4 The transverse momentum p ⊥ is orthogonal to p 1,2 , p ⊥ • p 1,2 = 0.The coefficients α 3,4 and β 3,4 in Eq.( 4) can be written as where we defined p 2 3 = m 2 3 , p 2 4 = m 2 4 and introduced standard Mandelstam variables s d) , i = 1, ..., 4, introduced in Eq.(3), depend on momenta and Lorentz indices.To make this dependence explicit, we decompose them into invariant form factors T i , i = 1, .., 17, We note that not all invariant form factors that appear in Eq.( 6) give independent contributions to physical amplitudes.This happens because we did not use the transversality condition for lepton currents p 3 • ǫ 3 = p 4 • ǫ 4 = 0.As we will see shortly, when physical amplitudes are computed, the number of relevant form factors will be reduced thanks to the transversality condition.
To calculate the physical amplitude in Eq.(2), we need to contract µν with polarization vectors of external vector bosons.As we already mentioned, they are given by matrix elements of the vector current between relevant leptonic states.Using the spinor-helicity formalism, we write Since, as we already pointed out, the helicity of the incoming quark fully determines the allowed helicity of the incoming antiquark, we need eight helicity amplitudes to fully describe production of two vector bosons.They are We note that, according to Eq.( 7), a change in lepton helicities can be obtained by interchanging momenta p 5 ↔ p 6 and p 7 ↔ p 8 , where necessary.Therefore, we obtain all required helicity amplitudes from M(q L , 5 L , 7 L ) and M(q R , 5 L , 7 L ) by simple permutations of lepton momenta.
For left-and right-handed incoming quarks we find where the following nine combinations of form factors enter We note that form factors only depend on s, t, u and m 2 3 , m 2 4 ; therefore, they are invariant under p 5 ↔ p 6 and p 7 ↔ p 8 permutations.For this reason, it is straightforward to obtain all the relevant amplitudes from Eq. (8).To find T 1,..,17 and F 1,2,..15 , we need to construct projection operators.To accomplish this, we write a generic amplitude as a matrix element of a string of Dirac matrices Multiplying it with ūp 1 Ôv p 2 and summing over spinor helicities, we find Choosing different operators O, we can project on individual T form factors or their combinations.
Below we list all the projection operators and the results that we get when these operators are convoluted with the amplitude Each of the G 1..17 projections can be calculated from Feynman diagrams by taking traces and integrating over loop momenta.These quantities depend on the Mandelstam invariants of the underlying 2 → 2 process q q′ → V 1 V 2 but not on the polarization vectors of electroweak bosons or on the momenta of final-state leptons.Thus, they can be expressed in terms of Feynman integrals of the type introduced in Refs.[20,21] and then reduced to master integrals using integration-by-parts identities.
We can use projections G 1,..,17 to derive the T form factors; the result reads The F form factors that are used for the evaluation of the amplitude can now be easily computed using Eq.( 9).We do not present these results here since they are not very illuminating.
One point worth emphasizing though is that in Eq.( 13) there are expressions (see e.g.equation for T 17 ) that appear to have spurious singularities in the limit d → 4. It would have been unfortunate if these singularities survive in the final formulas for form factors F since, in such a situation, computation of all pieces needed for the evaluation of G's, including master integrals, to higher orders in ǫ is required.It is therefore pleasing to observe that this does not happen and once results for the F form factors are written in terms of G 1,..,17 projectors, all spurious (d − 4) singularities disappear.
For future reference, we give results for leading order form factors that appear in the physical amplitudes Eq.( 8) We stress that these results are exact in a sense that no d → 4 limit was taken to obtain them; in other words, all the d-dependence cancels out completely once physical form factors are computed.
Before proceeding to the discussion of various couplings of vector bosons to fermions introduced in Eq.(2), we want to explain why A (d) and A (ng) amplitudes can be computed assuming that interactions of vector bosons with quarks are mediated by vector currents.The main reason is that we only consider here contributions of massless quarks and that, in such a case, the helicity is conserved.Therefore, we can always write couplings of vector bosons to fermions as linear combinations of left-and right-handed couplings and then move helicity projection operators to external lines.In this way the helicity-dependent couplings are generated and amplitudes that remain can be viewed as originating from pure vector current interactions of gauge bosons and quarks.
Of course, this discussion applies only to diagrams where vector bosons couple directly to external quark lines; in the notation of Eq.(2) such diagrams contribute to A (d) amplitudes.However, we will now argue that the same parametrization Eq.(3) can be used to compute amplitudes A (ng) , which receive contributions from diagrams where vector bosons couple to closed fermion loops.In fact, this would have been obvious provided that electroweak boson coupling to fermions is vector-like since in this case tracing over spin degrees of freedom in the internal quark loop gives us diagrams that are not very different from the ones that we already considered.Potential problems could be expected with axial couplings and, in particular, with terms where one axial and three vector couplings appear in ggV 1 V 2 Green's function.However, all such terms cancel because of C-parity conservation for massless fermions for any final state with two vector bosons [15].2Contributions of terms with two vector and two axial couplings are not anomalous and must be equal to those with four vector couplings.Hence, we conclude that it is sufficient to consider vector-current couplings of gauge bosons to internal fermion loops, to account for all non-vanishing contributions to the amplitude.The parametrization of the amplitude A (ng) is then taken from Eq.(3) and its expression through invariant form factors is taken from Eq. ( 8).
To complete our construction of the scattering amplitude for q q′ → V 1 V 2 in Eq.( 2), it remains to specify various helicity-dependent couplings that we introduced there.Below we present these couplings for various pairs of electroweak gauge bosons that can be produced in proton collisions.
A. γ * γ * production Photons are produced in the annihilation of a quark and an antiquark of the same flavor q = q ′ .Photon interactions with both quarks and leptons are pure vector-like and, therefore, helicityindependent.We find where Q l and Q q are electromagnetic quark and lepton charges in units of the positron charge.

Finally, C
γγ is given by the sum of up and down quark charges of the first two generations B.

ZZ production
Pairs of Z bosons are produced in annihilation of quarks and anti-quarks of the same flavor.Couplings of Z bosons to quarks and leptons depend on helicity and weak isospin of the corresponding particle.We find where The coefficient C ng ZZ is given by the sum of up and down quark vector and axial charges of the first two generations.We find C.
Zγ * production Production of Z and γ * occurs in annihilation of quarks and anti-quarks of the same flavor.The couplings to external fermions are given in Eqs.(15,17).The coefficient D.
W + γ * and for W − γ * production The Since W bosons interact with left-handed fermions, we have Coefficients C ng W + γ and C ng W − γ vanish identically due to electric charge conservation.

E.
W + Z and for W − Z production Production processes for W ± Z are identical to W ± γ * .The relevant couplings are Coefficients C ng W ± Z vanish identically due to electric charge conservation. F.
Finally, we display the necessary couplings for the pair production of two W bosons.This process can occur thanks to uū and d d annihilation.First, consider the process u(p 1 )ū(p 2 ) → W + (p 3 )W − (p 4 ).In this case, In case of d(p 1 ) d(p 2 ) → W + (p 3 )W − (p 4 ), the coupling constants read Couplings to leptons are given in Eq.( 22).The coefficient C (ng) W + W − receives contributions from two generations 3 and reads

III. CALCULATION OF THE FORM FACTORS AND THE AMPLITUDES
Once projection operators are established, one can compute the form factors and construct the scattering amplitudes for arbitrary di-boson final state.To this end, we note that all integrals that appear in the calculation of form factors can be associated with one of the six topologies, introduced in Refs.[20,21].Integrals that belong to each of these topologies are closed under the integrationby-parts identities [23].We use QGRAF [25] to generate Feynman diagrams and FORM [26] for algebraic manipulations and computation of projections G 1,..,17 .We use FIRE [27][28][29] to reduce all the integrals that appear in this calculation to master integrals.The master integrals were computed by us in Refs.[20,21].In principle, once all the ingredients are in place, computation of the amplitude becomes straightforward.In practice, however, it requires some effort to put all the pieces together primarily because algebraic expressions that appear e.g. in the course of the reduction to master integrals are quite large in size.The master integrals are expressed in terms of Goncharov polylogarithms up to weight four; we use GiNaC [30] implementation [31] to compute them.
The results for the amplitudes contain infra-red and ultraviolet divergences.The ultraviolet divergences are removed by the renormalization of the strong coupling constant.Since tree-level scattering amplitudes are independent of α s , we only need its renormalization through one loop.It reads where α (0) s is the bare and α s = α s (µ) is the renormalized coupling constant.In addition, β 0 = (11C A − 4T R n f )/6 = 11/2 − n f /3 is the QCD beta-function, n f is the number of massless flavors, C A = 3, T R = 1/2 and S ǫ = (4π) ǫ e −ǫγ E , with γ E being the Euler constant.
3 As we already mentioned several times, we do not consider the third quark generation in this paper.According to Catani [32], infra-red singularities of UV-renormalized amplitudes at next-to-next-toleading order are fully determined by leading and next-to-leading order amplitudes.This feature can be used as an important check of the correctness of the computation.To introduce Catani's result, we write a UV-renormalized amplitude as Leading order amplitude M 0 in this formula is finite and, in fact, d-independent.The next-toleading order amplitude contains infra-red divergences that, however, can be written in a factorized form For the process of vector boson pair production, the final state is neutral.Therefore, where L S = ln(−s/µ 2 − i0) = ln |s/µ 2 | − iπ and C F = 4/3.In variance with M 0 , M 1 and M 1,fin depend on the dimensional regularization parameter ǫ.This feature is important for proper comparison of Catani's formula with results of explicit computation.
The two-loop amplitude M 2 can be written in a similar way where The two constants that enter this formula read In Fig. 2 we show ratios of 1/ǫ singularities in the squared amplitude | Ã(d,2) LLL | 2 computed directly and using Catani's formula.Deviations of this ratio from one signal numerical instabilities.We ) in the bulk of the phase-space.Significant instabilities are observed for backward scattering (178 degrees), where the transverse momentum is close to 1 GeV.However, the situation improves considerably already for 176 degree scattering where the transverse momentum is 2 GeV.The forward scattering limit appears to be more stable; even at two degrees, the 1/ǫ contribution is computed properly to within a percent.
In the left pane of Fig. 3 we show the absolute value squared of the ratio of the finite part of the A (d,2) LLL amplitude and the leading order amplitude A (d,0) LLL .To understand numerical accuracy of these results, we compared the output obtained with the double-precision version of the Fortran code with the Mathematica implementation.The advantage of the latter is that it provides a possibility to compute amplitudes with arbitrary numerical precision thereby ameliorating the problem of numerical instability.For backward scattering, we find that up to 174 degrees, the finite part is computed to within a per mille.For forward direction, the situation is similar but, perhaps, slightly better.
Finally, in the right pane of Fig. 3  LLL .Numerical instabilities are apparent for backward scattering.In fact, at 170 degrees (p ⊥ ∼ 5 GeV), the agreement between double-precision Fortran code and the Mathematica code is about ten percent.In the forward direction, the situation is much better -the double-precision Fortran results agree with the results obtained using Mathematica implementation to better than a fraction of a percent for scattering angles as small as six degrees.
To summarize, while the two-loop amplitudes that we compute in this paper do exhibit numerical instabilities at small values of vector boson transverse momenta, we believe the stability is acceptable for phenomenological applications.Moreover, there are several ways to improve the situation.For example, it is possible to extend the Fortran code to provide results with quadruple precision.
Note that computation of master integrals with arbitrary precision is feasible since Goncharov polylogarithms implementation in GiNaC does provide this functionality [31].Moreover, it should also be possible to construct expansion of analytic expressions for scattering amplitudes that we obtained in this paper around singular limits, for example for forward or backward scattering and threshold production.If such expansions become available, computation of helicity amplitudes in singular limits will be significantly simplified.We leave these improvements for future work.

IV. CONCLUSIONS
In this paper we described computation of two-loop scattering amplitudes for the annihilation of a quark and an antiquark into four leptons, that occurs through the production of two electroweak gauge bosons.The invariant masses of gauge bosons are kept arbitrary.We have given explicit formulas for projection operators that allow one to compute contributions of individual Feynman diagrams to invariant form factors.We use these form factors to construct helicity amplitudes for vector boson pair production processes including all off-shell effects and leptonic decays of vector bosons.

Figure 1 :
Figure 1: Representative diagrams that contribute to two-loop amplitude for vector boson production in hadron collisions.
we show absolute value squared of the ratio of the finite part of the left amplitude A (ng ,2) LLL and the leading order LLL amplitude A (d,0)

Table I :
(34)887902 4.0059079 − i 3.0593301 0.5752861 + i 5.8440713 −9.6769949 + i 0.9775875 Ratios of selected one-loop helicity amplitudes and tree amplitudes, see Eq.(34).Momenta of external particles are given in the main text of the paper.Tree-amplitudes are A