Top quark mass effects in $gg\to ZZ$ at two loops and off-shell Higgs interference

TL;DR

This paper analyzes the effects of the top quark mass on the interference between Higgs production and background in gluon fusion to Z-boson pairs, crucial for Higgs width measurements.

Contribution

It introduces a method using conformal mapping and Padé approximants to accurately compute two-loop form factors with top mass effects in Higgs interference.

Findings

Accurate form factors for top-mass effects at two loops

Enhanced understanding of Higgs interference above top threshold

Improved theoretical predictions for Z-boson pair production

Abstract

We consider top-quark mass effects in the Higgs-interference contribution to $Z$-boson pair production in gluon fusion. While this production mechanism is formally of next-to-next-to leading order, its contribution is numerically important above the top threshold $M_{ZZ}^2=4m_t^2$. This region is essential to constrain the width of the Higgs boson and good control over the top-quark mass dependence is crucial. We determine the form factors that are relevant for the interference contribution at two-loop order using a method based on a conformal mapping and Pad\'e approximants constructed from the expansions of the amplitude for large top mass and around the top threshold.

Equations

$$|\mathcal{B}⟩=\frac{{\delta }^{AB}}{N_A}(p_1\cdot p_2{g}^{\mu \nu }-p_1^{\nu }{p}_2^{\mu })P_Z^{\alpha \rho }(p_3)P_{Z,\rho }^{\beta }(p_4)|B_{\mu \nu \alpha \beta }^{AB}⟩,$$

$$|\mathcal{B}⟩=\frac{i{g}_W^2}{4{\cos }^2{\theta }_W}\left(v_t^2|{\tilde{\mathcal{B}}}_{VV}⟩+a_t^2|{\tilde{\mathcal{B}}}_{AA}⟩\right),$$

$$|{\tilde{\mathcal{B}}}_i⟩=\frac{{\alpha }_s}{4\pi }|{\tilde{\mathcal{B}}}_i^{(1)}⟩+{\left(\frac{{\alpha }_s}{4\pi }\right)}^2|{\tilde{\mathcal{B}}}_i^{(2)}⟩+\dots ,$$

$$|{\tilde{\mathcal{F}}}_i^{(2)}⟩=|{\tilde{\mathcal{B}}}_i^{(2)}⟩+\frac{e^{ϵ{\gamma }_E}}{\Gamma (1-ϵ)}\left[\frac{2C_A}{{ϵ}^2}{\left(\frac{{\mu }^2}{-s}\right)}^{ϵ}+\frac{{\beta }_0}{ϵ}\right]|{\tilde{\mathcal{B}}}_i^{(1)}⟩,$$

$$|{\tilde{\mathcal{F}}}_i^{(1)}⟩$$

$$|{\tilde{\mathcal{F}}}_i^{(2)}⟩$$

$$|{\tilde{\mathcal{B}}}_{AA-VV}⟩\equiv |{\tilde{\mathcal{B}}}_{AA}⟩-|{\tilde{\mathcal{B}}}_{VV}⟩$$

$$z=\frac{4\omega }{(1+\omega {)}^2}$$

$$[n/m](\omega )=\frac{\sum _{i=0}^n{a}_i{\omega }^i}{1+\sum _{j=1}^m{b}_j{\omega }^j},$$

$$P_{AA-VV}^{(1)}(\omega )=$$

$$P_{AA-VV}^{(2)}(\omega )=$$

$$+s_{AA}^{(2)}(z(w))-s_{VV}^{(2)}(z(w)),$$

$$P_{VV}^{(1)}(\omega )=$$

$$P_{VV}^{(2)}(\omega )=$$

$$+s_{VV}^{(2)}(z(w)),$$

$$|{\tilde{F}}_{AA}^{(j)}(z(\omega ))⟩\simeq P_{AA-VV}^{(j)}(\omega )+P_{VV}^{(j)}(\omega ),$$

$$|{\tilde{F}}_{VV}^{(j)}(z(\omega ))⟩\simeq P_{VV}^{(j)}(\omega ).$$

$$m_Z=91.1876\text{ GeV},m_t=173\text{ GeV},$$

$$\text{Re}(z(\omega ))>0\text{and}|\omega |<1.2.$$

$$a_{R,i}\in [0.1,10],$$

$$|{\tilde{\mathcal{F}}}_i^{(1)}⟩$$

$$|{\tilde{\mathcal{F}}}_i^{(2)}⟩$$

$$[dl]=\frac{d^{d}l}{i{\pi }^{\frac{d}{2}}}{e}^{ϵ{\gamma }_E},$$

$$C_0=\int [dl]\frac{1}{l^2[(l+q{)}^2-1][(l+q-p_Z{)}^2-1]}$$

$$a_{AA}^{(3,0)}=$$

$$a_{AA}^{(5,0)}=$$

$$a_{AA}^{(7,0)}=$$

$$\displaystyle-788896\*r_{Z}^{6}+269568\*r_{Z}^{7}+16\*r_{p_{T}}\*(1-2\*r_{Z})^{2}\*(39-234\*r_{Z}+672\*r_{Z}^{2}-920\*r_{Z}^{3}+592\*r_{Z}^{4})\big{]}\,,$$

$$b_{AA}^{(2,1)}=$$

$$b_{AA}^{(3,0)}=$$

$$-64\text{*}{C}_0\text{*}(1-2\text{*}{r}_Z{)}^2\text{*}(1-4\text{*}{r}_Z{)}^2\text{*}(-1-7\text{*}{r}_Z+34\text{*}{r}_Z^2-44\text{*}{r}_Z^3+8\text{*}{r}_Z^4)$$

$$\displaystyle+64\*r_{Z}\*(1-4\*r_{Z})^{2}\*\sqrt{\frac{1-r_{Z}}{r_{Z}}}\*\arctan\bigg{(}\frac{2\*\sqrt{(1-r_{Z})\*r_{Z}}}{1-2\*r_{Z}}\bigg{)}\*(9-45\*r_{Z}+70\*r_{Z}^{2}-56\*r_{Z}^{3}+32\*r_{Z}^{4})$$

$$\displaystyle-(1-2\*r_{Z})\*(1-4\*r_{Z})\*\big{[}192+9\*\pi^{2}-56\*\ln(2)-1728\*r_{Z}-90\*\pi^{2}\*r_{Z}+560\*\ln(2)\*r_{Z}+6256\*r_{Z}^{2}+384\*\pi^{2}\*r_{Z}^{2}$$

$$\displaystyle\qquad-2016\*\ln(2)\*r_{Z}^{2}-12480\*r_{Z}^{3}-816\*\pi^{2}\*r_{Z}^{3}+3584\*\ln(2)\*r_{Z}^{3}+12032\*r_{Z}^{4}+576\*\pi^{2}\*r_{Z}^{4}-3584\*\ln(2)\*r_{Z}^{4}-2048\*r_{Z}^{5}\big{]}$$

$$\displaystyle+128\*(1-2\*r_{Z})\*\ln(2-4\*r_{Z})\*(2-23\*r_{Z}+98\*r_{Z}^{2}-184\*r_{Z}^{3}+152\*r_{Z}^{4}-112\*r_{Z}^{5}+96\*r_{Z}^{6})\bigg{]}\,,$$

$$b_{AA,\ln }^{(3,0)}=$$

$$b_{AA}^{(3,1)}=$$

$$b_{AA}^{(4,1)}=$$

$$b_{AA}^{(5,0)}=$$

$$+109350\text{*}{\pi }^2\text{*}{r}_Z+1881360\text{*}\ln (2)\text{*}{r}_Z+7528432\text{*}{r}_Z^2-855900\text{*}{\pi }^2\text{*}{r}_Z^2-14848320\text{*}\ln (2)\text{*}{r}_Z^2-44282176\text{*}{r}_Z^3$$

$$+3817800\text{*}{\pi }^2\text{*}{r}_Z^3+67172160\text{*}\ln (2)\text{*}{r}_Z^3+141881152\text{*}{r}_Z^4-10558080\text{*}{\pi }^2\text{*}{r}_Z^4-187708800\text{*}\ln (2)\text{*}{r}_Z^4$$

$$-245787136\text{*}{r}_Z^5+18135360\text{*}{\pi }^2\text{*}{r}_Z^5+319941120\text{*}\ln (2)\text{*}{r}_Z^5+202149888\text{*}{r}_Z^6-17763840\text{*}{\pi }^2\text{*}{r}_Z^6$$

$$\displaystyle\qquad-302008320\*\ln(2)\*r_{Z}^{6}-35688448\*r_{Z}^{7}+7326720\*\pi^{2}\*r_{Z}^{7}+117350400\*\ln(2)\*r_{Z}^{7}-25067520\*r_{Z}^{8}\big{]}$$

Full text

I Introduction

A direct measurement of the Higgs boson width $\Gamma_{H}$ is not possible at the LHC or even the envisioned next generation of collider experiments. However, indirect constraints can be obtained at the LHC by studying the process $pp\to H\to ZZ\,(\to 4l)$ on the Higgs boson peak where the cross section depends on the combination $g_{Hgg}^{2}g_{HZZ}^{2}/\Gamma_{H}$ and off the peak where the measurement of the cross section constrains the product $g_{Hgg}^{2}g_{HZZ}^{2}$ of the effective Higgs boson-gluon coupling $g_{Hgg}$ and the Higgs boson-$Z$ boson coupling $g_{HZZ}$, as proposed in Kauer and Passarino (2012); Caola and Melnikov (2013); Campbell et al. (2014a).444Note that the indirect way of constraining the Higgs width is not entirely model-independent Englert and Spannowsky (2014); Englert et al. (2015). The same strategy can be employed with $WW$ final states Campbell et al. (2014b). The latest studies from the LHC experiments give an upper limit of 14.4 MeV at 95% C.L. from the $ZZ$ final state at ATLAS Aaboud et al. (2018) and the value $3.2^{+2.8}_{-2.2}$ MeV from the combination of $VV$ final states in CMS Sirunyan et al. (2019), close to the SM prediction $\Gamma_{H}^{\text{SM}}=4.10\pm 0.06$ MeV de Florian et al. (2016). Measurements of the Higgs boson signal at large invariant mass can also be used to directly constrain physics beyond the Standard Model in the Higgs sector Gainer et al. (2015); Azatov et al. (2015); Buschmann et al. (2015); Azatov et al. (2016).

Here, we focus on the loop-induced continuum gluon fusion process $gg\to ZZ$ and in particular its interference with the off-shell Higgs contribution $gg\to H^{*}\to ZZ$. Despite the narrow width of the Higgs boson these interference effects are sizable with 10% of the Higgs signal stemming from the off-shell region where the invariant mass of the two decay products is greater than $2\,m_{Z}$ Kauer and Passarino (2012) and higher-order corrections are required to control the uncertainties. The Higgs-mediated amplitude only depends on two scales, the mass $m_{q}$ of the quark in the loop and the invariant mass $M_{ZZ}$ of the final state. Next-to-leading order (NLO) corrections with the full quark-mass dependence have been known for some time Spira et al. (1995); Harlander and Kant (2005); Anastasiou et al. (2007); Aglietti et al. (2007), and the top-quark mass dependence at next-to-next-to leading order (NNLO) has been reconstructed very recently Davies et al. (2019a) (see also Harlander et al. (2019)). On the other hand the continuum amplitude depends on four scales $m_{q}$, $m_{Z}$, $M_{ZZ}$ and the transverse momentum $p_{T}$ of one of the $Z$ bosons, and the exact result is only known at leading order (LO) Glover and van der Bij (1989) while an analytic NLO calculation appears extremely challenging. In the massless limit $m_{q}=0$ the two-loop amplitude has been determined in Caola et al. (2015a); von Manteuffel and Tancredi (2015); Caola et al. (2015b) and the NLO cross section in Caola et al. (2016). Recently, also the quark–gluon channel has been included Grazzini et al. (2019).

The contribution from top quarks at two-loop order has been computed in a large-mass expansion (LME) Melnikov and Dowling (2015); Campbell et al. (2016); Caola et al. (2016) and is known up to $1/m_{t}^{12}$. While the contribution from massless quarks dominates the interference correction at small invariant masses, the top-quark contribution is of the same size near the top threshold $M_{ZZ}=4m_{t}^{2}$ and dominates in the large invariant-mass regime. Since the LME ceases to provide a reliable description above the top threshold, the authors of Campbell et al. (2016) have improved their prediction by a conformal mapping and the construction of Padé approximants based on the available number of LME coefficients. In Gröber et al. (2018) we have extended this method by considering the expansion around the top threshold in addition to the LME and demonstrated that the top-mass effects can be reproduced correctly by comparing results for the two-loop amplitude for $gg\to HH$ with the numerical calculation from Borowka et al. (2016a, b); Heinrich et al. (2017).555Recently, an independent numerical calculation Baglio et al. (2019), several approximations Bonciani et al. (2018); Davies et al. (2018, 2019b); Xu and Yang (2019) which are consistent with the earlier results and a combined result Davies et al. (2019c) have appeared.

In this work we consider the form factors of the continuum $gg\to ZZ$ amplitude that are relevant for the interference contribution at one and two loops. The non-analytic terms in the expansion around the top threshold are computed up to at least order $(1-z)^{4}$, where $z=M_{ZZ}^{2}/(4m_{t}^{2})+i0$, and used to construct Padé approximants. Together with the exactly known real NLO top quark Campbell et al. (2014c); Campbell et al. (2016) and the massless quark corrections Caola et al. (2015a); von Manteuffel and Tancredi (2015); Caola et al. (2015b, 2016); Grazzini et al. (2019) this is sufficient to determine the full NLO interference contribution with realistic top-quark mass dependence.

II Form factors for interference

Up to the two loop level, the amplitude for the top-mediated non-resonant continuum production process $g(\mu,A,p_{1})+g(\nu,B,p_{2})\to Z(\alpha,p_{3})+Z(\beta,p_{4})$ receives contributions from both box and double-triangle diagrams, see figure 1. The latter are known for arbitrary quark masses Campbell et al. (2016, 2007) and will not be discussed in the following.

The box amplitude $\Ket{B_{\mu\nu\alpha\beta}^{AB}}$ has a complicated tensor structure Glover and van der Bij (1989); Costantini et al. (1971); Caola et al. (2015a); von Manteuffel and Tancredi (2015). However, the interference with the Higgs-mediated amplitude is described by a single form factor. Adopting the conventions of Campbell et al. (2016) it takes the form

[TABLE]

with $N_{A}=N_{c}^{2}-1$ and $P_{Z}^{\alpha\rho}(p)=-g^{\alpha\rho}+p^{\alpha}p^{\rho}/m_{Z}^{2}$. The form factor can be decomposed into a vector and axial-vector part

[TABLE]

where $a_{t}=1/2$ and $v_{t}=1/2-4/3\sin^{2}\theta_{W}$ denote the axial-vector and vector couplings for an up-type quark. Mixed $v_{t}a_{t}$ terms are forbidden by charge conjugation symmetry. The order in the strong coupling constant $\alpha_{s}$ is indicated as follows

[TABLE]

with $i=VV,AA$. At order $\alpha_{s}^{2}$ the renormalized form factors contain IR divergences, which cancel in the combination with real corrections, and we define the finite remainder by applying the subtraction Catani (1998)666Note that this subtraction differs at order $\epsilon^{0}$ from the one given in eq. (2.14) of Campbell et al. (2016).

[TABLE]

where $\beta_{0}=\tfrac{11}{3}C_{A}-\tfrac{4}{3}T_{f}n_{l},C_{A}=3,T_{f}=\tfrac{1}{2},n_{l}=5$, and the form factors $\Ket{\widetilde{\mathcal{B}}_{i}^{(1,2)}}$ are defined in $d=4-2\epsilon$ dimensions. The one-loop form-factors $\Ket{\widetilde{\mathcal{B}}_{i}^{(1)}}$ are already finite; we define $\Ket{\widetilde{\mathcal{F}}_{i}^{(1)}}=\Ket{\widetilde{\mathcal{B}}_{i}^{(1)}}$ for the sake of a consistent notation.

II.1 The amplitude near threshold

Above the top threshold at $z=1$ the top quarks in the loop can go on shell which manifests as non-analytic terms in the expansion of the form factors in $\bar{z}\equiv 1-z$, generating a sizable imaginary part. As shown in Gröber et al. (2018) the knowledge of these terms alone provides very valuable information for the determination of top-quark mass effects in our approach. The calculation of the non-analytic terms is significantly simpler than that of the analytic contributions and was described in detail in Gröber et al. (2018) for the three leading non-analytic expansion terms of the one and two-loop form factors for $gg\to HH$. For $gg\to ZZ$ we expand the amplitude up to high orders in $\bar{z}\equiv 1-z$ and therefore use the expansion by regions Beneke and Smirnov (1998); Jantzen (2011) to expand the full-theory diagrams instead of an EFT approach where a large number of effective vertices is required due to the deep expansion. We use QGRAF Nogueira (1993) to generate the Feynman diagrams which are processed and expanded using private FORM Vermaseren (2000) code. The IBP reduction Chetyrkin and Tkachov (1981) is performed with FIRE Smirnov (2015) which is based on the Laporta algorithm Laporta (2000).

Our results are given in Appendix A and an ancillary Mathematica file. They are of the form

[TABLE]

where $\bar{n}_{2}$ is $n$ modulo 2, the coefficients are functions of the dimensionless variables $r_{Z}=m_{Z}^{2}/M_{ZZ}^{2}$ and $\tilde{x}=(p_{T}^{2}+m_{Z}^{2})/M_{ZZ}^{2}$. We use the symbol $\asymp$ to indicate that terms which are analytic in $\bar{z}$ and currently unknown have been dropped on the right-hand side.

Threshold logarithms $\ln\bar{z}$ and logarithms $\ln(-4z)$ related to massless cuts in the amplitude first appear at two-loop order. While we generally compute the expansion coefficients up to $n=8$, i.e. expand up to $\bar{z}^{4}$, we find that for the massless-cut contribution proportional to $\ln(-4z)$ more input is required to achieve a reliable Padé approximation. We therefore compute the corresponding coefficients $b_{i,\ln}^{(n,m)}$ up to $n=9$.

As in Higgs pair production there is no S-wave contribution to the form factors relevant for the interference and the leading non-analytic terms involve the $\bar{z}$-suppressed P-wave Green function Beneke et al. (2014).

II.2 Behavior for $z\to\infty$

In addition to the LME and threshold expansions we can exploit scaling information in the small-mass limit $m_{t}\to 0$ which corresponds to $z\to\infty$. This does not require an additional calculation in this region but relies solely on the symmetries of QCD. The absence of infrared $1/m_{t}$ power divergences as $m_{t}\to 0$ implies that the form factors can only show logarithmic behavior as $z\to\infty$. Below we show that the difference

[TABLE]

vanishes as $z\to\infty$. To prove this we note that chirality is conserved in massless QCD and hence the four-point correlator of two vector currents, a left-handed and a right-handed current, which we denote in short by [V,V,V-A,V+A], vanishes in the limit of zero quark masses.777To make the double-triangle contribution shown in Fig. 1 anomaly free we have to consider doublets of quarks and not just a single (top) quark, but we omit this technicality here since the double-triangle contribution is known and not considered below. Using that the correlator [V,V,V,A] vanishes due to charge conjugation we immediately conclude that [V,V,V,V]-[V,V,A,A]$\to 0$ as $z\to\infty$. We exploit this below and reconstruct the top-mass dependence of $\Ket{\widetilde{\mathcal{F}}_{VV}^{(i)}}$ and $\Ket{\widetilde{\mathcal{F}}_{AA}^{(i)}}-\Ket{\widetilde{\mathcal{F}}_{VV}^{(i)}}$ where we have one additional condition for the latter.

III The method

We approximate the box form factors (2) using our approach from Gröber et al. (2018). First, we introduce subtraction functions $s_{VV}^{(2)},s_{AA}^{(2)}$ in such a way that the combinations $\Ket{\widetilde{\mathcal{F}}_{i}^{(2)}}-s_{i}^{(2)}$ retain their analytic structure for $|z|<1$ but have threshold expansions which are free of logarithms $\ln(\bar{z})$ up to the highest known order, i.e. up to $\bar{z}^{4}$. The construction of such subtraction functions is detailed in Gröber et al. (2018) and we give the ones we explicitly need in Appendix B. Note that even after this subtraction the threshold and large mass expansions of the two-loop form factors still receive contributions proportional to a single logarithm $L_{s}\equiv\ln(-4z)$ from diagrams with massless cuts. We therefore split the subtracted two-loop form factors into a constant and a logarithmic part and construct separate approximants for each part.

The top mass dependence is contained in the variable $z$ and the conformal transformation Fleischer and Tarasov (1994)

[TABLE]

is used to map the entire complex $z$ plane onto the unit disc $|\omega|\leq 1$ with the branch cut for $z\geq 1$ corresponding to the perimeter. Thus, the top-mass dependence is encoded by a function that is analytic in the region $|\omega|<1$ and can be reconstructed using Padé approximants

[TABLE]

where the $n+m+1$ coefficients $a_{i},b_{j}$ can be fixed by imposing the condition that the expansion of eq. (8) in the LME and threshold region must reproduce the known coefficients for given, fixed values of $r_{Z}$ and $\tilde{x}$. The small-mass behavior discussed in sec. II.2 is not used to further constrain the Padé coefficients, but is taken into account by a rescaling of the Padé ansatz. Hence, we use approximation functions of the form

[TABLE]

where $P_{AA-VV}^{(j)}$ is used to approximate the difference between the axial-vector and vector form factors, whereas the vector form factors in isolation are approximated using $P_{VV}^{(j)}$. The limit $z\to\infty$ corresponds to $\omega\to-1$ where the approximants in eq. (8) approach a constant value. Thus, the rescaling eq. (9) enforces the correct asymptotic behavior for $z\to\infty$ discussed in sec. II.2 and provides us with free parameters $a_{R,i}$ that can be varied in addition to the polynomial degrees $n,m,k$ and $l$ to assess the stability of the approximation. We note that these variations are performed independently for all the terms in eq. (9). Our final ansätze for the form factor approximation are then

[TABLE]

IV Results

Before showing our results at NLO for the form factors, we can compare the LO form factors constructed as discussed in the previous sections with the full analytic result. We choose as input for the on-shell $Z$-boson and top quark masses

[TABLE]

and show results for two different values of $\tilde{x}$ in fig. 2. The plots contain the maximum information we have available from the LME at LO (see Campbell et al. (2016)) and our threshold expansion. By construction the Padé ansatz in eq. (8) contains poles in the complex $\omega$ plane whereas the functions it approximates are analytic in $z$ implying the absence of poles in the unit disc $|\omega|\leq 1$. Furthermore poles in the vicinity of the unit disc can cause unphysical behavior in the reconstructed form factors. We were not able to construct Padé approximants without poles inside a larger disc $|\omega|\leq 1.2$. Therefore we focus on the time-like region of the form factors and construct only Padé approximants which do not contain poles for

[TABLE]

We obtain an uncertainty estimate for our results in the following way. For every phase space point, we calculate the mean and standard deviation for each contributing Padé approximant in eq. (9). To this end, we vary the rescaling parameters $a_{R,i}$ in the region

[TABLE]

and vary $[n/m]$ within $|n-m|\leq 3$, where $n+m+1$ is the number of available constraints. We construct 100 variants for each Padé approximant. Our final prediction then follows from the sum of the mean values of the Padé approximants, with an uncertainty obtained by adding the individual errors in quadrature.

Fig. 2 shows the Padé approximants from eq. (10) and eq. (11) for the LO form factors $\Ket{\widetilde{F}_{VV}^{(1)}}$ and $\Ket{\widetilde{F}_{AA}^{(1)}}$ including our uncertainty estimate as points with error bars. We observe good agreement with the full results, which are indicated by the solid lines, up to large values of the invariant mass $M_{ZZ}$ of the $Z$-boson pair. The error remains small throughout the whole invariant mass range, increasing somewhat towards large $M_{ZZ}$. The behavior for different values of $\tilde{x}$ is similar. To demonstrate the importance of including the threshold expansion we also show an approximation based solely on the LME as shaded regions. For this we adopt the prescription given in ref. Campbell et al. (2016) and show the envelope of the [2/2], [2/3], [3/2] and [3/3] Padé approximants which we have constructed without applying the rescaling of eq. (9) or the pole criterion eq. (13). We note that the resonant structure near $z=1$ in the upper right plot showing the vector form factor for maximal transverse momentum is caused by a pole near $w=1$ in the [3/3] Padé approximant. In our full results from eq. (10) and eq. (11) we apply the criterion eq. (13) to exclude approximants which feature such resonances in the time-like region $z\geq 0$. We conclude that the threshold expansion is essential for the reconstruction of the full top mass dependence above the top quark threshold.

We now turn to the NLO form factors. In fig. 3 we show the results for the virtual corrections to the form factors $\Ket{\widetilde{F}_{VV}^{(2)}}$ (upper panel) and $\Ket{\widetilde{F}_{AA}^{(2)}}$ (lower panel) for two values of $\tilde{x}$. Note that we do not include the double-triangle contribution to the form factors, as they have been computed analytically in Campbell et al. (2016). As at LO, we include only the top quark contributions. The uncertainty associated with the Padé construction increases with $M_{ZZ}$. Since we input information mainly at low $M_{ZZ}$ this behavior is expected. With the exception of the vector form factor $\Ket{\widetilde{F}_{VV}^{(2)}}$ for small transverse momenta (upper left panel in fig. 3) we find that the Padé approximation based on the LME alone does not yield a realistic reconstruction of the top-quark mass effects of the form factors. In particular, the important axial-vector form factor suffers from very large uncertainties. We remark though that in Campbell et al. (2016) for the NLO cross section the Padé prediction was improved by a reweighting with the full LO cross section.

We note that $\Ket{\widetilde{F}_{VV}^{(2)}}$ shows a small oscillation in the region of large $M_{ZZ}$ when the transverse momentum of the $Z$ bosons is small as is evident from the upper left plot in fig. 3. We trace the appearance of the second peak back to the contribution proportional to $L_{s}$ stemming from diagrams with massless cuts. In general, we find that this contribution shows worse convergence behavior than the non-logarithmic terms when including more and more terms in the LME and the threshold expansion. This is shown in fig. 4 where we compare our results from fig. 3 to the Padé approximants obtained with the same procedure but only using threshold input up to the order $\bar{z}^{2}$ and $\bar{z}^{3}$. We observe good convergence in the case of the axial-vector form factor. On the other hand, the $\mathcal{O}(\bar{z}^{2})$ approximation for the vector form factor does not feature the oscillatory behavior described above and there is no overlap with the full approximation in a significant part of the phase space. However, the $\mathcal{O}(\bar{z}^{3})$ and $\mathcal{O}(\bar{z}^{4})$ results are in good agreement with the full approximation where we have also included the $\mathcal{O}(\bar{z}^{5})$ term in the coefficient of the logarithm $L_{s}$ to verify that this stabilization persists with the addition of higher orders in the threshold expansion. We conclude from this discussion, that the Padé approximation can be improved systematically when including higher orders in the various expansions. Nevertheless, we believe that the prediction for $\Ket{\widetilde{F}_{VV}^{(2)}}$ should be taken with a grain of salt above $M_{ZZ}\geq 500\text{ GeV}$ because of the slower convergence.

In fig. 5 we show the virtual corrections to the form factor $v_{f}^{2}\Ket{\widetilde{F}_{VV}^{(2)}}+a_{f}^{2}\Ket{\widetilde{F}_{VV}^{(2)}}$ as it enters in the interference term with the Higgs boson exchange. The dashed lines show the form factor $v_{f}^{2}\Ket{\widetilde{F}_{VV}^{(2)}}$ increased by a factor of 300. This clearly demonstrates that the interference term will be dominated by $\Ket{\widetilde{F}_{AA}^{(2)}}$ and we therefore choose not to modify the uncertainty estimate for the vector form factor. The fact that $\Ket{\widetilde{F}_{VV}^{(2)}}$ is negligible compared to $\Ket{\widetilde{F}_{AA}^{(2)}}$ allows us to make trustworthy predictions for the interference with the Higgs production with subsequent decay to $Z$ bosons up to $M_{ZZ}\to\infty$, even though as stated above we trust our results for $\Ket{\widetilde{F}_{VV}^{(2)}}$ only for $M_{ZZ}\leq 500\text{ GeV}$.

The numerical implementation of the form factors is available as a FORTRAN routine on request and can be combined with existing computations of the massless loop contributions and the real corrections for the interference of the Higgs exchange with decay to $ZZ$ with the continuum background.

V Conclusions and outlook

We have considered top-quark mass effects in the continuum process $gg\to ZZ$, focusing on the form factors relevant for the NLO interference with the production of a Higgs boson and its subsequent decay into two $Z$ bosons. We have presented a Padé-based approximation using information from an expansion around a large top quark mass and an expansion around the top quark pair production threshold.

At LO, we have shown that our Padé construction approximates very well the full top mass dependence of the form factors for the whole range of the invariant mass $M_{ZZ}$ of the $Z$ bosons. At NLO, we provide a new prediction with very small uncertainties at small and moderate $M_{ZZ}$, with an increased uncertainty towards large $M_{ZZ}$. We expect that adding more information into the Padé construction at large $M_{ZZ}$ would improve the description also in this region.

Our results can be combined both with virtual corrections mediated by massless loops and the real corrections. The latter constitute a one-loop process and can therefore be computed with well-established techniques. The Padé construction can also be applied to the remaining form factors contributing to $gg\to ZZ$, which do not interfere with the Higgs signal.

We note also that while in this work we have applied our method to the production of on-shell $Z$ bosons, there is no obstruction for applying it also to off-shell $Z$ boson production. Indeed, the LME for off-shell $Z$ boson production is already known up to the order $z^{4}$ Caola et al. (2016). While a calculation of the full top mass dependence for on-shell $Z$ bosons with numerical methods seems to be feasible with current techniques in a reasonable time-frame (see Agarwal (27th May 2019); Agarwal and von Manteuffel ) a computation of the off-shell form factors appears to be beyond the current state-of-the-art.

Acknowledgements

We thank B. Agarwal, A. von Manteuffel and N. Kauer for their comments on the manuscript. This work has received funding from the European Union’s Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement No. 764850, SAGEX. RG is supported by the “Berliner Chancengleichheitsprogramm”.

Appendix A Threshold expansion of form factors

In the following we give explicit expressions for the coefficients in the threshold expansions of the form factors. For convenience, we quote the definition already given in eq. (5):

[TABLE]

where $i\in\{VV,AA\}$ and $\bar{n}_{2}$ is $n$ modulo $2$. The coefficients $a,b$ are most conveniently written in terms of the two dimensionless ratios $r_{Z}=\tfrac{m^{2}_{Z}}{M^{2}_{ZZ}}$ and $r_{p_{T}}=\tfrac{p_{T}^{2}}{M^{2}_{ZZ}}=\tilde{x}-r_{Z}$. We define the loop integral measure as

[TABLE]

and use the short-hand notation

[TABLE]

with $q^{2}=1,p_{Z}^{2}=4\*r_{Z}^{2},q\cdot p_{Z}=1$. The coefficients $b_{i,\ln}^{(n,1)}$ and $b_{i,\ln}^{(2n,m)}$ vanish. Furthermore, coefficients with $m=0$ and even $n$ do not contribute to the imaginary part and are therefore not listed here. We have calculated the remaining coefficients $a_{i}^{n,0},b_{i}^{n,m}$ up to $n=8$ and the coefficients $b_{i,\ln}^{(n,0)}$ up to $n=9$, obtaining the following results:

[TABLE]

for the expansion of the axial-vector component. The corresponding coefficients in the expansion of the vector part read

[TABLE]

Appendix B Subtractions

In this appendix, we give the functions $s_{i}$ with $i\in\{VV,AA\}$ used to subtract the threshold logarithms. We write them in terms of auxiliary subtraction functions $s_{n},n\in\mathbb{N}$, i.e.

[TABLE]

where the coefficients $C_{i,n}$ are constants and $s_{n}\mathop{\asymp}\limits^{z\to 1}\bar{z}^{\frac{n}{2}}\ln(\bar{z})+\mathcal{O}(\bar{z}^{\frac{n+1}{2}})$ in the threshold region. We construct these auxiliary functions based on the known analytical results for the vacuum polarization function. The subtraction functions and their threshold expansions are

[TABLE]

where we have used the symbol $\asymp$ to denote that terms analytical in $(1-z)$ have been dropped on the right-hand side, and we only use subtractions for the logarithmic terms, hence no subtraction functions $s_{1}(z)$ and $s_{3}(z)$ are necessary. We have used

[TABLE]

and $\Pi^{(1),v}$ is the well-known two-loop correction to the vacuum polarization Källén and Sabry (1955) in the convention of Kiyo et al. (2009). The functions $s_{i}$ in eqs. (54) – (59) are constant as $z\to 0$ and only diverge logarithmically as $z\to\infty$.