Electron-positron annihilation to photons at $O(\alpha^3)$ revisited

We apply the modern multiloop methods to the calculation of the total cross sections of electron-positron annihilation to 2 and 3 photons exactly in $s/m^2$ with the accuracy $O(\alpha^3)$. Examining the asymptotics of our results, we find agreement with Ref. [Andreassi_et_al_1962] and discover mistakes in the results of Refs. [Eidelman&Kuraev_1978,Berends&Kleiss_1981]. This mistake is due to the terms, omitted in differential cross section in Refs. \[Eidelman&Kuraev_1978,Berends&Kleiss_1981], which are peaked in the kinematic region with all three photons being quasi-parallel to the collision axis. After restoring these terms, we find an agreement of the corrected result of Ref. [Berends&Kleiss_1981] with our result.


Introduction
Modern methods of multiloop calculations crucially reduce the efforts required to check and improve the available results on radiative corrections. In this work we use this fortunate circumstance in order to calculate the total cross sections of the processes e + e − → 2γ and e + e − → 3γ with accuracy O(α 3 ) for arbitrary energies. Surprisingly, we find that several results available in the ultrarelativistic limit contain errors. In particular, there is no correct result for the total cross section of e + e − → 3γ in the center-of-mass frame 1 . Our technique is based on the Cutkosky rule which allows one to represent the phase-space integrals via the loop integrals with cut propagators. We apply the differential equations method to calculate the emerging two-loop integrals. We use the dimensional regularization d = 4 − 2 to treat both infrared and ultraviolet divergences.
The paper is organized as follows. In the next section we present our results and discuss important issues related to them. Other sections contain details of the calculation. Th conclusion is presented in the last section.

Results
Let us present our results. Below we use the units = c = m = 1, where m is the electron mass. Since the total cross sections σ e + e − →2γ and σ e + e − →3γ are both infrared divergent, we define σ e + e − →2γ (ω 0 ) and σ e + e − →3γ (ω 0 ) which depend on the soft cut-off ω 0 . The quantity σ e + e − →3γ (ω 0 ) is the cross section of the process e + e − → 3γ integrated over the kinematic region where the energy of any photon is greater than ω 0 . The contribution of the complementary kinematic region (when the energy of one of the three photons is less than ω 0 ) is then added to σ e + e − →2γ to form the finite quantity σ e + e − →2γ (ω 0 ). Note that the restriction of the integration region introduces the dependence of the cross section on the frame, which we denote by the upper superscript f, as in σ f e + e − →3γ (ω 0 ). In the center-of-mass frame we have Here is the Born cross section of the process e + e − → 2γ, and we use the symmetrization symbol It worth noting that βσ cmf e + e − →3γ (ω 0 ) is an analytical function of β 2 (or, equivalently, of s − 4) in the vicinity of β 2 = 0.
The cross section of the process e + e − → 2γ with the account of the first radiative correction has the form The first term here, 1 + π α v σ 0 , is nothing but the Born cross section σ 0 , multiplied by the expansion of the Sommerfeld-Sakharov factor 2πα/v 1−e −2πα/v with v = 2β 1+β 2 being the relative velocity. It is remarkable that, apart from the contribution of term π α v σ 0 , the cross section σ cmf e + e − →2γ (ω 0 ), multiplied by β, is again an analytic function of β 2 in the vicinity of β 2 = 0.
The corresponding cross sections in the rest frame of the electron read where Note that the sum σ e + e − →2γ(γ) = σ e + e − →3γ + σ e + e − →3γ is independent of ω 0 and, hence, of the frame.

Asymptotics
Let us now discuss the asymptotics of the presented results.
Threshold asymptotics.. We start from the threshold asymptotics β 1. The threshold asymptotics of σ cmf,rf e + e − →3γ (ω 0 ) reads The first term in braces is well known and determines the orthopositronium decay width. The threshold asymptotics of σ cmf,rf e + e − →2γ (ω 0 ) reads The first term in braces is known for a long time, see, e.g., Ref. [4]. In particular, this term determines the radiative correction to the parapositronium decay width. Note that the threshold expansion of δσ in Eq. (7) starts from O(β 4 ), so Eqs. (8) and (9) hold both for center-of-mass frame and electron rest frame.
Ultrarelativistic limit.. Let us now discuss the high-energy asymptotics s 1. For the cross section = σ e + e − →3γ we have The asymptotics of σ rf e + e − →3γ (ω 0 ) in the electron rest frame exactly coincides with the corresponding result of Ref. [1]. However, the asymptotics of σ cmf e + e − →3γ (ω 0 ) in the center-of-mass frame does not coincide with the two available results [2,3]. Moreover, these two results differ from each other: σ cmf, Ref. [3] e + e − →3γ (ω 0 ) ≈ We have been able to trace the origin of discrepancy of our result with that of Ref. [3]. Namely, it appeared that Refs. [2,3] have overlooked in the differential cross section the terms that contribute to the total cross section in triply collinear kinematic region, see Appendix.
The ultrarelativistic asymptotics of σ e + e − →2γ reads dotted, and dash-dotted curves correspond to consecutive approximations of the threshold (truncation at β 1 , β 2 , β 3 , left side of the graph) and high-energy (truncation at ( 1 These two asymptotics coincide with the corresponding results of Refs. [3] and [1], respectively. The comparison of the exact cross section with the asymptotic expansions is demonstrated in Fig. 1. Let us present a few terms of high-energy expansion of the cross section σ e + e − →2γ(γ) : It is remarkable that if we diminish by a factor of 2 the term on the last line 2 , we will obtain an extremely good approximation for the exact cross section σ e + e − →2γ(γ) with the largest deviation about 2% taking place at the threshold point.

Calculation of
We start with the calculation of the total Born cross section of the 3-photon annihilation 3 . The diagrams are shown in Fig. 2. We define two LiteRed bases, pdb and xdb, corresponding to the denominators of diagrams iv, v in Fig. 2, respectively. These two bases are sufficient for the IBP reduction of all scalar integrals appearing in the cross section of the process e + e − → 3γ. There are 7 distinct master integrals which we choose as shown in Fig. 3. We use Libra 4 package [5] to reduce the system to -form [6,7]. The  new set of functions, J 1 , . . . , J 7 is related to j 1 , . . . , j 7 via where β = 1 − 4/s and α n = α . . . (α + n − 1) is the Pochhammer symbol. They satisfy the differential system in -form where We fix the boundary conditions at β → 0 (s → 4) by evaluating the following coefficient in the asymptotics: where [j k ] β µ denotes to coefficient in front of β µ in small-β asymptotics of j k and we have explicitly indicated all coefficients which are obvious zeros. Thus, we are left with three nontrivial coefficients, [j 1,4,6 ] β 0 , which are nothing but the naive values of the corresponding integrals at the threshold, j th 1,4,6 = j 1,4,6 (s = 4). Performing the IBP reduction we find that The two remaining integrals j th 1,4 can be calculated exactly in in terms of hypergeometric function, however we choose to follow the same approach as in Ref. [8] when calculating the parapositronium decay width to 4γ. We choose the constant (but -dependent) overall normalization so that j th 1 = 1. Then we have This fully fixes the boundary conditions. Since the total cross section is infrared divergent at = 0, we have to be careful with the overall normalization, namely, we should pay attention to the factors which tend to unity as → 0. We choose the following n-particle phase-space definition in d = 4 − 2 dimensions: where γ E = 0.577 . . . is the Euler constant. The factor e γ E 4π (n−1) conveniently removes γ E and ln 4π in our intermediate formulae 5 . At = 0 the definition turns into the usual definition of phase-space. We also normalize the trace of Dirac matrices by the condition Substituting the results for the master integrals J k , we obtain where z = 1−β 1+β . Note that the cross section contains −1 term, which is due to the infrared divergent contribution of the region where the energy of one of the outgoing photons is small. Thus, in order to obtain the finite quantity, we have to subtract the contribution of this region. We derive the corresponding formulae in Section 4.

Soft-photon contribution
The probability to emit soft photon is usually regulated by the fictitious photon mass. However, within our approach, we must stick to the dimensional regularization. As, to the best of our knowledge the relevant expressions are not in the literature, we derive them here with some details.

Radiation probability.
We start from the following formula 6 for the probability of soft photon radiation: Here k = (ω, k) = (|k|, k) is the photon momentum, n∈i∪f denotes the sum over initial and final particles, e n = q n |e| and p n = m n u n = m n (γ n , γ n β n ) are their charges and momenta (with m n = p 2 n being the mass), and σ n = +1 (σ n = −1) when n ∈ i (n ∈ f ). Note that we have again introduced a factor e γ E 4π for consistency with our previous definitions. The integration over ω will be restricted from above by the infrared cut parameter ω 0 and can be trivially performed. Thus, we have where .

(33)
They satisfy the following differential system in Pfaff form where we have used the notationã = a − x 1 x 2 x 3 /a. The physical region is defined by the inequalities We fix the boundary conditions at the point x 1 = x 2 = x 3 = 1 and travel to the generic point (x 1 , x 2 , x 3 ) in the physical region along the contour γ(0 τ 2) defined piece-wise as The boundary conditions appear to be trivial with the only nonzero constant being We finally obtain where Finally, we obtain where
Let us now derive the cross section of e + e − → 3γ integrated over the kinematic region where the energies of all photons are restricted from below by some experimental cut-off ω 0 . This restriction obviously introduces the frame dependence, and we will specialize our formulae to two physically relevant frames: the center-of-mass frame and the rest frame of the initial electron.
The cross section σ e + e − →3γ (ω i > ω 0 ) is obtained by subtracting from σ e + e − →3γ the contribution of the soft region: where f = cmf and f = rf for the center-of-mass frame and the electron rest frame, respectively. We have The two-photon annihilation Born cross section σ e + e − →2γ should also be calculated with 1 terms retained: We finally arrive at Eq. (1).
Let us now briefly describe the calculation of the virtual correction to the total cross section of e + e − → 2γ. We calculate the contribution of the diagrams depicted in Fig. 4. The IBP reduction of the two-loop diagrams reveals 14 master integrals depicted in Fig. 5. We use Libra to reduce the differential system forj 1−14 to -form. The 'canonical' master integralsJ 1−14 are defined as follows Figure 5: Two-loop (one loop cut) master integrals for the virtual correction to the total cross section of e + e − → 2γ.
They satisfy the differential system We fix the boundary conditions by considering the asymptotic coefficients at s → 4. Most of the nontrivial boundary constants correspond to the naive values of the integrals at the threshold. The only exception is the leading threshold asymptotics of j 3 , proportional to β 1−2 . Namely, we explicitly calculate the following constants: Here we again have chosen the overall factor so thatj th 1 = 1. These boundary conditions, are sufficient to fix the specific solutions for the integralsJ . Using these solutions, we obtain for the "bare" cross section where we have neglected terms suppressed in . Note that the renormalized cross section σ e + e − →2γ still contains −1 terms due to infrared divergence. In order to obtain the observable cross section σ f e + e − →2γ (ω 0 ) we have to add the soft-photon contribution W f (ω 0 )σ 0 , where W f (ω 0 ) is defined in Eq. (41) for f = cmf, rf. Finally, we obtain Eqs. (4) and (6) for the cross sections σ cmf e + e − →2γ (ω 0 ) and σ rf e + e − →2γ (ω 0 ), respectively.

Conclusion
In the present paper we have calculated the total cross sections σ f e + e − →3γ (ω 0 ) and σ f e + e − →2γ (ω 0 ) for arbitrary energies with O(α 3 ) accuracy. The energy cut ω 0 for soft photons has been applied in the centerof-mass frame (f = cmf) and in the rest frame of the electron (f = rf). We have found errors in the high-energy results available in the literature for σ cmf e + e − →3γ (ω 0 ). As an additional check of our results for the 3-photon annihilation cross section, we have performed a numerical integration of the differential cross section using the Cuba library [11]. Table 1 (1), for σ cmf e + e − →3γ with that of Monte-Carlo integration performed using Cuba library [11]. The errors in last column are those provided by Cuba. minimal photon energy ω 0 is of the order of electron mass in the rest frame of the initial electron (marked with wavy lines in the table). This is, of course, quite expected as the photon with energy of the order of electron mass in electron rest frame can have energy of the order of √ s when one passes to the center-of-mass frame (thus, the soft-photon approximation breaks).
expression for s|M | 2 indeed leads to the result (13) of Ref. [3] 7 . So, the origin of the discrepancy of our result with that of Ref. [3] can be only in the initial expression for the differential cross section. Indeed, a thorough inspection of the exact expression for the differential cross section from Ref. [12] has revealed the overlooked in Refs. [2,3] terms which contribute to the total cross section. Namely, in s|M | 2 one has to take into account also the terms These terms contribute in the kinematic region where simultaneously two photons have small scattering angles. Then the third photon necessarily has scattering angle close to π and we have the following power counting: We have checked that these terms, overlooked in Refs. [2,3], give exactly the contribution − 4α 3 π 2 3s to the total cross section, in agreement with our asymptotics (11).