Electroweak Effects in e⁺e⁻→ZH Process

Abstract

Electroweak radiative corrections to the cross-section of the process $e^{+}{e}^{-}\to ZH$ are considered. The complete one-loop electroweak radiative corrections are evaluated with the help of the SANC system. Higher-order contributions of the initial-state radiation are computed in the QED structure function formalism. Numerical results are produced by a new version of the ReneSANCe event generator and MCSANCee integrator for the conditions of future electron-positron colliders. The resulting theoretical uncertainty in the description of this process is estimated.

1. Introduction

The Standard Model (SM) is extremely successful in describing various phenomena in particle physics. Despite this fact, there are many reasons to consider the SM as an effective model, i.e., a low-energy approximation of a more general theory. Looking for the limits of the SM applicability domain is one of the most valuable problems in modern fundamental physics. On the other hand, a deep investigation of the SM properties at the quantum level is still an important task since this model is relevant for many applications in high-energy physics as well as astrophysics and cosmology. In this context, exploring the Higgs boson sector of the SM is crucial for checking the mechanism of spontaneous symmetry breaking and finalizing the verification of the model within the energy range achieved at modern accelerators.

To perform an in-depth verification of the Standard Model and define the energy region of its applicability, we certainly need a new high-energy accelerator. An electron-positron collider with energies of a few hundred GeV looks now to be the best option. Several projects of this kind of machine are being under consideration, e.g., ILC [1], CLIC [2], CEPC [3], and FCC-ee [4]. Programs of all these colliders, except CLIC, necessarily include the option to run as Higgs factories with center-of-mass system (c.m.s.) energies of about 240 GeV. At this energy, the maximal count rate of events of the $e^{+}{e}^{-}\to ZH$ processes can be reached. Collecting several million such events will substantially increase the precision of the Higgs boson mass and the determination of the partial decay widths [5,6].

The expected high statistics of events with Higgs bosons challenges the theory to provide very accurate SM predictions for the corresponding processes with uncertainty at the permille level. So, we need to take into account radiative corrections in the first and higher orders of perturbation theory (Non-perturbative effects due to strong interactions are also relevant in running EW couplings and in producing extra meson pairs). The status of high-precision calculations for FCC-ee (and other future $e^{+}{e}^{-}$ colliders in general) is described in [7].

In this work, we analyze QED and electroweak (EW) radiative corrections to the higgsstrahlung process

$$e^{+}+e^{-}\to Z+H.$$

(1)

This process is the most promising one in studying the Higgs boson properties. So, the accuracy of its theoretical description should be higher than the experimental precision so that the combined uncertainty in the results of data analysis would not be spoiled by the theory. The uncertainty estimate should be as complete as possible.

In this paper, we evaluate the complete one-loop corrections supplemented by higher-order (HO) QED contributions in the leading logarithmic approximation (LLA) [8]. Our aim is to analyze the size of different HO contributions, estimate the resulting theoretical uncertainty, and verify the necessity to include other HO corrections. Please note that in this work we do not consider decays of Z and H, which are left for further study.

The complete one-loop electroweak radiative corrections to the process under consideration were computed with the help of the SANC computer system and reported in [9]. Here we will concentrate on the analysis of the HO QED effects. Recently, effects due to higher-order initial-state radiation (ISR) of photons in the process $e^{+}{e}^{-}\to {\mu }^{+}{\mu }^{-}b\overline{b}$ were considered in [10]. The channel with Z and H bosons in this process was included. It was claimed that the third order leading logarithmic contribution is numerically important and should be included. Here we study QED ISR corrections in more detail with taking into account radiation of light pairs, photonic leading logarithmic contributions up to the fourth order, and the complete one-loop (electro)weak effects.

Two-loop QED corrections due to the initial-state radiation for a general process of high-energy electron-positron annihilation through a virtual photon or Z boson were calculated in [11] and recently corrected in [12]. Higher-order QED ISR contributions in the leading and next-to-leading logarithmic approximations up to the order $\mathcal{O}\left({\alpha }^6{L}^5\right)$ were given in [13]. These results are performed in an inclusive set-up where only the distribution in the invariant mass of the final state particles is available. So, they provide a benchmark for comparisons while for practical applications one needs a Monte Carlo simulation with complete kinematics.

The paper is organized as follows. In Section 2, we outline the contributions due to the higher-order QED initial-state radiation order by order. In Section 3, we present the numerical results for the cross-section of associated $ZH$ production in the energy region $\sqrt{s}=200-500$ GeV. Our conclusions are given in Section 4.

2. ISR Corrections in LLA Approximation

2.1. General Notes

Let us consider ISR corrections to processes of high-energy electron-positron annihilation within the LLA. They can be evaluated with the help of the QED structure function formalism [8]. For ISR corrections in the annihilation channel the large logarithm is $L=ln(s/m_e^2)$ where the total c.m.s. energy $\sqrt{s}$ is chosen as factorization scale.

The master formula for a general $e^{+}{e}^{-}$ annihilation cross-section with ISR QED corrections in the leading logarithmic approximation has the same structure as the one for the Drell–Yan process, it reads

$$\begin{array}{c}\hfill {\sigma }^{\mathrm{LLA}}=\underset{0}{\overset{1}{\int }}\mathrm{d}{x}_1\underset{0}{\overset{1}{\int }}\mathrm{d}{x}_2{\mathcal{D}}_{ee}\left(x_1\right){\mathcal{D}}_{ee}\left(x_2\right){\sigma }_0(x_1,x_2,s)\Theta \left(\mathrm{cuts}\right),\end{array}$$

(2)

where ${\sigma }_0(x_1,x_2,s)$ is the Born-level cross-section of the annihilation process with reduced energies of the incoming particles. Here we do not take into account “photon induced” contributions, since the corresponding kernel cross-sections $\sigma (\gamma e\to eZH)$ and $\sigma (\gamma \gamma \to ZH)$ are very much suppressed by extra powers of the fine structure constant $\alpha$.

The electron structure functions, or better to say electron parton density functions, ${\mathcal{D}}_{ee}\left(x\right)$ describe the density of probability to find an electron with an energy fraction x in the initial electron beam [8,14,15]. Please note that the electron parton distribution functions in QED are completely analogous to proton PDF in QCD except the possibility to compute them in the QED case.

In the LLA approximation we can separate the pure photonic corrections (marked “$\gamma$”) and the remaining ones which include the pure pair and mixed photon-pair effects (marked as “pair”) as follows:

$${\mathcal{D}}_{ee}\left(x\right)={\mathcal{D}}_{ee}^{\gamma }\left(x\right)+{\mathcal{D}}_{ee}^{\mathrm{pair}}\left(x\right),$$

(3)

$$\begin{array}{cc}\hfill {\mathcal{D}}_{ee}^{\gamma }\left(x\right)& =\delta (1-x)+\frac{\alpha }{2\pi }(L-1)P^{\left(1\right)}\left(x\right)+{\left(\frac{\alpha }{2\pi }(L-1)\right)}^2\frac{1}{2!}{P}^{\left(2\right)}\left(x\right)\hfill \\ & +{\left(\frac{\alpha }{2\pi }(L-1)\right)}^3\frac{1}{3!}{P}^{\left(3\right)}\left(x\right)+{\left(\frac{\alpha }{2\pi }(L-1)\right)}^4\frac{1}{4!}{P}^{\left(4\right)}\left(x\right)+\mathcal{O}\left({\alpha }^5{L}^5\right),\hfill \end{array}$$

(4)

$$\begin{array}{cc}\hfill {\mathcal{D}}_{ee}^{\mathrm{pair}}\left(x\right)& ={\left(\frac{\alpha }{2\pi }L\right)}^2\left[\frac{1}{3}{P}^{\left(1\right)}\left(x\right)+\frac{1}{2}{R}_s\left(x\right)\right]\hfill \\ & +{\left(\frac{\alpha }{2\pi }L\right)}^3\left[\frac{1}{3}{P}^{\left(2\right)}\left(x\right)+\frac{4}{27}{P}^{\left(1\right)}\left(x\right)+\frac{1}{3}{R}_p\left(x\right)-\frac{1}{9}{R}_s\left(x\right)\right]+\mathcal{O}\left({\alpha }^4{L}^4\right).\hfill \end{array}$$

(5)

Pair corrections can be separated into singlet $(\sim$$R_{s,p})$ and non-singlet ones $(\sim$$P^{\left(n\right)})$. We take into account both by default.

Non-singlet splitting functions can be represented in the form

$$\begin{array}{c}\hfill P^{\left(n\right)}\left(x\right)=\underset{\Delta \to 0}{lim}\left\{\delta (1-x)P_{\Delta }^{\left(n\right)}+P_{\Theta }^{\left(n\right)}\left(x\right)\Theta (1-\Delta -x)\right\}\end{array}$$

(6)

with $\Delta \ll 1$ being the soft-hard separator. For example,

$$\begin{array}{c}\hfill P_{\Delta }^{\left(1\right)}=2ln\Delta +\frac{3}{2},P_{\Theta }^{\left(1\right)}\left(x\right)=\frac{1+x^2}{1-x}.\end{array}$$

(7)

Higher-order non-singlet pure photonic splitting functions are obtained by iterations of convolution

$$\begin{array}{c}\hfill P^{(n+1)}\left(x\right)=\underset{0}{\overset{1}{\int }}\mathrm{d}{x}_1\underset{0}{\overset{1}{\int }}\mathrm{d}{x}_2\delta (x-x_1{x}_2)P^{\left(n\right)}\left(x_1\right)P^{\left(1\right)}\left(x_2\right),\end{array}$$

(8)

see the details in Ref. [15].

The singlet splitting functions $R_s$ and $R_p$ are not singular at $x\to 1$, so they do not contain $\Delta$ parts. Explicit expressions for all relevant splitting functions are given in Refs. [15,16].

The Born-level partonic cross-section ${\sigma }_0(x_1,x_2,s)$ is known in the partonic c.m.s. as ${\sigma }_{\mathrm{Born}}\left(\widehat{s}\right)$, where $\widehat{s}=x_1{x}_{2}s$. The transition from the partonic c.m.s. into the laboratory reference frame is required if $x_1{x}_2\ne 1$.

Let us classify contributions with different kinematics:

I.

${\left(SV\right)}_1\times {\left(SV\right)}_2$ The Born kinematics: additional contributions to Born+Soft+Virt.

II.

$H_1\times {\left(SV\right)}_2$ One hard photon collinear to the first initial particle with possible soft and/or virtual (Soft+Virt) corrections to the second one. Hereafter “One hard photon” means “at least one hard photon in the same direction”.

III.

${\left(SV\right)}_1\times H_2$ Soft+Virt to the first initial particle and one hard photon along the second initial particle.

IV.

$H_1\times H_2$ One hard photon along the first initial particle and one along the second one.

Separation of hard and soft photon emission is provided by the dimensionless parameter $\Delta \ll 1$ with typical values ${10}^{-3}$, ${10}^{-4}$. Under all integrals relevant (process dependent) cuts on the lower values of $x_1$ and $x_2$ variables should be applied.

Application of representation (6) in structure functions (3) and their substitution into the master Equation (2) gives the corrected cross-section in the LLA approximation. We expand the result in $\alpha$ and look at the second, third, and fourth order contributions.

A few general comments are in order:

•

The first lower index below denotes the order in $\alpha L$.

•

Factorials and coefficients are given explicitly in order to see their structure.

•

For pure photonic LLA corrections the traditional shift $L⟶(L-1)$ is carried out, it takes into account part of the next-to-leading (NLO) corrections. However, for pair corrections such a shift is not well justified, and we keep the large log unchanged.

2.2. First Order LLA Contributions

There are only photonic corrections in $\mathcal{O}\left(\alpha \right)$. Below we list the contributions of different kinematics.

I. Born kinematics

$$\begin{array}{ccc}& & \delta {\sigma }_{1,\gamma }^{\left(I\right)}=\left(\frac{\alpha }{2\pi }(L-1)\right){\sigma }_0(1,1,s)\left\{2P_{\Delta }^{\left(1\right)}\right\}.\hfill \end{array}$$

(9)

II. Emission only along the first particle

$$\begin{array}{ccc}& & \delta {\sigma }_{1,\gamma }^{\left(II\right)}=\left(\frac{\alpha }{2\pi }(L-1)\right)\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_1{\sigma }_0(x_1,1,s)\left\{P_{\Theta }^{\left(1\right)}\left(x_1\right)\right\}.\hfill \end{array}$$

(10)

III. Emission only along the second particle

$$\begin{array}{ccc}& & \delta {\sigma }_{1,\gamma }^{\left(III\right)}=\left(\frac{\alpha }{2\pi }(L-1)\right)\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_2{\sigma }_0(1,x_2,s)\left\{P_{\Theta }^{\left(1\right)}\left(x_2\right)\right\}.\hfill \end{array}$$

(11)

IV. Emission along both initial particles

$$\begin{array}{c}\hfill \delta {\sigma }_{1,\gamma }^{\left(IV\right)}=0.\end{array}$$

(12)

2.3. Second Order LLA Contributions

I. Born kinematics

$$\begin{array}{ccc}& & \delta {\sigma }_{2,\gamma }^{\left(I\right)}={\left(\frac{\alpha }{2\pi }(L-1)\right)}^2{\sigma }_0(1,1,s)\left\{2\frac{1}{2!}{P}_{\Delta }^{\left(2\right)}+P_{\Delta }^{\left(1\right)}{P}_{\Delta }^{\left(1\right)}\right\},\hfill \end{array}$$

(13)

$$\begin{array}{ccc}& & \delta {\sigma }_{2,\mathrm{pair}}^{\left(I\right)}={\left(\frac{\alpha }{2\pi }L\right)}^2{\sigma }_0(1,1,s)\left\{2\frac{1}{3}{P}_{\Delta }^{\left(1\right)}\right\}.\hfill \end{array}$$

(14)

II. Emission only along the first particle

$$\begin{array}{ccc}& & \delta {\sigma }_{2,\gamma }^{\left(II\right)}={\left(\frac{\alpha }{2\pi }(L-1)\right)}^2\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_1{\sigma }_0(x_1,1,s)\left\{\frac{1}{2!}{P}_{\Theta }^{\left(2\right)}\left(x_1\right)+P_{\Theta }^{\left(1\right)}\left(x_1\right)P_{\Delta }^{\left(1\right)}\right\},\hfill \end{array}$$

(15)

$$\begin{array}{ccc}& & \delta {\sigma }_{2,\mathrm{pair}}^{\left(II\right)}={\left(\frac{\alpha }{2\pi }L\right)}^2\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_1{\sigma }_0(x_1,1,s)\left\{\frac{1}{3}{P}_{\Theta }^{\left(1\right)}\left(x_1\right)+\frac{1}{2}{R}_s\left(x_1\right)\right\}.\hfill \end{array}$$

(16)

III. Emission only along the second particle

$$\begin{array}{ccc}& & \delta {\sigma }_{2,\gamma }^{\left(III\right)}={\left(\frac{\alpha }{2\pi }(L-1)\right)}^2\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_2{\sigma }_0(1,x_2,s)\left\{\frac{1}{2!}{P}_{\Theta }^{\left(2\right)}\left(x_2\right)+P_{\Theta }^{\left(1\right)}\left(x_2\right)P_{\Delta }^{\left(1\right)}\right\},\hfill \end{array}$$

(17)

$$\begin{array}{ccc}& & \delta {\sigma }_{2,\mathrm{pair}}^{\left(III\right)}={\left(\frac{\alpha }{2\pi }L\right)}^2\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_2{\sigma }_0(1,x_2,s)\left\{\frac{1}{3}{P}_{\Theta }^{\left(1\right)}\left(x_2\right)+\frac{1}{2}{R}_s\left(x_2\right)\right\}.\hfill \end{array}$$

(18)

IV. Emission along both initial particles

$$\begin{array}{ccc}& & \delta {\sigma }_{2,\gamma }^{\left(IV\right)}={\left(\frac{\alpha }{2\pi }(L-1)\right)}^2\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_1\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_2{\sigma }_0(x_1,x_2,s)\left\{P_{\Theta }^{\left(1\right)}\left(x_1\right)P_{\Theta }^{\left(1\right)}\left(x_2\right)\right\},\hfill \end{array}$$

(19)

$$\begin{array}{ccc}& & \delta {\sigma }_{2,\mathrm{pair}}^{\left(IV\right)}=0.\hfill \end{array}$$

(20)

2.4. Third Order LLA Contributions

I. Born kinematics

$$\begin{array}{ccc}& & \delta {\sigma }_{3,\gamma }^{\left(I\right)}={\left(\frac{\alpha }{2\pi }(L-1)\right)}^3{\sigma }_0(1,1,s)\left\{2\frac{1}{3!}{P}_{\Delta }^{\left(3\right)}+2P_{\Delta }^{\left(1\right)}\frac{1}{2!}{P}_{\Delta }^{\left(2\right)}\right\},\hfill \end{array}$$

(21)

$$\begin{array}{ccc}& & \delta {\sigma }_{3,\mathrm{pair}}^{\left(I\right)}={\left(\frac{\alpha }{2\pi }L\right)}^3{\sigma }_0(1,1,s)\left\{2\frac{1}{3}{P}_{\Delta }^{\left(2\right)}+2\frac{4}{27}{P}_{\Delta }^{\left(1\right)}+2\frac{1}{3}{P}_{\Delta }^{\left(1\right)}{P}_{\Delta }^{\left(1\right)}\right\}.\hfill \end{array}$$

(22)

II. Emission only along the first particle

$$\delta {\sigma }_{3,\gamma }^{\left(II\right)}={\left(\frac{\alpha }{2\pi }(L-1)\right)}^3\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_1{\sigma }_0(x_1,1,s)\{\frac{1}{3!}{P}_{\Theta }^{\left(3\right)}\left(x_1\right)+\frac{1}{2!}{P}_{\Theta }^{\left(2\right)}\left(x_1\right)P_{\Delta }^{\left(1\right)} + P_{\Theta }^{\left(1\right)}\left(x_1\right)\frac{1}{2!}{P}_{\Delta }^{\left(2\right)}\},$$

(23)

$$\begin{array}{cc}\hfill \delta {\sigma }_{3,\mathrm{pair}}^{\left(II\right)}& ={\left(\frac{\alpha }{2\pi }L\right)}^3\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_1{\sigma }_0(x_1,1,s)\{\frac{1}{3}{P}_{\Theta }^{\left(2\right)}\left(x_1\right)+\frac{4}{27}{P}_{\Theta }^{\left(1\right)}\left(x_1\right)+\frac{1}{3}{R}_p\left(x_1\right)\hfill \\ & -\frac{1}{9}{R}_s\left(x_1\right)+2\frac{1}{3}{P}_{\Theta }^{\left(1\right)}\left(x_1\right)P_{\Delta }^{\left(1\right)}+\frac{1}{2}{R}_s\left(x_1\right)P_{\Delta }^{\left(1\right)}\}.\hfill \end{array}$$

(24)

III. Emission only along the second particle

$$\delta {\sigma }_{3,\gamma }^{\left(II\right)}={\left(\frac{\alpha }{2\pi }(L-1)\right)}^3\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_2{\sigma }_0(1,x_2,s)\{\frac{1}{3!}{P}_{\Theta }^{\left(3\right)}\left(x_2\right)+\frac{1}{2!}{P}_{\Theta }^{\left(2\right)}\left(x_2\right)P_{\Delta }^{\left(1\right)}+P_{\Theta }^{\left(1\right)}\left(x_2\right)\frac{1}{2!}{P}_{\Delta }^{\left(2\right)}\},$$

(25)

$$\begin{array}{cc}\hfill \delta {\sigma }_{3,\mathrm{pair}}^{\left(II\right)}& ={\left(\frac{\alpha }{2\pi }L\right)}^3\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_2{\sigma }_0(1,x_2,s)\{\frac{1}{3}{P}_{\Theta }^{\left(2\right)}\left(x_2\right)+\frac{4}{27}{P}_{\Theta }^{\left(1\right)}\left(x_2\right)+\frac{1}{3}{R}_p\left(x_2\right)\hfill \\ & -\frac{1}{9}{R}_s\left(x_2\right)+2\frac{1}{3}{P}_{\Theta }^{\left(1\right)}\left(x_2\right)P_{\Delta }^{\left(1\right)}+\frac{1}{2}{R}_s\left(x_2\right)P_{\Delta }^{\left(1\right)}\}.\hfill \end{array}$$

(26)

IV. Emission along both initial particles

$$\delta {\sigma }_{3,\gamma }^{\left(IV\right)}={\left(\frac{\alpha }{2\pi }(L-1)\right)}^3\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_1\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_2{\sigma }_0(x_1,x_2,s)\{\frac{1}{2!}{P}_{\Theta }^{\left(2\right)}\left(x_1\right)P_{\Theta }^{\left(1\right)}\left(x_2\right)+P_{\Theta }^{\left(1\right)}\left(x_1\right)\frac{1}{2!}{P}_{\Theta }^{\left(2\right)}\left(x_2\right)\},$$

(27)

$$\begin{array}{cc}\hfill \delta {\sigma }_{3,\mathrm{pair}}^{\left(IV\right)}& ={\left(\frac{\alpha }{2\pi }L\right)}^3\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_1\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_2{\sigma }_0(x_1,x_2,s)\{\frac{1}{3}{P}_{\Theta }^{\left(1\right)}\left(x_1\right)P_{\Theta }^{\left(1\right)}\left(x_2\right)+\frac{1}{2}{R}_s\left(x_1\right)P_{\Theta }^{\left(1\right)}\left(x_2\right)\hfill \\ & +P_{\Theta }^{\left(1\right)}\left(x_1\right)\frac{1}{3}{P}_{\Theta }^{\left(1\right)}\left(x_2\right)+P_{\Theta }^{\left(1\right)}\left(x_1\right)\frac{1}{2}{R}_s\left(x_2\right)\}.\hfill \end{array}$$

(28)

2.5. Fourth Order LLA Contributions

Here we list only pure photonic contributions due to the smallness of pair corrections in the fourth order.

I. Born kinematics

$$\delta {\sigma }_{4,\gamma }^{\left(I\right)}={\left(\frac{\alpha }{2\pi }(L-1)\right)}^4{\sigma }_0(1,1,s)\left\{2\frac{1}{4!}{P}_{\Delta }^{\left(4\right)}+2P_{\Delta }^{\left(1\right)}\frac{1}{3!}{P}_{\Delta }^{\left(3\right)}+{\left(\frac{1}{2!}{P}_{\Delta }^{\left(2\right)}\right)}^2\right\}.$$

(29)

II. Emission only along the first particle

$$\begin{array}{cc}\hfill \delta {\sigma }_{4,\gamma }^{\left(II\right)}& ={\left(\frac{\alpha }{2\pi }(L-1)\right)}^4\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_1{\sigma }_0(x_1,1,s)\{\frac{1}{4!}{P}_{\Theta }^{\left(4\right)}\left(x_1\right)+\frac{1}{3!}{P}_{\Theta }^{\left(3\right)}\left(x_1\right)P_{\Delta }^{\left(1\right)}\hfill \\ & +\frac{1}{2!}{P}_{\Theta }^{\left(2\right)}\left(x_1\right)\frac{1}{2!}{P}_{\Delta }^{\left(2\right)}+P_{\Theta }^{\left(1\right)}\left(x_1\right)\frac{1}{3!}{P}_{\Delta }^{\left(3\right)}\}.\hfill \end{array}$$

(30)

III. Emission only along the second particle

$$\begin{array}{cc}\hfill \delta {\sigma }_{4,\gamma }^{\left(II\right)}& ={\left(\frac{\alpha }{2\pi }(L-1)\right)}^4\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_2{\sigma }_0(1,x_2,s)\{\frac{1}{4!}{P}_{\Theta }^{\left(4\right)}\left(x_2\right)+\frac{1}{3!}{P}_{\Theta }^{\left(3\right)}\left(x_2\right)P_{\Delta }^{\left(1\right)}\hfill \\ & +\frac{1}{2!}{P}_{\Theta }^{\left(2\right)}\left(x_2\right)\frac{1}{2!}{P}_{\Delta }^{\left(2\right)}+P_{\Theta }^{\left(1\right)}\left(x_2\right)\frac{1}{3!}{P}_{\Delta }^{\left(3\right)}\}.\hfill \end{array}$$

(31)

IV. Emission along both initial particles

$$\begin{array}{cc}\hfill \delta {\sigma }_{4,\gamma }^{\left(IV\right)}& ={\left(\frac{\alpha }{2\pi }(L-1)\right)}^4\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_1\underset{0}{\overset{1-\Delta }{\int }}\mathrm{d}{x}_2{\sigma }_0(x_1,x_2,s)\{\frac{1}{2!}{P}_{\Theta }^{\left(2\right)}\left(x_1\right)\frac{1}{2!}{P}_{\Theta }^{\left(2\right)}\left(x_2\right)\hfill \\ & +P_{\Theta }^{\left(1\right)}\left(x_1\right)\frac{1}{3!}{P}_{\Theta }^{\left(3\right)}\left(x_2\right)+\frac{1}{3!}{P}_{\Theta }^{\left(3\right)}\left(x_1\right)P_{\Theta }^{\left(1\right)}\left(x_2\right)\}.\hfill \end{array}$$

(32)

2.6. LLA Contributions for Helicity States

The leading-order (LO) splitting function $P_{ee}\left(x\right)$ given by Equation (7) preserves helicity [17], i.e.,

$$\begin{array}{c}\hfill P_{e_{+}{e}_{+}}^{\left(n\right)}\left(x\right)=P_{e_{-}{e}_{-}}^{\left(n\right)}\left(x\right)=P_{ee}^{\left(n\right)}\left(x\right),P_{e_{+}{e}_{-}}^{\left(n\right)}\left(x\right)=P_{e_{-}{e}_{+}}^{\left(n\right)}\left(x\right)=0.\end{array}$$

(33)

However, singlet contributions of pair corrections can be separated for different helicities:

$$\begin{array}{ccc}& & R_s\left(x\right)=R_{s,e_{+}{e}_{+}}\left(x\right)+R_{s,e_{-}{e}_{+}}\left(x\right),R_p\left(x\right)=R_{p,e_{+}{e}_{+}}\left(x\right)+R_{p,e_{-}{e}_{+}}\left(x\right),\hfill \end{array}$$

(34)

$$\begin{array}{ccc}& & R_{s,e_{+}{e}_{+}}\left(x\right)=3(1-x)+2(1+x)lnx+\frac{2(1-x^3)}{3x},\hfill \end{array}$$

(35)

$$R_{s,e_{-}{e}_{+}}\left(x\right)=-2(1-x)+\frac{2(1-x^3)}{3x},$$

(36)

$$\begin{array}{cc}\hfill R_{p,e_{+}{e}_{+}}\left(x\right)& =-\frac{31}{6}(1-x)+6(1-x)ln(1-x)-2lnx+4xlnx+\frac{4(1-x^3)}{3x}ln(1-x)\hfill \\ & +\frac{4}{3}{x}^{2}lnx+4(1+x)lnxln(1-x)-(1+x){ln}^{2}x+4(1+x){\mathrm{Li}}_2\left(1-x\right),\hfill \end{array}$$

(37)

$$\begin{array}{cc}\hfill R_{p,e_{-}{e}_{+}}\left(x\right)& =\frac{4}{3}(1-x)+\frac{4(1-x^3)}{3x}ln(1-x)-4(1-x)ln(1-x)\hfill \\ & +2(1-x)lnx+\frac{4}{3}{x}^{2}lnx.\hfill \end{array}$$

(38)

Since QED preserves parity,

$$\begin{array}{c}\hfill R_{i,e_{-}{e}_{+}}\left(x\right)=R_{i,e_{+}{e}_{-}}\left(x\right),R_{i,e_{-}{e}_{-}}\left(x\right)=R_{i,e_{+}{e}_{+}}\left(x\right).\end{array}$$

(39)

2.7. Scheme with Exponentiation

In the master Equation (2), the electron structure functions can be taken in the exponentiated form [18,19]. That would mean continuous integration over $x_{1,2}$ without the auxiliary parameter $\Delta$. The result of [19] contains only the pure photonic corrections and corresponds to the inclusion of exact leading logs up to the order $\mathcal{O}\left({\alpha }^5{L}^5\right)$ together with approximate (incomplete) HO LLA effects. Please note that the HO exponentiated effects become exact (in LLA) in the soft photon limit.

•

The pair LLA corrections can be added to the result of [19] as, e.g., in Ref. [20] with a possible update for higher-order pair corrections listed above.

•

The exponentiated structure functions include the LLA part of one-loop QED radiative corrections. To avoid double counting with complete one-loop corrections, we need to subtract the first order leading logarithmic terms from one-loop corrections. So, the final result reads

$$\begin{array}{cc}\hfill {\sigma }^{\mathrm{corr}.}& =\underset{0}{\overset{1}{\int }}\mathrm{d}{x}_1\underset{0}{\overset{1}{\int }}\mathrm{d}{x}_2{\mathcal{D}}_{ee}^{\exp }\left(x_1\right){\mathcal{D}}_{ee}^{\exp }\left(x_2\right){\sigma }_0(x_1,x_2,s)\Theta \left(\mathrm{cuts}\right)\hfill \\ & +\left[{\sigma }^{\mathrm{Soft}+\mathrm{Virt}}-{\sigma }_{\mathrm{LLA}}^{\mathrm{Soft}+\mathrm{Virt}}\right]+\left[{\sigma }^{\mathrm{Hard}}-{\sigma }_{\mathrm{LLA}}^{\mathrm{Hard}}\right].\hfill \end{array}$$

(40)

The “Soft+Virt” part has the Born-like kinematics:

$$\begin{array}{c}\hfill {\sigma }_{\mathrm{LLA}}^{\mathrm{Soft}+\mathrm{Virt}}={\sigma }^{\mathrm{Born}}2\frac{\alpha }{2\pi }(L-1)\left[2ln\omega +\frac{3}{2}\right],\end{array}$$

(41)

where $\omega$ is the dimensionless soft-hard separator from the original complete one-loop formulae implemented in the Monte Carlo generator.

The “Hard LLA” term is rewritten to match “Hard” kinematics:

$$\begin{array}{c}\hfill {\sigma }_{\mathrm{LLA}}^{\mathrm{Hard}}=\frac{\alpha }{2\pi }\int \frac{{\mathrm{d}}^{3}k}{k_0^{2}2\pi }\left[\frac{E^2}{k{p}_1}{\sigma }^{\mathrm{Born}}(x_1,1,s)+\frac{E^2}{k{p}_2}{\sigma }^{\mathrm{Born}}(1,x_2,s)\right]\left(2-2\frac{k_0}{E}+\frac{k_0^2}{E^2}\right).\end{array}$$

(42)

3. Numerical Results

In this section, we show numerical results for one-loop EW and HO QED radiative corrections to the $e^{+}{e}^{-}\to ZH$ process. The input parameters can be found in [9]. The results are obtained without any angular cuts. The relative correction $\delta$ (in %) is defined as

$$\delta =\frac{{\sigma }_{\mathtt{LLA}}}{{\sigma }^{\mathrm{Born}}}-1.$$

(43)

To illustrate the trends of the ISR contribution behavior, we present separate distributions for each $\mathcal{O}\left({\alpha }^n{L}^n\right)$, $n=2-4$ term and their sum as a function of the c.m.s. energy.

Figure 1 shows the values of the dominant contribution $\mathcal{O}\left({\alpha }^2{L}^2\right)\gamma$ and the sum of all considered orders of the ISR terms ${\sum }_{n=2}^4\mathcal{O}\left({\alpha }^n{L}^n\right)$ vs. c.m.s. energy. The dominant contribution $\mathcal{O}\left({\alpha }^2{L}^2\right)\gamma$ is about 3% at the c.m.s. energy $\sqrt{s}=220$ GeV, then it crosses the zero line approximately at $\sqrt{s}=250$ GeV, reaches the minimum value about $-0.9\%$ at 325 GeV and goes to $-0.5\%$ at $\sqrt{s}=500$ GeV. The sum is mainly determined by the $\mathcal{O}\left({\alpha }^2{L}^2\right)\gamma$ term in the region of c.m.s energies $\sqrt{s}=220-290$ GeV and becomes close to zero at $\sqrt{s}=500$ GeV.

Figure 2 demonstrates the values of the contribution of relative correction (in %): (a) $\mathcal{O}\left({\alpha }^2{L}^2\right)e^{+}{e}^{-}$, (b) $\mathcal{O}\left({\alpha }^2{L}^2\right){\mu }^{+}{\mu }^{-}$, (c) $\mathcal{O}\left({\alpha }^3{L}^3\right)\gamma$ vs. c.m.s. energy. Figure 3 shows the values of the contribution of relative correction (in %): (a) $\mathcal{O}\left({\alpha }^3{L}^3\right)e^{+}{e}^{-}$, (b) $\mathcal{O}\left({\alpha }^3{L}^3\right){\mu }^{+}{\mu }^{-}$, (c) $\mathcal{O}\left({\alpha }^4{L}^4\right)\gamma$ vs. c.m.s. energy. The second order contributions due to light pair emission are smaller than second order photonic corrections. The suppression of pair corrections with respect to photonic ones in the same order is typical in annihilation processes at LEP energies [21]. However, there is a kinematical region where they are comparable. See, for example regions near 240 GeV (

Table 1. ISR corrections in the LLA approximation for the $e^{+}{e}^{-}\to ZH$ process at $\sqrt{s}=240$ GeV. No cuts are imposed. Here ${\delta }_{\mathrm{ISR}\mathrm{LLA}}\equiv \delta {\sigma }_{\mathrm{ISR}\mathrm{LLA}}/{\sigma }_0$. The Born cross-section is ${\sigma }_0=225.74\left(1\right)$ fb.

	$\mathcal{O}\left({\mathit{\alpha }}^2{\mathit{L}}^2\right)$	$\mathcal{O}\left({\mathit{\alpha }}^2{\mathit{L}}^2\right)$	$\mathcal{O}\left({\mathit{\alpha }}^2{\mathit{L}}^2\right)$	$\mathcal{O}\left({\mathit{\alpha }}^3{\mathit{L}}^3\right)$	$\mathcal{O}\left({\mathit{\alpha }}^3{\mathit{L}}^3\right)$	$\mathcal{O}\left({\mathit{\alpha }}^3{\mathit{L}}^3\right)$	$\mathcal{O}\left({\mathit{\alpha }}^4{\mathit{L}}^4\right)$	${\displaystyle {\sum }_{\mathit{n}=2}^4}\mathcal{O}\left({\mathit{\alpha }}^{\mathit{n}}{\mathit{L}}^{\mathit{n}}\right)$
	$\mathbf{\gamma }$	${\mathbf{e}}^{+}{\mathbf{e}}^{-}$	${\mathbf{\mu }}^{+}{\mathbf{\mu }}^{-}$	$\mathbf{\gamma }$	${\mathbf{e}}^{+}{\mathbf{e}}^{-}$	${\mathbf{\mu }}^{+}{\mathbf{\mu }}^{-}$	$\mathbf{\gamma }$
$\delta {\sigma }_{\mathrm{LLA}}$, fb	$1.128\left(1\right)$	$-0.368\left(1\right)$	$-0.218\left(1\right)$	$0.176\left(1\right)$	$0.019\left(1\right)$	$0.011\left(1\right)$	$-0.023\left(1\right)$	$0.727\left(1\right)$
${\delta }_{\mathrm{LLA}}$, %	$0.500\left(1\right)$	$-0.163\left(1\right)$	$-0.097\left(1\right)$	$0.078\left(1\right)$	$0.008\left(1\right)$	$0.005\left(1\right)$	$-0.010\left(1\right)$	$0.322\left(1\right)$

) and 250 GeV (

Table 2. ISR corrections in the LLA approximation for the $e^{+}{e}^{-}\to ZH$ process at $\sqrt{s}=250$ GeV. No cuts are imposed. Here ${\delta }_{\mathrm{ISR}\mathrm{LLA}}\equiv \delta {\sigma }_{\mathrm{ISR}\mathrm{LLA}}/{\sigma }_0$. The Born cross-section is ${\sigma }_0=225.59\left(1\right)$ fb.

	$\mathcal{O}\left({\mathit{\alpha }}^2{\mathit{L}}^2\right)$	$\mathcal{O}\left({\mathit{\alpha }}^2{\mathit{L}}^2\right)$	$\mathcal{O}\left({\mathit{\alpha }}^2{\mathit{L}}^2\right)$	$\mathcal{O}\left({\mathit{\alpha }}^3{\mathit{L}}^3\right)$	$\mathcal{O}\left({\mathit{\alpha }}^3{\mathit{L}}^3\right)$	$\mathcal{O}\left({\mathit{\alpha }}^3{\mathit{L}}^3\right)$	$\mathcal{O}\left({\mathit{\alpha }}^4{\mathit{L}}^4\right)$	${\displaystyle {\sum }_{\mathit{n}=2}^4}\mathcal{O}\left({\mathit{\alpha }}^{\mathit{n}}{\mathit{L}}^{\mathit{n}}\right)$
	$\mathbf{\gamma }$	${\mathbf{e}}^{+}{\mathbf{e}}^{-}$	${\mathbf{\mu }}^{+}{\mathbf{\mu }}^{-}$	$\mathbf{\gamma }$	${\mathbf{e}}^{+}{\mathbf{e}}^{-}$	${\mathbf{\mu }}^{+}{\mathbf{\mu }}^{-}$	$\mathbf{\gamma }$
$\delta {\sigma }_{\mathrm{LLA}}$, fb	$-0.223\left(1\right)$	$-0.268\left(1\right)$	$-0.159\left(1\right)$	$0.211\left(1\right)$	$-0.010\left(1\right)$	$-0.006\left(1\right)$	$-0.016\left(1\right)$	$-0.468\left(1\right)$
${\delta }_{\mathrm{LLA}}$, %	$-0.099\left(1\right)$	$-0.119\left(1\right)$	$-0.070\left(1\right)$	$0.094\left(1\right)$	$-0.004\left(1\right)$	$-0.003\left(1\right)$	$-0.007\left(1\right)$	$-0.207\left(1\right)$

) and also near 500 GeV. The third (fourth) order photonic corrections are approximately 10 (100) times smaller than $\mathcal{O}\left({\alpha }^2{L}^2\right)\gamma$.

One can see that in the threshold energy region there are several competing contributions with different behavior. This confirms the necessity to take into account HO QED ISR contributions in the studies of the higgsstrahlung process at future colliders. The $\mathcal{O}\left({\alpha }^3{L}^3\right)$ pair and $\mathcal{O}\left({\alpha }^4{L}^4\right)$ photon contributions should be taken into account for the ${10}^{-4}$ accuracy goal (see Figure 3).

Figure 4 illustrates the behavior of the cross-sections with respect to the c.m.s. energy. It is seen that at the peak in the threshold region, the one-loop QED corrections change not only the height of the peak but also its shape and position.

Figure 5 complements Figure 4. It shows the size of the relative RC in different approximations. One can see the difference between the first order $\mathcal{O}\left(\alpha L\right)$ term (a) and exact $\mathcal{O}\left(\alpha \right)$ QED corrections (b) and that the HO ISR LLA contributions provide a small but visible shift (the difference between lines (c) and (d)) from the complete one-loop EW correction. Moreover, this shift changes its sign.

In Table 1 and Table 2 we show the ISR corrections of different order of $\mathcal{O}\left({\alpha }^n{L}^n\right),n$ = 2–4 in the LLA approximation for the c.m.s. energies $\sqrt{s}=240$ GeV and 250 GeV in the $\alpha \left(0\right)$ EW scheme. We provide numbers in two points because these energies are particular for the process under consideration. Specifically, the Born-level cross-section has a peak at ∼240 GeV while the present plans of future $e^{+}{e}^{-}$ colliders envisage operation at 250 GeV where the counting rate of the signal is higher [3,4].

It is seen in Table 1 and Table 2 that the corrections for the sum of all considered orders of the ISR terms ${\sum }_{n=2}^4\mathcal{O}\left({\alpha }^n{L}^n\right)$ are about 0.322% for the c.m.s. energy $\sqrt{s}=240$ GeV and about −0.207% for the c.m.s. energy $\sqrt{s}=250$ GeV. For the c.m.s. energy $\sqrt{s}=240$ GeV the most significant HO contribution is of course the photonic one of the order $\mathcal{O}{\left(\alpha L\right)}^2$. It composes about half a percent while from pairs we obtain about −(0.1–0.2)%. For the c.m.s. energy $\sqrt{s}=250$ GeV the dominant contributions of the second order are about $-0.099\%$ for $\gamma$ and $-0.119\%$ for $e^{+}{e}^{-}$-pairs ($-0.070\%$ for ${\mu }^{+}{\mu }^{-}$-pairs), respectively. When considering HO corrections, we see that it is certainly sufficient to take into account corrections up to the fourth order.

Variation of the factorization scale in the argument of the large logarithm can simulate the next-to-leading logarithmic corrections, e.g., $\mathcal{O}\left({\alpha }^{2}L\right)$. In the same manner as in estimates of scale variation uncertainties in QCD, we apply factors 2 and 1/2 to the energy scale in the argument of the large logarithm. This leads to the following values of the HO LLA corrections at $\sqrt{s}=240$ GeV: ${\delta }_{\mathrm{LLA}}\left(2\sqrt{s}\right)=0.361\left(1\right)\%$ and ${\delta }_{\mathrm{LLA}}(\sqrt{s}/2)=0.286\left(1\right)\%$, respectively. And for $\sqrt{s}=250$ GeV we obtain ${\delta }_{\mathrm{LLA}}\left(2\sqrt{s}\right)=-0.228\left(1\right)\%$ and ${\delta }_{\mathrm{LLA}}(\sqrt{s}/2)=-0.187\left(1\right)\%$.

In

Table 3. The Born and pure weak corrections in different EW schemes at the c.m.s. energy $\sqrt{s}=240$ GeV.

EW Scheme	$\mathit{\alpha }\left(0\right)$	${\mathit{G}}_{\mathit{\mu }}$	$\mathit{\alpha }\left({\mathit{M}}_{\mathit{Z}}^2\right)$
${\sigma }_{\mathrm{Born}}$, fb	$225.74\left(1\right)$	$240.43\left(1\right)$	$254.65$(1)
${\sigma }_{\mathrm{Born}+\mathrm{PW}}$, fb	$231.88\left(1\right)$	$233.25\left(1\right)$	$231.80\left(1\right)$
$\delta {\sigma }_{\mathrm{PW}}$, fb	$6.15\left(1\right)$	$-7.18\left(1\right)$	$-22.85\left(1\right)$
${\delta }_{\mathrm{PW}}$, %	$2.72\left(1\right)$	$-2.99\left(1\right)$	$-8.97\left(1\right)$

and

Table 4. The Born and pure weak corrections in different EW schemes at the c.m.s. energy $\sqrt{s}=250$ GeV.

EW Scheme	$\mathit{\alpha }\left(0\right)$	${\mathit{G}}_{\mathit{\mu }}$	$\mathit{\alpha }\left({\mathit{M}}_{\mathbf{Z}}^2\right)$
${\sigma }_{\mathrm{Born}}$, fb	$225.59\left(1\right)$	$240.27\left(1\right)$	$254.49\left(1\right)$
${\sigma }_{\mathrm{Born}+\mathrm{PW}}$, fb	$231.17\left(1\right)$	$232.49\left(1\right)$	$231.01\left(1\right)$
$\delta {\sigma }_{\mathrm{PW}}$, fb	$5.58\left(1\right)$	$-7.78\left(1\right)$	$-23.48\left(1\right)$
${\delta }_{\mathrm{PW}}$, %	$2.47\left(1\right)$	$-3.24\left(1\right)$	$-9.22\left(1\right)$

, the results of the Born cross-sections, the sum of the Born and pure weak (PW) contributions as well as the relative corrections $\delta$ (%) for the c.m.s. energies $\sqrt{s}=240$ and 250 GeV in the $\alpha \left(0\right)$, $G_{\mu }$ and $\alpha \left(M_Z^2\right)$ EW schemes are presented (We define the pure weak contribution as the difference between the complete one-loop electroweak correction and the pure QED part of it). The $\alpha \left(M_{\mathrm{Z}}^2\right)$ = 1/129.02 value was used in the calculations. As seen the corrections in the $\alpha \left(0\right)$ scheme are positive and equal to 2.72% for the c.m.s. energy $\sqrt{s}=240$ GeV and 2.47% for $\sqrt{s}=250$ GeV. The calculations in the $G_{\mu }$ scheme reduce RC to about 5–6 %, they become negative and equal to −2.99% for the c.m.s. energy $\sqrt{s}=240$ GeV and −3.24% for the c.m.s energy $\sqrt{s}=250$ GeV. In the case of the $\alpha \left(M_Z^2\right)$ scheme, RCs get even more negative and achieve the value −8.97% and −9.22% for the c.m.s. energies $\sqrt{s}=240$ GeV and $\sqrt{s}=250$ GeV, respectively. These results show that there is no most suitable EW scheme of calculations for minimizing the value of the pure weak corrections for the $e^{+}{e}^{-}\to ZH$ reaction. However, the sensitivity to the choice of input EW scheme is reduced for the Born+PW cross-sections compared to the Born one. In [22,23], the mixed QCD and EW NNLO corrections were considered and a further reduction of the EW scheme dependence was observed.

In

Table 5. Comparison between results with order-by-order and exponentiated structure functions. Only pure photonic corrections are taken in account. Here $R_i={\sigma }_i/{\sigma }_{\mathrm{Exp}\mathrm{Mul}}^{\left(3\right)}$, $i=(\mathrm{LLA},\mathrm{Exp}\mathrm{Add},\mathrm{Exp}\mathrm{Mul})$, ${\sigma }_{\mathrm{Exp}\mathrm{Add}}$ calculated with the electron structure functions taken in the additive exponentiated form [18] and ${\sigma }_{\mathrm{Exp}\mathrm{Mul}}$ in the multiplicative exponentiated form [19].

	${\displaystyle \sum _{\mathit{n}=1}^{\mathit{N}}}\mathcal{O}\left({\mathit{\alpha }}^{\mathit{n}}{\mathit{L}}^{\mathit{n}}\right)$
$\mathbf{N}$	1	2	3	4
$\sqrt{s}=240$ GeV
$R_{\mathrm{LLA}}$	$0.9934$	$0.9993$	$1.0002$	$1.0001$
$R_{\mathrm{Exp}\mathrm{Add}}$	$0.9975$	$1.0054$	$1.0000$
$R_{\mathrm{Exp}\mathrm{Mul}}$	$0.9996$	$1.0000$	$1.0000$
$\sqrt{s}=250$ GeV
$R_{\mathrm{LLA}}$	$1.0001$	$0.9990$	$1.0000$	$0.9999$
$R_{\mathrm{Exp}\mathrm{Add}}$	$0.9969$	$1.0002$	$1.0000$
$R_{\mathrm{Exp}\mathrm{Mul}}$	$0.9994$	$1.0000$	$1.0000$

, we verified the difference between order-by-order and exponentiated (“additive” according to the prescription of Kuraev and Fadin [8] and “multiplicative” proposed by Jadach and Ward [24]) realizations of the electron structure function. Results are shown up to $\mathcal{O}\left({\alpha }^3{L}^3\right)$ finite terms for exponentiated forms and up to $\mathcal{O}\left({\alpha }^4{L}^4\right)$ for order-by-order calculation. It can be seen that result using multiplicative exponentiated form converges faster. However, taking into account four orders in the order-by-order calculation is enough to reach the ${10}^{-4}$ accuracy.

4. Conclusions

We considered the contributions due to the QED initial-state radiation (photons and pairs) to the higgsstrahlung process. Their impact has been analyzed order by order. The complete one-loop electroweak one-loop corrections were presented. Higher-order ISR QED contributions were calculated within the leading logarithmic approximation. The known expressions for contributions of the collinear electron structure function of the orders $\mathcal{O}\left({\alpha }^n{L}^n\right),n=$ 2–4 for photons and pairs were used. These corrections are known to be very important in the case of resonances, e.g., at the Z-boson peak studied at LEP. We would like to emphasize that higher-order QED ISR corrections can be large not only at resonances but also near the reaction thresholds. Please note that the cross-section of this process has a peak at the threshold.

By looking at the magnitude of the complete one-loop electroweak and higher-order LLA QED corrections, we can estimate the theoretical uncertainty and define what other contributions should be taken into account. Specifically, a safe estimate of the theoretical uncertainty in EW and LLA RC can be derived by variation of the EW scheme and the factorization scale, respectively. One can see that to meet the high precision of future experiments, we need to go beyond the approximation explored here. At least the next-to-leading QED ISR logarithmic corrections should also be taken into account. One needs to improve the uncertainty in pure weak contributions. That can be done by taking into account higher-order EW and mixed QCD and EW effects in the Z boson propagator and vertices. Note also that corrections for the whole processes with different decay modes of Z and Higgs bosons should be evaluated.

For the one permille precision tag relevant for future studies of the higgsstrahlung process, we see that there is a good agreement between the order-by-order results and the known exponentiated QED LLA corrections [18,19]. So either approach can be used. Presumably, the exponentiated one is more suitable for Monte Carlo simulations, while the order-by-order one can be used for benchmarks and cross-checks.

The numerical results presented here were obtained by means of the Monte Carlo generator ReneSANCe [25] and MCSANCee integrator [26] which allow evaluation of arbitrary differential cross-sections. These computer codes can be downloaded from the SANC project homepage sanc.jinr.ru, (accessed on 29 June 2021) and the ReneSANCe HEPForge page renesance.hepforge.org, (accessed on 29 June 2021).