Impact of NNLO QED corrections on lepton-proton scattering at MUSE

Abstract

We present the complete next-to-next-to-leading order (NNLO) pure pointlike QED corrections to lepton-proton scattering, including three-photon-exchange contributions, and investigate their impact in the case of the MUSE experiment. These corrections are computed with no approximation regarding the energy of the emitted photons and taking into account lepton-mass effects. We contrast the NNLO QED corrections to known next-to-leading order corrections, where we include the elastic two-photon exchange (TPE) through a simple hadronic model calculation with a dipole ansatz for the proton electromagnetic form factors. We show that, in the low-momentum-transfer region accessed by the MUSE experiment, the improvement due to more sophisticated treatments of the TPE, including inelastic TPE, is of similar if not smaller size than some of the NNLO QED corrections. Hence, the latter have to be included in a precision determination of the low-energy proton structure from scattering data, in particular for electron-proton scattering. For muon-proton scattering, the NNLO QED corrections are considerably smaller.

1 Introduction

The scattering of electrons and muons off protons has been used for decades to obtain information on the structure of the proton. Still, in the regime of low energies, where the quark content of the proton is not yet resolved and the scattering is described with the help of form factors, there are several open questions and discrepancies, see [1] for a recent review. In view of this unsatisfactory situation it is important to revisit the theoretical aspects related to the extraction of form factors of the proton, with careful consideration of all effects that influence the differential distributions of the final-state particles. In addition to the uncertainty budget of radiative corrections due to hadronic contributions, dominated by the two-photon-exchange (TPE), this also includes standard QED corrections. The latter can lead to additional real photons in the final state, and a precise confrontation of theory with experiment needs to specify how such radiative events are treated.

The analyses carried out so far have taken into account QED corrections at next-to-leading order (NLO), often with additional approximations [2,3,4,5,6,7,8]. However, perturbative calculations of QED corrections with pointlike particles to fully differential cross sections have now reached a maturity that allows to obtain complete next-to-next-to-leading order (NNLO) corrections to $2\rightarrow {2}$ processes [9,10,11,12,13,14,15]. In the following, we refer to these types of corrections as pure QED corrections. These computations can be done including mass effects and without making any approximation on the energy range of the emitted photons. This provides an opportunity to obtain unprecedented accuracy for the pure pointlike QED part of the low-energy lepton-scattering processes.

The presence of non-pointlike hadrons poses an additional challenge. While the emission of a single photon from an on-shell proton line can be described by two electromagnetic form factors, Feynman diagrams with more complicated topologies, e.g. involving hadronic intermediate states, are more difficult to describe. Experimental analyses used to include the TPE as evaluated in the article by Mo-Tsai [2] or the later article by Maximon-Tjon [3]. That is the elastic TPE, which has a proton in the intermediate state, as well as real radiation (bremsstrahlung), both in the limit of soft photons. The precision of modern scattering experiments – take for example the A1 [16] and initial-state-radiation (ISR) [17, 18] experiments at MAMI or the PRad [19] experiment at JLab – required to go beyond that approximation and consider a more complete treatment of TPE and real radiation [6, 20]. Corrections beyond the soft-photon approximation are sometimes referred to as “hard TPE” and hard-photon radiation, respectively [1]. In the following, we refer to corrections from diagrams with exchange of two virtual photons, shown in (8a), as virtual TPE, and to the interference of one-photon-exchange (OPE) diagrams with a single bremsstrahlung photon radiated from the lepton and proton line, respectively, shown in (8b), as real TPE.

At low energies, the TPE contributions cannot be computed in perturbative QCD directly. They need to be modeled or, preferably, evaluated without model dependence in an effective-field theory framework [21,22,23] or through the use of dispersion relations. The latter require further experimental input, see [24,25,26] for a selection of recent data-driven evaluations. Since the first works suggested an insufficiently precise description of the hard TPE as the origin of the discrepancy between form factor extractions from unpolarised and polarisation-transfer measurements [27,28,29,30], a vast literature on how to obtain and improve virtual TPE contributions appeared, see [31,32,33,34] for reviews focusing solely on virtual TPE in lepton-proton scattering.

The main aim of this investigation is not to improve the predictions for TPE as such but, rather, to critically assess the impact of TPE corrections available in the literature, relative to other corrections to lepton-proton scattering. To this end, we combine a simplified implementation of the TPE corrections with state-of-the-art NNLO QED corrections.

We focus our application on the high-precision muon-scattering experiment MUSE [35, 36], which uses a beam of electrons and muons of both charges ($e^+$ and $\mu ^+$ as well as $e^-$ and $\mu ^-$), with three different beam momenta¹$p_\text {beam}=115, 153, 210$ MeV. Its aim is to compare extractions of the proton charge radius from electron and muon scattering, respectively, obtained with the same experimental setup, and to experimentally determine TPE corrections making use of both beam polarities. The MUSE kinematics is limited to the low momentum-transfer region (0.08 GeV$^2$ for $p_\text {beam}=210$ MeV), where the TPE corrections are dominated by the elastic TPE, while the inelastic TPE is smaller than the anticipated accuracy of the MUSE cross-section measurements [38]. Therefore, as a reasonable first approximation in the MUSE kinematics, we implement a simple model for the elastic TPE contribution and neglect the inelastic part.

All considered corrections are implemented in the McMule framework [10]. This goes beyond the NLO radiative corrections from [6] applied in the recent MUSE analysis of instrumental uncertainties [37]. In particular, we adapt the recent NNLO computation for muon-electron scattering [14] to obtain the NNLO QED corrections for lepton-proton scattering for pointlike protons in a fully differential way. It is the first time that these corrections (including three-photon exchange contributions) are taken into account, and we assess their relevance relative to variations in the treatment of TPE corrections.

In Sect. 2 we will give a detailed description of the contributions that are included in our calculation. This allows us to present in Sect. 3 results for MUSE with $p_\text {beam}=210$ MeV and study the impact of NNLO QED corrections. Our conclusions and an outlook towards further work will be presented in Sect. 4.

2 Calculation

In order to maximally exploit the technical progress in the computation of QED corrections, we take as the starting point a pointlike interaction $ie \,q_p\,\gamma ^\mu$ of the photon with the proton, and we will call this the pure QED contribution. We introduce the charge of the proton $q_p=1$ in units of e for bookkeeping purposes. The non-pointlike structure of the proton will be taken into account by additional contributions, denoted by $F\times \delta ^\mu$. Again F has been introduced for bookkeeping purposes. Thus, for the photon–proton interaction we write

(1)

with the four-momentum of the photon q, the spacelike virtuality of the photon $q^2=-Q^2<0$, the Dirac and Pauli form factors of the proton, $F_1(Q^2)$ and $F_2(Q^2)$, the proton mass M, and our notation for the antisymmetric combination of Dirac matrices $\sigma ^{\mu \nu }=\frac{i}{2}\left( \gamma ^\mu \gamma ^\nu -\gamma ^\nu \gamma ^\mu \right)$. Note that for real photons, the form factors are normalised through their charge and anomalous magnetic moment $\kappa$: $F_1(0)=q_p=1$ and $F_2(0)=\kappa$. In what follows we will describe in detail which contributions we include, up to and including NNLO.

Starting at LO we obtain the matrix element (squared)

(2)

by computing the tree-level amplitude of the $2\rightarrow {2}$ lepton-proton process with the full photon-proton vertex, depicted as grey blobs in accordance with (1). The subscript n indicates the number of final-state particles, i.e. $n=2$ for the process considered in this paper. In the argument of the amplitude $\mathcal {A}^{(0)}_n$ we indicate that there is a single power of the coupling of the photon to the lepton, where $q_\ell =\pm {1}$ is the charge of the positron or electron in units of e, and that the full photon-proton vertex (1) with arbitrary form factors $F_1$ and $F_2$ is included. We suppress the dependence on the external momenta but use the convention

$$\begin{aligned} \ell (p_1)\, p(p_2) \rightarrow \ell (p_3)\, p(p_4)\ + \{\gamma _i(k_i)\} \end{aligned}$$

(3)

with up to two additional photons in the final state and either lepton polarity.

Thus, the unpolarised tree-level cross section is obtained by integrating (2) over the two-body phase space $d\Phi _n$ and including the standard flux and spin average factors

$$\begin{aligned}&d\sigma ^{(0)}(q^2_\ell ,F^2)= \frac{1}{2s} \frac{1}{4} \int d\Phi _n \mathcal {M}^{(0)}_n(q^2_\ell ,F^2) S(p_3,p_4) \, . \nonumber \\ \end{aligned}$$

(4)

The differential nature of the computation is coded in the measurement function $S(p_3,p_4)$ that allows to include arbitrary cuts on the final states.

At NLO, virtual and real corrections contribute to the cross section. This leads to ultraviolet (UV) and infrared (IR) divergences. Both types of singularities are regularised dimensionally. For the UV divergences, we use the on-shell renormalisation scheme for the masses and the coupling. The IR singularities cancel when combining real and virtual corrections, as discussed in [1]. We perform the phase-space integrations numerically, using the FKS$^\ell$ subtraction method [39] to achieve this cancellation for arbitrary IR-safe observables.

The so-called leptonic corrections consist of Feynman diagrams with additional photons solely attached to the lepton line. This is a gauge-invariant subset of the complete NLO correction, and corresponds to the OPE approximation. In the case of the electron, these corrections are expected to dominate due to collinear emission. This results in large logarithms of the form $\log (m_e^2/E^2)$ where the energy scale of the process E is much larger than the electron mass $m_e$. Representative diagrams for the leptonic NLO corrections are

(5a)

(5b)

The virtual corrections $d\sigma ^{(1)}_v(q^4_\ell ,F^2)$ are obtained by integration of the one-loop matrix element $\mathcal {M}^{(1)}_n(q^4_\ell ,F^2)$ over the two-body phase space, whereas for the real corrections $d\sigma ^{(1)}_r(q^4_\ell ,F^2)$ we have to integrate the matrix element $\mathcal {M}^{(0)}_{n+1}(q^4_\ell ,F^2)$ describing real radiation over the three-body phase space $d\Phi _{n+1}$. All these calculations can be performed with arbitrary form factors [2, 3] using standard methods, and we obtain the leptonic NLO corrections

$$\begin{aligned} d\sigma ^{(1)}(q^4_\ell ,F^2)&= d\sigma ^{(1)}_v(q^4_\ell ,F^2) + d\sigma ^{(1)}_r(q^4_\ell ,F^2) \end{aligned}$$

(6)

to the cross section.

In analogy to the leptonic corrections, the protonic corrections include emission solely from the proton line. Technically speaking, they are also OPE contributions. Note that we do not absorb the pure QED virtual corrections, i.e. the vertex corrections, into the form factors, since they are IR divergent. In general, we consider all QED reducible contributions as independent from the form factors. As we will see, in practice all this has very limited impact, since these corrections are very small compared to the leptonic OPE, even for $\ell =\mu$. This is due to the lack of logarithmic enhancements since collinear radiation only results in logarithms of the form $\log (M^2/E^2)$ with $M^2 \approx E^2$ for the energy scale considered here. Thus, in a standard OPE approach, these corrections are often neglected. In the following, we take the protonic corrections into account in the pointlike proton approximation. Again, we have virtual and real corrections

(7a)

(7b)

where we only show a single diagram for illustrative purposes.

The NLO corrections that involve additional photon couplings to both, the lepton and proton line, we call mixed corrections or simply TPE. As a first approximation to the TPE at NLO, we consider the elastic contribution which is due to an intermediate proton. For the virtual TPE correction this results in box (and crossed box) diagrams whereas for the real TPE corrections the intermediate proton is between the exchange photon and the real photon. Concretely, we take into account

(8a)

(8b)

We include the real TPE contribution $d\sigma _r^{(1)}(q^3_\ell ,F^3)$ by integrating $\mathcal {M}^{(0)}_{n+1}(q^3_\ell ,F^3)$ over the three-body phase space $d\Phi _{n+1}$. Combining this with the corresponding virtual TPE correction $d\sigma _v^{(1)}(q^3_\ell ,F^3)$ leads to an IR finite result for any IR safe observable.

Besides the elastic TPE with a proton intermediate state, there is the so-called inelastic TPE contribution with inelastic intermediate states

(9)

The latter are denoted by the oval blobs in the above Feynman diagrams, and may contain anything: pions, resonances like the $\Delta (1232)$, etc. In the low-$Q^2$ region, relevant for the MUSE experiment, the elastic TPE dominates, while the inelastic TPE can be neglected. Therefore, the inelastic TPE is presently not included in McMule. Again, we leave this for a future update.

In this work, the virtual TPE has been implemented through a simple hadronic model calculation of the box and crossed-box diagrams in (8a), assuming on-shell proton form factors. The same approach has been used for instance in [28, 49] for electron and muon scattering, respectively. In [50], the hadronic model calculation has been compared to a dispersive evaluation with one subtraction function. The latter involves an s-channel cut through the box diagram, thus, only needs input from on-shell form factors and does not require off-shell form factors, as the box calculation would in principle. Both approaches agree with $5\times 10^{-4}$ relative accuracy for the kinematics of muon scattering in the MUSE experiment [50], and the same quality of the approximation can be assumed for electron scattering at MUSE. Thus, this model dependence in our approach can be safely neglected. For the proton electric and magnetic Sachs form factors

$$\begin{aligned} G_E(Q^2)= & {} F_1(Q^2)-\frac{Q^2}{4\,M^2} F_2(Q^2), \end{aligned}$$

(10a)

$$\begin{aligned} G_M(Q^2)= & {} F_1(Q^2)+ F_2(Q^2), \end{aligned}$$

(10b)

we use a dipole ansatz

$$\begin{aligned} G_E(Q^2)= & {} G_D(Q^2)=\frac{G_M(Q^2)}{1+\kappa }, \qquad \text {with} \quad \nonumber \\{} & {} G_D(Q^2) =\left( \frac{\Lambda ^2}{\Lambda ^2+Q^2}\right) ^2. \end{aligned}$$

(10c)

Note that in the limit $\Lambda \rightarrow \infty$ the pure (pointlike) QED vertex is recovered. Of course, describing both the normalised electric and magnetic Sachs form factor through one single parameter $\Lambda$ is a rather naive ansatz. Furthermore, the simple dipole form can only ever be a rough approximation to any form factor. Nevertheless, the standard dipole with $\Lambda ^2=0.71$ GeV$^2$ is widely used as a reasonable first approximation to the proton form factors, and serves well our purpose to examine the relative importance of TPE and NNLO QED corrections. While we leave the implementation of the elastic TPE correction with input from modern form-factor parametrisations to a future version of McMule, we want to illustrate the impact of the form factors and the uncertainties in their description by considering various values for $\Lambda$. To motivate our choice, we consider the slopes of the Sachs form factors at $Q^2=0$, which are related to the charge and magnetisation radii of the proton

$$\begin{aligned} R_{E,M}=\sqrt{-\frac{6}{G_{E,M}(0)}\frac{dG_{E,M}(Q^2)}{dQ^2}}\Bigg \vert _{Q^2=0} \end{aligned}$$

(11)

shown in Fig. 1. Besides the standard dipole, we use $\Lambda ^2=0.86$ GeV$^2$ reproducing the small $R_M$ from [40], $\Lambda ^2=0.66$ GeV$^2$ reproducing $R_E$ as extracted with unprecedented precision from the Lamb shift in muonic hydrogen [43], and $\Lambda ^2=0.60$ GeV$^2$ reproducing the large $R_E$ from [41]. The resulting impact on the virtual TPE correction, defined usually as

$$\begin{aligned} \delta _{2\gamma }(\text {IR}) = \frac{\mathcal {M}^{(1)}_{n}(q^3_\ell ,F^3)\big |_\text {IR}}{\mathcal {M}^{(0)}_n(q^2_\ell ,F^2)}, \end{aligned}$$

(12)

is shown in Fig. 2 for electron and muon scattering. The “IR” label indicates that the omission of $\mathcal {M}^{(0)}_{n+1}(q^3_\ell ,F^3)$ in (12) requires an unphysical subtraction of the IR singularity in $\mathcal {M}^{(1)}_{n}(q^3_\ell ,F^3)$. For the purpose of Fig. 2 we have used the Maximon-Tjon subtraction [1, 3].

Fig. 1

Comparison of extractions of the proton electric and magnetic radii, $R_E$ and $R_M$, from different parametrisations of the proton Sachs form factors [16, 40,41,42], including the standard dipole with $\Lambda ^2=0.71$ GeV$^2$, to the extraction from the muonic-hydrogen Lamb shift [43], and the CODATA recommended values for $R_E$, before (CODATA ’14 [44]) and after (CODATA ’18 [45]) inclusion of the data from muonic hydrogen. Note that the displayed Bernauer results [16] are including hard TPE corrections from [46] (solid) and [32, 47] (dashed), respectively. On the additional axes, we show which value of $\Lambda _{E,M}^2$ would reproduce the radii if a dipole ansatz was assumed for the form factors

Going towards larger values of $Q^2$, excited intermediate states, e.g., resonances [26], eventually do lead to sizeable contributions. In [38], the inelastic TPE correction has been estimated through a dispersive approach for near-forward kinematics, relating it to forward doubly-virtual Compton scattering amplitudes, which are in turn reconstructed dispersively with empirical input for the unpolarised proton structure functions. In their estimate, the inelastic TPE contribution, $\delta _{2\gamma }\sim 5\times 10^{-4}$, is smaller than the anticipated $1 \,\%$ accuracy of the cross-section measurements at MUSE. Coincidentally, for ep scattering, our evaluation of the elastic TPE with a dipole parameter of $\Lambda ^2=0.86$ GeV$^2$ is very close to their prediction for the total TPE [38]. This can be seen from Fig. 2, where we compare our spread of results approximating the elastic TPE to the dispersive evaluation (solid cyan line) [38] and an empirical extraction (solid violet line and error band) [16] of the total TPE. Note that the elastic TPE included in [38] agrees with our result for the standard dipole (short-dashed blue line). The uncertainty on their inelastic TPE is small on the scale of the total TPE, and thus, omitted from the figure.

Fig. 2

Comparison of virtual TPE corrections, $\delta _{2\gamma }$, to $e^- p$ (left column) and $\mu ^- p$ (right column) scattering for three different beam momenta envisaged by the MUSE experiment [48]: $p_\text {beam}=115, 153, 210$ MeV. The soft singularities are subtracted using the Maximon-Tjon prescription [3]. Shown are the elastic TPE from our box model calculation with proton dipole form factors and different values of $\Lambda ^2=0.60, 0.66, 0.71$ and 0.86 GeV$^2$ (orange dotted, red dot-dashed, blue short-dashed and pink long-dashed lines), compared to the theoretical prediction for the total TPE from [38] (solid cyan line), and the empirical extraction of the total TPE from [16] (solid violet line with error band)

All contributions discussed so far were NLO. Moving towards NNLO, we start again with the leptonic or OPE corrections. This gauge invariant subset of NNLO corrections has been computed [9, 10] for any choice of form factors. It consists of double-virtual, real-virtual, and double-real corrections

$$\begin{aligned} d\sigma ^{(2)}(q^6_\ell ,F^2)&= d\sigma ^{(2)}_{vv}(q^6_\ell ,F^2) + d\sigma ^{(2)}_{rv}(q^6_\ell ,F^2)\nonumber \\ {}&\quad + d\sigma ^{(2)}_{rr}(q^6_\ell ,F^2) \end{aligned}$$

(13)

which are obtained by integrating the two-loop matrix element $\mathcal {M}^{(2)}_n(q^6_\ell ,F^2)$ over $d\Phi _n$, the one-loop matrix element $\mathcal {M}^{(1)}_{n+1}(q^6_\ell ,F^2)$ over $d\Phi _{n+1}$, and the tree-level matrix element $\mathcal {M}^{(0)}_{n+2}(q^6_\ell ,F^2)$ over $d\Phi _{n+2}$, respectively. For any IR-safe observable, the IR singularities of the individual parts in (13) cancel in the sum. Representative Feynman diagrams of the various matrix elements are

(14a)

(14b)

(14c)

The one-loop amplitude squared, cf. the second line in (14a), is included in the two-loop matrix element. We note that some NLO diagrams for the process of $\ell p \rightarrow \ell p \gamma$, corresponding to an IR-finite subset of the $\ell p \rightarrow \ell p$ process at NNLO, have been previously included in [20] in approximate ways. The full set of leptonic NNLO corrections depicted in (14a), (14b), and (14c) has been computed in [9] with a slicing approach and later with the McMule framework in [10]. The two results disagree substantially and a corresponding discussion can be found in [10].

The leptonic corrections are expected to be dominant, at least for the case $\ell =e$, since they contain hard collinear emission from the electron line. This leads to large logarithms. As we will see, the size of these corrections depends crucially on the precise definition of the observable. More concretely, the way additional photon radiation is treated in the experiment will have a decisive impact. Thus, these corrections have to be under control for empirical extractions of form factors and TPE effects.

The leptonic corrections are technically the most simple NNLO corrections. Going beyond OPE, we have to consider one-loop pentagon diagrams (for the real-virtual corrections) and, more challenging, a set of topologically non-trivial two-loop diagrams, including (crossed) double-box diagrams. In the language of lepton-proton scattering, they correspond to three-photon exchange contributions. With current techniques, it is not possible to do such a computation including form factors. Hence, for all NNLO corrections beyond OPE we use the approximation of a pointlike proton $\Lambda \rightarrow \infty$. In this case, the NNLO corrections to lepton-proton scattering can be obtained from those of muon-electron scattering, with adapted masses. The latter have been computed [14, 51] using the two-loop integrals of [52], as well as OpenLoops [53] and Package-X [54] for the one-loop amplitudes. They can be split into gauge-invariant parts consisting of terms with a fixed power of $q_\ell ^n$ and $q_p^m$ with $n+m=8$. As an example, we illustrate the virtual contributions to $d\sigma ^{(2)}(q_\ell ^4,q_p^4)$ which requires the two-loop matrix element

(15a)

Of course, real-virtual and double-real corrections have to be considered as well. Representative examples for $d\sigma ^{(2)}_{rv}(q_\ell ^5,q_p^3)$ and $d\sigma ^{(2)}_{rr}(q_\ell ^5,q_p^3)$ are

(16a)

(16b)

All contributions considered so far are collectively called photonic corrections. In addition, there are vacuum polarisation contributions. We include electron, muon, tau loops in the vacuum polarisation $\Pi$, and collectively refer to them as fermionic corrections. Note that we also include hadronic contributions in $\Pi$, however, they are about a factor 100 smaller than the fermionic contributions. Vacuum polarization starts to contribute at NLO through virtual effects. At NNLO, there are virtual and real fermionic corrections to be included. In analogy to the other corrections, we use a form factor for the OPE contributions and a pointlike proton interaction for TPE. Sample diagrams for the virtual corrections are

(17a)

(17b)

(17c)

There are also contributions with either two one-loop insertions of $\Pi$ or a single insertion of a two-loop vacuum polarisation. They contribute to virtual corrections only and we denote them collectively by $d\sigma ^{(2)}(q_\ell ^2,\Pi ^2,F^2)$. In the case of leptons, the analytic form of the two-loop vacuum polarisation from [55] is used, whereas for hadronic loops, the Fortran library alphaQED [56] is employed. We do not include real corrections with an additional $e^+\,e^-$ pair in the final state [57]. This is a measurably different process. However, depending on the details of the experimental analysis, this process can contribute to lepton-proton cross sections at NNLO.

To summarise, our results for the NNLO cross section

$$\begin{aligned} d\sigma _2&=d\sigma ^{(0)}+d\sigma ^{(1)}+d\sigma ^{(2)} \end{aligned}$$

(18)

include

$$\begin{aligned} d\sigma _0&=d\sigma ^{(0)}(q_\ell ^2,F^2), \end{aligned}$$

(19a)

$$\begin{aligned} d\sigma ^{(1)}&=d\sigma ^{(1)}(q_\ell ^4,F^2) +d\sigma ^{(1)}(q_\ell ^3,F^3) \nonumber \\&\quad + d\sigma ^{(1)}(q_\ell ^2,q_p^4) +d\sigma ^{(1)}(q_\ell ^2,\Pi ,F^2), \end{aligned}$$

(19b)

$$\begin{aligned} d\sigma ^{(2)}&=d\sigma ^{(2)}(q_\ell ^6,F^2) \nonumber \\ {}&\quad + \Bigg (\sum _{j=3}^{5} d\sigma ^{(2)}(q_\ell ^j,q_p^{8-j}) \Bigg ) + d\sigma ^{(2)}(q_\ell ^2,q_p^6)\nonumber \\ {}&\quad +\Big (d\sigma ^{(2)}(q_\ell ^4,\Pi ,F^2) + d\sigma ^{(2)}(q_\ell ^2,\Pi ^2,F^2)\Big )\nonumber \\ {}&\quad +d\sigma ^{(2)}(q_\ell ^3,\Pi ,q_p^3) +d\sigma ^{(2)}(q_\ell ^2,\Pi ,q_p^4), \end{aligned}$$

(19c)

where all parts are individually gauge independent. All contributions proportional to $F^2$ can easily be computed with arbitrary form factors. The term $d\sigma ^{(1)}(q_\ell ^3,F^3)$ has been computed using a dipole ansatz for the electromagnetic form factors. For the remaining terms we use the pointlike proton approximation.

In Sect. 3 we will present results for different values of the dipole parameter $\Lambda$ in (10c), including $\Lambda =\infty$ for pointlike protons. In order to indicate the dependence on $\Lambda$, we will use the compact notation

$$\begin{aligned} d\sigma _0^\Lambda&=d\sigma ^{(0)}(q_\ell ^2,F(\Lambda )^2), \end{aligned}$$

(20a)

$$\begin{aligned} d\sigma ^{(1)\Lambda }&= d\sigma _\ell ^{(1)\Lambda } +d\sigma _{x}^{(1)\Lambda } +d\sigma _{p}^{(1)\infty } +d\sigma _{\Pi }^{(1)\Lambda }, \end{aligned}$$

(20b)

$$\begin{aligned} d\sigma ^{(2)\Lambda }&= d\sigma _\ell ^{(2)\Lambda } +d\sigma _x^{(2)\infty } +d\sigma _p^{(2)\infty } +d\sigma _{\ell \Pi }^{(2)\Lambda } \nonumber \\&\quad +d\sigma _{x\Pi }^{(2)\infty } +d\sigma _{p\Pi }^{(2)\infty }, \end{aligned}$$

(20c)

where the terms in (20) are in one-to-one correspondence with those of (19). The labels $\ell \in \{e,\mu \}$, p, and x stand for leptonic (i.e. electronic or muonic), protonic, and mixed corrections. The terms $\sim q_\ell ^2\Pi ^2 F^2$ have been absorbed into $d\sigma _{\ell \Pi }^{(2)\Lambda }$ as a matter of convention. If the proton is treated pointlike, we set $\Lambda =\infty$ in the notation. Otherwise, we use the label $\Lambda \in \{60,71,86\}$ to indicate the value of the dipole parameter that has been used, where e.g. the label $\Lambda =60$ corresponds to $\Lambda ^2=0.60\,\textrm{GeV}^2$.

3 Results and discussion

This section presents results for lepton-proton scattering, tailored to the characteristics of the MUSE experiment [35, 36]. Particular emphasis is given to the impact of NNLO pure QED corrections compared to hadronic effects at NLO, focusing mainly on the elastic TPE discussed in Sect. 2. The code employed for this study is publicly available at https://gitlab.com/mule-tools/mcmule and the whole set of results can be found in the relevant directory of the McMule user library https://mule-tools.gitlab.io/user-library/ along with user, menu and configuration files, and the Python code that generates the plots in the paper [58]. The production runs employed version v0.5.0 of the McMule public release.

The kinematics of the process is defined by the momenta in (3) together with the lepton and the proton mass, $m_\ell$ and M. Both polarities of the lepton are considered. For the purpose of illustration, we consider an incoming lepton of momentum

$$\begin{aligned} p_{\textrm{beam}}=|\vec {p}_1|=210\,\textrm{MeV} \end{aligned}$$

(21)

scattering off a proton at rest. This is consistent with one of the MUSE setups, and corresponds to a centre-of-mass energy of $\sqrt{s}\approx 1.2$ GeV. The experimental setup defines a window for the lepton scattering angle, $20^\circ< \theta _\ell < 100^\circ$, and for the lepton final-state momentum

$$\begin{aligned} |\vec {p}_3| > 15\,\textrm{MeV} \equiv p_{\textrm{min}}\,. \end{aligned}$$

(22)

All the results use the input parameters [59]

$$\begin{aligned} \alpha&= 1/137.035999084 \,,&m_e&= 0.510998950 \,\,\textrm{MeV} \,, \\ m_\mu&= 105.658375 \,\,\textrm{MeV} \,,&M&= 938.2720813 \,\,\textrm{MeV} \,. \end{aligned}$$

The most recent version alphaQEDc19 of the Fortran library alphaQED is used for the evaluation of diagrams with hadronic loop insertions.

In order to measure elastic scattering, kinematical cuts have to be applied to suppress radiative events. In the case of the MUSE experiment this is done by means of a forward-angle calorimeter [37]. In the following, results are presented for two different scenarios, S0 and S1, depending on whether an additional inelasticity cut is enforced. In our simulation, the energy of photon(s) emitted in the lab frame within a 100 mrad angle is cumulated into $E_\gamma$. In scenario S1, if $E_\gamma > 0.4\,p_{\textrm{beam}}$ the corresponding event is removed from the analysis. The kinematical details discussed above and used in our scenarios are summarised in Table 1.

Table 1 Kinematical scenarios analysed in the McMule prediction

Fig. 3

Complete NLO+NNLO corrections for ep (left panel) and $\mu p$ (right panel) scattering, normalised to the Born cross sections. Shown are both kinematic scenarios S0 and S1. Since corrections due to real-photon emissions are included, the distributions differ whether they are plotted as functions of $Q^2_\ell$ or $Q^2_p$

Figure 3 illustrates the impact of this forward-angle inelasticity cut on $E_\gamma$. It depicts the sum of the NLO and NNLO corrections normalised to the Born cross section as a function of the leptonic momentum transfer $Q_\ell ^2 = -(p_1-p_3)^2$ as well as its protonic counterpart $Q_p^2=-(p_2-p_4)^2$. Results for both scenarios S0 and S1 are shown. In the presence of radiation we have $Q_\ell ^2 \ne Q_p^2$ and the deviation of the two curves can be taken as a measure of inelastic effects. In the case of $\mu p$ scattering hard radiation is not collinearly enhanced due to the larger lepton mass. For ep scattering, on the other hand, sizable deviations can be observed. We therefore restrict the discussion to this more interesting case.

For small momentum transfer, the S1 distributions w.r.t. $Q_e^2$ and $Q_p^2$ are almost identical. Hence, the $E_\gamma$ cut is able to remove most hard photon radiation in this region. This is not the case for larger momentum transfer where the $Q_e^2$ and $Q_p^2$ curves start to deviate. This behaviour can be understood as follows. Small momentum transfer corresponds to forward scattering of the lepton. In this case both initial-state as well as final-state collinear radiation is emitted in forward direction and therefore removed by the $E_\gamma$ cut. For larger $Q^2$, on the other hand, final-state collinear radiation is not forwardly directed and is thus not vetoed. Nevertheless, the inelasticity cut still removes initial-state collinear radiation. However, this effect seems to be very small since the S0 and S1 scenarios approach each other in this region. Thus, final-state radiation dominates the inelastic effects for larger momentum transfer irrespective of the cut on $E_\gamma$.

The order-by-order contributions, $\sigma ^{(i)}$, to the LO, NLO, and NNLO integrated cross section are shown in Tables 2 and 3 for ep and $\mu p$ scattering, respectively. The various contributions are denoted as in (20). The results are presented for both kinematical scenarios and for different values of $\Lambda$ entering the dipole ansatz of the proton form factors (10c). In particular, $\Lambda ^2\in \{0.60,\,0.71,\,0.86\}$ GeV$^2$, while the label $\Lambda \infty$ stands for a pure pointlike QED photon-proton vertex. The NNLO corrections have been evaluated only for the point-like proton ($\Lambda \infty$) and the standard dipole form factors ($\Lambda 71$). We reiterate that inelastic contributions to TPE are not taken into account in the calculation.

Table 2 Integrated cross sections for ep scattering, for both S0 and S1 scenarios, at LO, NLO, and NNLO. All digits shown are significant. $\Lambda \infty$ denotes pure QED contributions with a pointlike proton, $\{\Lambda 60,\,\Lambda 71,\,\Lambda 86\}$ stand for proton finite-size corrections with proton form factors modeled through a dipole ansatz and $\Lambda ^2$ set to $\{0.60,\,0.71,\,0.86\}$ GeV$^2$, respectively. When applicable, the different contributions with positive and negative electrons are shown

Table 3 Same as Table 2 but for $\mu p$ scattering

In order to better grasp the impact of different contributions to ep and $\mu p$ scattering, Figs. 4 and 5 show the same corrections presented in Tables 2 and 3 as bar plots, for both kinematical scenarios. For each contribution, photonic and fermionic corrections are plotted with a different colour. Each contribution with a superscript “71” indicates the corresponding pure QED contribution with the additional inclusion of proton from factors, using the dipole ansatz in (10c) with $\Lambda ^2=0.71$ GeV$^2$. For the LO and NLO contributions, a black band represents the variation obtained with $0.60\,\textrm{GeV}^2< \Lambda ^2 < 0.86\,\textrm{GeV}^2$. When a black band refers to a sum of photonic and fermionic corrections, this is plotted as the square root of the quadrature sum of the two contributions. In the case of ep scattering, the black bands for the NLO electronic correction are particularly small compared to the scale of the plot, for both S0 and S1.

Fig. 4

Integrated cross sections for ep scattering, for both S0 (left panel) and S1 (right panel), at LO, NLO, and NNLO. Yellow bars indicate contributions with fermionic loop insertions, blue bars indicate photonic contributions. The fermionic contributions $\sigma _{\Pi }^{(1)\Lambda }$, $\sigma _{e\Pi }^{(2)\Lambda }$, $\sigma _{x\Pi }^{(2)\infty }$, and $\sigma _{p\Pi }^{(2)\infty }$ are individually combined with $\sigma _{e}^{(1)\Lambda }$, $\sigma _{e}^{(2)\Lambda }$, $\sigma _{x}^{(2)\infty }$, and $\sigma _{p}^{(2)\infty }$, respectively. $\Lambda ^2=0.71$ GeV$^2$ is taken as the reference value for the dipole parameter. Black bands denote variations due to $\Lambda ^2\in [0.60,\, 0.86]\,\textrm{GeV}^2$. Note that the black bands should not be interpreted as the uncertainty estimate of our theory prediction. They merely illustrate the impact of the form factors and their possible uncertainties, assuming our naive dipole ansatz. Further uncertainties, e.g., due to higher-order corrections, model dependence or missing inelastic TPE, are not shown

Fig. 5

Same as Fig. 4 but for $\mu p$ scattering

The discussion about ep and $\mu p$ scattering somewhat differs because of the different lepton mass. We start our analysis with the former. As touched upon in Sect. 2, collinear photon emission introduces logarithms of the form $\log (m_i^2/E^2)$ where $m_i \in \{m_e,M\}$ depending on whether the photon is emitted from the electron or the proton line. In the former case, the logarithm is large and the corresponding contribution enhanced. This explains the hierarchy between the different photonic NLO contributions in Table 2: the electronic OPE corrections are dominant compared to TPE, and even more compared to OPE protonic corrections. If a cut that limits hard forward-angle radiation (S1) is applied, the collinear enhancement is reduced and the hierarchy is less pronounced.

In order to illustrate the relative importance of various corrections in $e^- p$ scattering w.r.t. the TPE, we will use $\Lambda ^2=0.71$ GeV$^2$ for the reference form factor and $\sigma _{x}^{(1){71}}$ as normalisation. Electronic and fermionic NLO corrections by far outweigh the TPE,

$$\begin{aligned} \frac{\{\sigma _e^{(1)71}, \sigma _\Pi ^{(1)71}, |\sigma _p^{(1)}|\}}{\sigma _{x}^{(1){71}}} \approx {\left\{ \begin{array}{ll} \{40, 3.5, 0.01\} &{} [\texttt{S0}]\\ \{6, 3.5, 0.01\} &{} [\texttt{S1}] \end{array}\right. } , \quad \end{aligned}$$

(23)

even in the kinematical scenario S1 in which the electronic contribution is reduced by almost a factor 7 due to the cut on $E_\gamma$. Protonic corrections, on the other hand, are much smaller, justifying their neglect or approximate (pointlike proton) inclusion.

The impact of form factor insertions at LO and NLO can be quantified, for both kinematical scenarios, as

$$\begin{aligned} \begin{aligned}&\frac{|\sigma _0^{\infty }-\sigma _0^{71}|}{\sigma _{x}^{(1){71}}} \approx 10 \,, \\ {}&\quad \frac{\{|\sigma _e^{(1)\infty }-\sigma _e^{(1){71}}|, |\sigma _x^{(1)\infty }-\sigma _x^{(1){71}}|, |\sigma _\Pi ^{(1)\infty }-\sigma _\Pi ^{(1){71}}|\}}{\sigma _{x}^{(1){71}}}\\&\quad \approx {\left\{ \begin{array}{ll} \{0.06, 0.20, 0.15\} &{} [\texttt{S0}]\\ \{0.15 , 0.20, 0.15 \} &{} [\texttt{S1}] \end{array}\right. } \, . \end{aligned} \end{aligned}$$

(24)

The first relation simply states that the impact of the form factor at tree level is clearly dominating any TPE effect, as expected. Comparing the other relations in (24) we note that the form factor insertion turns out to be more relevant for the mixed and fermionic corrections and less relevant for the electronic OPE correction. With the same normalisation, the impact of varying $\Lambda$ at LO and NLO can be quantified as

$$\begin{aligned} \begin{aligned}&\frac{|\sigma _0^{86}-\sigma _0^{60}|}{\sigma _{x}^{(1){71}}} \approx 6 \,, \\&\quad \frac{\{|\sigma _e^{(1)86}-\sigma _e^{(1){60}}|, |\sigma _x^{(1)86}-\sigma _x^{(1){60}}|, |\sigma _\Pi ^{(1)86}-\sigma _\Pi ^{(1){60}}|\}}{\sigma _{x}^{(1){71}}}\\&\quad \approx {\left\{ \begin{array}{ll} \{0.02, 0.11, 0.10\} &{} [\texttt{S0}]\\ \{0.08, 0.11, 0.10\} &{} [\texttt{S1}] \end{array}\right. } \, . \end{aligned} \end{aligned}$$

(25)

Not surprisingly, variations in $\Lambda$ have the largest impact on the LO result. At NLO, the effect on the electronic, mixed, and fermionic corrections is roughly the same, and still more than half the size of the effect of the form factor inclusion itself, shown in (24). Thus, in the context of ep scattering, a calculation of NLO corrections necessarily requires a precise inclusion of the proton form factors.

Comparing the NNLO pure QED corrections (involving also pointlike three-photon exchanges) to the TPE effects we find

$$\begin{aligned} \begin{aligned}&\frac{\{ |\sigma _e^{(2)\infty }|, \sigma _{x}^{(2)\infty }, \sigma _{e\Pi }^{(2)\infty }, \sigma _{x\Pi }^{(2)\infty } \}}{\sigma _{x}^{(1){71}}}\\&\quad \approx {\left\{ \begin{array}{ll} \{0.08, 0.17, 0.34, 0.03\} &{} [\texttt{S0}]\\ \{0.06, 0.03, 0.06, 0.03\} &{} [\texttt{S1}] \end{array}\right. } \,. \end{aligned} \end{aligned}$$

(26)

Their relative size is of the same order (if not bigger) as the impact of adding form factors to the NLO TPE corrections, (24), or considering uncertainties of the TPE implementation, (25). This is particularly the case if forward energetic photons are not restricted. Hence, a detailed effort to improve the description of TPE contributions needs to be combined with NNLO QED corrections.

The case of $\mu p$ scattering behaves differently in some respects. The difference between muon and proton mass is much smaller than the difference between electron and proton mass. This is why we do not observe any collinear enhancements, and the NLO muonic corrections are smaller than the electronic corrections in (23),

$$\begin{aligned}{} & {} \frac{\{|\sigma _\mu ^{(1)71}|, \sigma _\Pi ^{(1)71}, |\sigma _p^{(1)}|\}}{\sigma _{x}^{(1){71}}}\approx \{0.4, 4.2,0.02\} . \quad \end{aligned}$$

(27)

Again, the protonic corrections are small enough to justify their approximate (pointlike proton) inclusion. Here and in the following, we only discuss the scenario S0, because the cut on energetic (forward) photons has a marginal effect only.

Considering the same quantities as for $e^- p$ scattering, (24), but now for $\mu ^- p$ scattering, again normalising by the TPE corrections $\sigma _{x}^{(1){71}}$, we find

$$\begin{aligned}&\frac{|\sigma _0^{\infty }-\sigma _0^{71}|}{\sigma _{x}^{(1){71}}} \approx 15 \,, \nonumber \\ {}&\quad \frac{\{|\sigma _\mu ^{(1)\infty }-\sigma _\mu ^{(1){71}}|, |\sigma _x^{(1)\infty }-\sigma _x^{(1){71}}|, |\sigma _\Pi ^{(1)\infty }-\sigma _\Pi ^{(1){71}}|\}}{\sigma _{x}^{(1){71}}}\nonumber \\&\quad \approx \{0.05, 0.30, 0.20\} \, . \end{aligned}$$

(28)

Thus, the impact of adding the form factor is slightly larger in the $\mu p$ case. Proceeding in analogy with the ep case, we next consider the impact of varying $\Lambda$ for $\mu p$

$$\begin{aligned}{} & {} \frac{|\sigma _0^{86}-\sigma _0^{60}|}{\sigma _{x}^{(1){71}}} \approx 8 , \quad \nonumber \\ {}{} & {} \frac{\{|\sigma _\mu ^{(1){86}}-\sigma _\mu ^{(1){60}}|, |\sigma _{x}^{(1){86}}-\sigma _{x}^{(1){60}}|, |\sigma _{\Pi }^{(1){86}}-\sigma _{\Pi }^{(1){60}}| \} }{\sigma _{x}^{(1){71}}} \nonumber \\{} & {} \quad \approx \{0.02,0.13,0.12\} \end{aligned}$$

(29)

and note that the results are similar as in (25). Finally, comparing TPE to pure NNLO QED, we get

$$\begin{aligned} \begin{aligned}&\frac{\{ |\sigma _\mu ^{(2)\infty }|, \sigma _{x}^{(2)\infty }, \sigma _{\mu \Pi }^{(2)\infty }, \sigma _{x\Pi }^{(2)\infty } \}}{\sigma _{x}^{(1){71}}} \\ {}&\quad \approx \{0.00, 0.02, 0.05, 0.03\} \,. \end{aligned} \end{aligned}$$

(30)

Here we see a clear difference between (26) and (30). As expected, higher-order QED radiative corrections are less relevant for $\mu p$ scattering. The impact of variation in the TPE evaluation is larger than the pure NNLO QED corrections. From this perspective, it is thus advantageous to study TPE in $\mu p$ scattering. However, the pure NNLO QED corrections add up to 10% of the TPE. Therefore, a precision study still benefits from inclusion of state-of-the-art QED corrections.

Fig. 6

Differential cross section for ep scattering w.r.t. with respect to the electron scattering angle, for S0 (left panel) and S1 (right panel). The cross section is split into different contributions at LO (gray), NLO (yellow and green) and NNLO (dashed cyan, dotted dark cyan and dark blue). The impact of the proton finite size on the OPE at LO is shown in red, and the impact on the elastic TPE at NLO is shown in pink. Some contributions are presented with their absolute value as the scale is logarithmic

Fig. 7

Same as Fig. 6 but for $\mu p$ scattering

Fig. 8

Difference between the $\theta _e$ differential cross sections for $e^-p$ and $e^+p$ scattering, for S0 (left panel) and S1 (right panel). The corrections to the cross section are split into different contributions at NLO (blue) and NNLO (yellow). The only non-zero corrections are those with odd powers of the formal charge $q_\ell$

Fig. 9

Same as Fig. 8 but for $\mu p$ scattering

We complement our discussion with results at differential level, considering differential distributions w.r.t. the lepton scattering angle. Figures 6 and 7 present such distributions for both kinematical scenarios, for ep and $\mu p$ scattering, respectively. In each plot, pure NLO leptonic and fermionic (yellow curve) as well as mixed corrections (green curve) are compared to the difference of LO effects with and without inclusion of the proton form factors (red curve). Furthermore, pure NNLO QED corrections (blue curves) are compared to the difference of NLO corrections with and without inclusion of the proton form factors (pink curve). Thus, the impact of the proton form factor inclusion at LO and NLO can be contrasted at differential level to NLO and NNLO pointlike corrections, respectively. The statements previously made for the integrated cross section are confirmed at the differential level. In the case of ep scattering, NLO pure QED corrections outweigh the effect of the form factor inclusion at LO, and NNLO pure QED corrections are comparable to the impact of the form factor inclusion in the NLO TPE correction, particularly for S0. The impact of pure QED corrections on $\mu p$ scattering, which can be read from the same set of curves, is smaller than in the ep case but still not negligible.

The cross-section difference between $\ell ^-p$ and $\ell ^+p$ allows to cancel radiative contributions with an even power of the formal charge $q_\ell$. Thus, up to and including NNLO, this leaves only the $\sigma _x^{(1)}$ and $\sigma _{x\Pi }^{(2)}$ contributions, and a subset of the $\sigma _{x}^{(2)}$ contribution. Figures 8 and 9 show the latter in both kinematical scenarios for ep and $\mu p$ scattering, respectively. One can see again that, in the case of S1 for ep scattering and in general for $\mu p$ scattering, NNLO corrections are more suppressed than in S0 for ep. Nevertheless, when extracting the TPE effect empirically from the cross-section difference measured at MUSE, i.e. in scenario S1 for ep or $\mu p$ scattering, it is important to take higher-order radiative corrections into account. In both cases, NNLO contributions will lead to a $10\%$ correction on the extraction, thus, cannot be neglected.

4 Conclusions and outlook

We have presented an update of the McMule framework for the process of lepton-proton scattering [10] with inclusion of additional proton-structure effects from elastic TPE, and complete pointlike QED corrections at NNLO, with lepton-mass effects. In Sect. 2, we have given a detailed description of the contributions that are included in the latest version of McMule. Our notation for individual contributions has been introduced in (18)–(20). In Sect. 3, we have studied the impact of higher-order QED radiative corrections on the unpolarised cross section for lepton-proton scattering at MUSE, focusing on one particular choice of beam momentum ($p_\text {beam}=210$ MeV). The availability of both electrons and muons, with both polarities, is a remarkable advantage for the MUSE experiment, as it allows to analyse a diversified phenomenology and to keep under control QED radiative corrections, if needed. This is achieved either with physical cuts on hard forward photons or by using muons, which are less inclined to irradiate.

Hadronic corrections are usually known less precisely than pure QED corrections (with a pointlike proton). In this work, our main aim has been to assess the relative size of NNLO pure QED corrections, as compared to the LO and NLO corrections with inclusion of the proton form factors and their uncertainties. A particular focus has been on TPE effects, referred to also as the NLO mixed corrections. Since the MUSE kinematics is limited to the low momentum-transfer region ($Q^2<0.08$ GeV$^2$), the inelastic TPE is small enough to be neglected in view of the anticipated $1\%$ accuracy of the cross section measurement. Therefore, only the elastic TPE has been implemented through a simple hadronic model assuming on-shell proton form factors described by a dipole ansatz (10). The dipole parameter (10c) has been varied around the standard dipole $\Lambda ^2=0.71$ GeV$^2$ within a broad range $0.60\,\textrm{GeV}^2< \Lambda ^2 < 0.86\,\textrm{GeV}^2$ to illustrate the impact of form factors uncertainties.

We conclude that while it is sufficient to evaluate the protonic NLO corrections in pure QED with a pointlike proton, all other NLO corrections, in particular the mixed and fermionic, necessarily require a precise inclusion of the proton form factors. Furthermore, we haven shown that NNLO pure QED corrections can be almost as sizeable as the NLO TPE corrections. Even for ep scattering with cuts on hard forward photon emission (scenario S1), see (26), or for $\mu p$ scattering, see (30), where higher-order radiative corrections are more suppressed, NNLO QED corrections should always be included together with an improved description of TPE effects. Equivalently, NNLO pure QED corrections need to be included when extracting the TPE effect empirically to better than $10\%$ accuracy from the cross-section difference between $\ell ^-p$ and $\ell ^+p$ scattering.

The same analysis can be readily repeated within the McMule framework for different kinematical scenarios or other observables, and also broadened to cover further experiments with different $Q^2$ ranges. To this end, the implementation of the elastic TPE correction with input from modern form-factor parametrisations, and the implementation of inelastic TPE corrections, are planned for a future version of McMule.