Rare decay $\pi^0\to e^+e^-$: on corrections beyond the leading order

The preceding experimental and theoretical results on the rare decay $\pi^0 \to e^+e^-$ are briefly summarized. Already computed two-loop QED corrections are reviewed and the bremsstrahlung contribution beyond the soft-photon approximation is analytically calculated. The possible further contribution of QCD loop corrections is estimated using the leading logarithm approximation. The complete result can be used to fit the value of the contact interaction coupling $\chi^r$ to the recent KTeV experiment with the result $\chi^r(M_\rho) = 4.5 \pm 1.0$.


Motivation
Experimental measurements of the rare decay of a neutral pseudoscalar meson to a lepton pair and its comparison with theoretical predictions offer an interesting way to study low-energy (long-distance) dynamics in the Standard Model (SM) [1,2,3]. Systematical theoretical treatment of the process dates back to 1959, when the first prediction of the decay rate was published by Drell [4]. While the possible contributions of the weak sector of the SM are small enough to be neglected, the leading order QED contribution is described by two virtual photon exchange triangle diagram. That is why the double off-shell pion transition form factor F π 0 γ * γ * , which is not known from the first principles, plays essential role.
Because of this one-loop structure for the leading order, the process is very rare and suppressed in the comparison to two photons decay (π 0 → γγ) by a factor of 2(αm e /M π 0 ) 2 due to the approximate helicity conservation of the interaction and thus may be sensitive to possible effects of the physics beyond the SM (expected branching ratio from the pure SM calculation is about 10 −7 ).
Recently, this decay has attracted attention of the theorists again in connection with a new precise branching ratio measurement. The KTeV-E799-II experiment at Fermilab [5] has observed π 0 → e + e − events (altogether 794 candidates), where K L → 3π 0 decay was used as a source of (1) Here we have introduced the Dalitz variable where p, q and k are four-momenta of electron, positron and photon, respectively, P = (p + q + k) is the fourmomentum of neutral pion π 0 with a mass M and E k is the energy of the real outgoing photon in the pion CMS. The lower bound of the Dalitz variable x is used to suppress the contribution of the Dalitz decay π 0 → e + e − γ, which naturally arises with lower x. By means of extrapolating the Dalitz branching ratio in (1) to the full range of x, the branching ratio of the neutral pion decay into an electron-positron pair was determined to be equal to Here the first error is from data statistics alone and the second is the total systematic error. For the matter of interest, current PDG average value (6.46±0.33)×10 −8 [6] is mainly based on this new result. The KTeV collaboration used the result (3) for further calculations. They used the early calculation of Bergström [7] to extrapolate the full radiative tail beyond x > 0.95 and to scale the result back up by the overall radiative corrections of 3.4 % to get the lowest order rate (with the final state radiation removed) for π 0 → e + e − process. The final result is B no-rad KTeV (π 0 → e + e − ) = (7.48 ± 0.29 ± 0.25) × 10 −8 . (4) Subsequent comparison with theoretical predictions of the SM was made in [1,2] using pion transition form factor data from CELLO [8] and CLEO [9] experiments. Finally, it has been found, that according to SM the result should be This can be interpreted as a 3.3 σ discrepancy between the theory and the experiment. Of course, the discrepancy initiated further theoretical investigation of its possible sources [10,11]. Aside from the attempts to find the corresponding mechanism within the physics beyond the SM, also the possible revision of the SM predictions has been taken into account. Many corrections of this kind have been already made, but so far with no such a significant influence on the final result.

Leading order
According to the Lorentz symmetry the on-shell invariant matrix element of the π 0 → e + e − process can be generally written in terms of just one pseudoscalar form factor and, as a consequence, the total decay rate is given by where m stands for electron mass and ν ≡ 2m/M . The leading order in the QED expansion is depicted as the left hand side of the graphical equation in the Fig. 1. Here the shaded blob corresponds to the off-shell pion transition form factor F π 0 γ * γ * (l 2 , (P − l) 2 ) where l is the loop momentum. This form factor serves as an effective UV cut-off due to its 1/l 2 asymptotics governed by OPE (see e.g. [12]) and the loop integral over d 4 l is therefore convergent. It is convenient to pick up explicitly the non-analytic contribution of the two-photon intermediate state (the imaginary part 1 is determined uniquely up to the normalization given by the on-shell value of F π 0 γ * γ * (0, 0) ≡ F π 0 γγ ) and express the form factor in the following way (cf. [13]) 1 Imaginary part of this contribution is given by Cutkosky rules cutting the two virtual photon lines in the Fig. 1. Figure 1: Leading order contribution in the QED expansion and its representation in terms of the leading order of the chiral perturbation theory.
Here, Li 2 is the dilogarithm, and µ represents the intrinsic scale connected with the form factor 2 F π 0 γ * γ * . The function χ P 2 /µ 2 , m 2 /µ 2 represents the remainder which collects the contributions of higher intermediate states and is real and analytic 3 for P 2 /µ 2 < 1.
The leading order terms in the chiral expansion of the form factor P LO are depicted as the right hand side of the graphical equation in Fig. 1. The π 0 γγ vertex in the loop graph is local and corresponds to the leading order term of the chiral expansion of the form factor F π 0 γ * γ * . Therefore the loop integration is no more UV finite and a counterterm (represented by the tree graph in the Fig. 1) is necessary. The sum of these two terms can be written in the form (8), where the transition form factor F π 0 γγ and the remainder χ P 2 /µ 2 , m 2 /µ 2 are replaced by their leading orders in the chiral expansion where χ (r) (µ) is the finite part of the above mentioned counterterm renormalized at scale µ. The graphical equation in the Fig. 1 can be understood as the matching condition for χ (r) (µ) at the leading order in the chiral expansion. It enables to determine χ (r) (µ) once the form factor F π 0 γ * γ * is known. The latter can be theoretically modeled e.g. by the lowest meson dominance (LMD) approximation to the large-N C spectrum of vector meson resonances yielding [13] where M ρ = 770 MeV is the mass of the ρ meson. For other alternative estimates cf. Tab. 1 and for the complete discussion see [1]. Using the value (11) we get for the π 0 → e + e − branching ratio numerically Model CLEO+OPE QCDsr LMD+V NχQM Table 1: Numerical values of χ (r) in different models according to [1,3]. The first two columns denoted as CLEO+OPE and QCDsr correspond to various treatments of CLEO data. LMD+V is an improvement of the LMD ansatz and NχQM stands for the nonlocal chiral quark model.

Two-loop virtual radiative corrections
The full two-loop virtual radiative (pure QED) corrections of order O(α 3 p 2 ) were calculated in [3]. In this section we will present a short review of the main results.
The relevant contributions to the amplitude are shown in Fig. 2. There are six two-loop diagrams. Listed sequentially, we have two vertex corrections (a, b), electron self-energy insertion (c), box-type correction (d) and two vacuum polarization insertions (e, f). Of course, for every such diagram a one-loop graph with corresponding counterterm must be added to renormalize the subdivergences. The relevant finite parts of these counterterms can be fixed by the requirement that the parameters m and α coincide with their physical values. After the subdivergences are canceled, the remaining superficial divergences has to be renormalized by another additional tree counter-term with coupling ξ. The finite part ξ (r) (µ) of this coupling has been estimated in [3] using its running with the renormalization scale as Besides the UV divergences, the graph (d) in the Fig. 2 is also IR divergent. It is therefore necessary to consider IR-safe decay width of the inclusive process π 0 → e + e − (γ) with additional real photon in the final state. In [3] the real photon bremsstrahlung has been taken into account using the soft-photon approximation. The final result depends on the experimental upper bound on the soft photon energy which can be expressed in terms of the lower bound x cut on the Dalitz variable x (see (2)). The result can be expressed in terms of the correction factor δ(x cut ) defined as where Γ LO is the leading order width and Γ NLO is the next-to leading O(α 3 p 2 ) correction. The x cut dependent overall correction δ(x cut ) has various sources and to emphasize the origin of its constituents, we will use the same symbol decorated with appropriate indices. For the complete QED two-loop correction δ (2) including soft-photon bremsstrahlung and KTeV cut x cut =0.95, in [3] it was obtained where only the uncertainties of χ (r) and ξ (r) were taken as the source of the error. This result differs significantly from the previous approximate calculations done by Bergström [7] or Dorokhov et al. [10], where for δ (2) (0.95) we would get −13.8 % and −13.3 %, respectively. There is a simple interrelation of this partial result of the QED radiative corrections and the branching ratio (3) obtained by KTeV experiment (for the details see [3]). We can write the theoretical prediction for the branching ratio measured by KTeV as where the only experimental input is the precise branching ratio B(π 0 → γγ) = (98.823 ± 0.034) %. In the above formula, corresponds to the unsubtracted fraction of the Dalitz de-cay background 4 omitted in the KTeV analysis and discussed in [14,3]. In what follows we will concentrate on the last missing ingredient of the formula (16), namely which is the difference between the exact bremsstrahlung and its soft photon approximation. This difference has been only roughly estimated in [3] and this estimate has been taken as a source of the error. Our aim is to calculate ∆ BS exactly and test the adequacy of the soft photon approximation for the cut x cut = 0.95 used in the KTeV analysis.

Bremsstrahlung
In this section, we discuss the above mentioned exact bremsstrahlung (BS), i.e. the real radiative correction corresponding to the process π 0 → e + e − (γ) beyond the softphoton approximation. As a consequence of the gauge invariance, the invariant amplitude for the BS correction (where k and ε * ρ (λ) (k) is the photon momentum and polarization vector, respectively) has to satisfy the Ward identity for on-shell k and thus it can be generally expressed in the form [14] iM BS ρ (p, q, k) = in terms of scalar form factors P , A and T . These are functions of two independent kinematic variables (x, y), defined as As mentioned above, x is the Dalitz variable (i.e. a normalized square of the total energy of e + e − pair in their CMS) and y has the meaning of a rescaled cosine of the angle included by the directions of outgoing photon and positron in the e + e − CMS. The modulus squared of the amplitude has the form [14] |M and using the variables x, y the differential decay rate is To the amplitude M (λ) (p, q, k) five Feynman diagrams contribute (cf. Fig. 3). Four of them correspond to the photon emission from the outgoing fermion lines (see Fig. 3a-3d). Naively, one would expect that only these four diagrams are necessary to consider since only they include IR divergences which are needed to cancel the IR divergences stemming from the virtual corrections (see graph (d) in the Fig. 2 and the corresponding one-loop diagram with counterterm). However, this result would not be complete. The reason is that the Ward identity (20) would be violated 5 . Thus it is necessary to add the third (box) diagram (Fig. 3e, photon emitted from the inner fermion line) to fulfill this relation. In the graphs (3a) and (3b) the πγγ vertex stems from the Wess-Zumino-Witten action [15,16] and the remaining vertices correspond to standard QED Feynman rules. These graphs are UV divergent by power counting and have to be regularized. In what follows, we use the dimensional regularization. In order to bypass the problems with intrinsically four-dimensional objects like γ 5 and the Levi-Civita pseudo-tensor ε µναβ , we use its variant known as Dimensional Reduction 6 (cf. [17]), which keeps the algebra of γ-matrices four-dimensional while the loop tensor integrals are regularized dimensionally and expressed in terms of the scalar one-loop integrals using the Passarino-Veltman reduction [18]. Within this framework we first get rid of the Levi-Civita tensor using the four-dimensional identities, e.g.
and then contract the reduced tensor integrals with the γmatrix structures 7 . The contributions of the box diagram Fig. 3e turns out to be finite while the triangle diagrams Fig. 3a and Fig. 3b contain subdivergences which have to be renormalized by means of the tree graphs with counterterm corresponding to the coupling χ (see Fig. 3c and 3d). Summing all the relevant contributions and using the fourdimensional Dirac algebra, we get finally the form factors P , A, and T , the explicit form of which is summarized in Appendix A. The differential decay rate dΓ BS (x, y) (cf. (24)) give rise to IR divergences when integrated over the phase space.  The divergences originate from the soft-photon region which is defined in terms of the variables (x, y) by means of the cut on the Dalitz variable x > x cut . These divergences are exactly the same as those stemming from an analogous integral of the differential decay rate dΓ BS soft (x, y) calculated within the soft-photon approximation. The latter is already included in the two-loop result [3], we therefore present our result for the exact BS as a difference the integral of which is IR finite. The result for dΓ BS diff (x, y) is shown in Fig. 4 and (integrated over the allowed region of y given by (22)) in Fig. 5. For ∆ BS (x cut ) we get finally The dependence of ∆ BS (x cut ) on x cut is shown in Fig. 6. For x cut = 0.95 and for χ (r) given by (11) we get numerically where the error stems from the uncertainty in χ (r) (M ρ ). In other words, using this cut of Dalitz variable in KTeV experiment, the soft-photon approximation is a very good approach to the exact result. The dependence of ∆ BS (0.95) on χ (r) is shown in Fig. 7. Now we have all ingredients needed in formula (16) under control and we can thus fit the value of the coupling χ (r) to meet the experiment with the result The error is dominated by the experimental uncertainty, while the theoretical error corresponding to the estimate (13) is negligible. To compare, some previously estimated values, which were considered as relevant, are shown in Tab. 1.

Estimate of the theoretical uncertainty of χ (r)
The above determination of χ (r) represents an effective LO value of this coupling and includes therefore implicitly higher order chiral contributions. The corrections to the LO value of χ (r) start at the NLO and stem from the two-loop graphs which correspond to a substitution of the one-loop subgraphs (and corresponding counterterms) for the shaded blob on the left hand side of the graphical equation depicted in Fig. 1. The relative size of such corrections is set by the factor (M/4πF ) 2 ∼ 10 −2 and can be naively treated as negligible, however it can be significantly numerically enhanced by the large double logarithm terms like log 2 (µ 2 /m 2 ) ∼ 10 2 for µ ∼ M ρ . Complete calculation of the NLO corrections is beyond the scope of the present article. In this section, we will only restrict ourselves to the rough estimate based on explicit calculation of the above mentioned leading (double) logarithms, which are expected to represent a numerically relevant part of the full NLO contribution. According to the Weinberg consistency relation [19], this can be achieved by means of evaluation of infinite parts of one-loop graphs only. In what follows, we will adapt this relation to our case.
Let us write the contribution of the above mentioned two-loop graphs as where the first three terms correspond to one-particle irreducible (1PI) contributions (including two-loop graphs, one-loop graphs with counterterms and tree counterterm graphs) and the last term represents the renormalization of the external pion line by means of the one-loop Z-factor. The contributions of the 1PI loop graphs P 2-loop can be written schematically 8 as an expansion in ε = 2 − d 2 In the same way, for P 1-loop CT we get (see Fig. 8) and the one-loop ingredients of the term Z 1-loop 1/2 P LO are then in the same way (see Fig. 9) Here l (r) (µ) and χ (r) (µ) are finite parts of the one-loop counterterms. We use the standard notation for the two-flavour Chiral Perturbation theory (ChPT) both in the even [20,21] and in the odd sector [22]. The coefficient β χ can be obtained from (8) and (10) and β 4 will be discussed below. The Weinberg condition is based on absence of nonlocal divergences of the form log µ 2 /ε. It can be expressed as the following constraint (36) 8 Because we are interested only in the singular parts we ignore the difference between M S, M S and M Sχ subtraction schemes in what follows. Such an omission can affect only the finite parts which are irrelevant for the leading log calculation.
The contribution of the leading double logs P LL is (37) Using the constraint (36), we get finally Let us now discuss the ingredients of the formula (38). The infinite parts of the couplings χ and l 4 are From the finiteness of P LO , it follows For the couplings c W i , the infinite parts depend on the form of the l 4 term in the chiral Lagrangian (see (see [23,24,25]) for details). For the standard choice we get (in this case, β 4 = 0), while for equivalent case, which differs by terms proportional to the LO equation of motion we get β 4 = −(M/F ) 2 and Because both choices have to lead to the same result, we get the following relation The Z factor is not a physical observable therefore is both sensitive to the field redefinition and in principle infinite. To calculate it we will use the exponential parametrization U = exp (iφ/F ) (see e.g. [24]): 7 The only missing ingredients are then P 1-loop 11,−1 , P 1-loop 13,−1 and P 1-loop χ,−1 , which correspond to the one-loop graphs depicted in the Fig. 8. Explicitly, we get Putting all these ingredients together, we find that i=7,11 (50) and we get finally which implies the following leading log correction, which has to be subtracted from the experimentally determined coupling (30) Numerically which is well below the uncertainty of χ (r) in (30). This can be taken as an indication of the robustness of our determination of χ (r) with respect to the NLO chiral corrections.

Conclusion
In this article we have revisited the decay π 0 → e + e − . It has attracted a lot of attention since its recent precise measurement by KTeV Collaboration at Fermilab due to the discrepancy with the theoretical predictions. Provided that the measured quantity is in agreement with the future experiments one can attribute the existing discrepancy to the quantum corrections, correct modeling of the double off-shell pion transition form factor F π 0 γ * γ * and/or possible contribution of the new physics. Our focus here was on the first part, i.e. standard model corrections to the leading order calculation. We have first briefly summarized recent precise theoretical works dealing with the two-loop QED corrections. The missing bremsstrahlung contribution to this process has been calculated. We have shown that the soft-photon approximation is an adequate approach in the region of KTeV experiment. Besides the electromagnetic corrections we have also studied possible stability in the strong sector. It is best modeled using the higher pion-loop contributions for example in the framework of SU (2) ChPT. It is often the case that in the twoflavour ChPT the order of these corrections can be estimated by the size of the chiral logarithms. In fact they represent the potential enhancement of the usual counting. We have explicitly calculated the coefficient of the leading logarithm and due to the large suppression factor 1/72 (see (51)) it turns out to be very small. This might be an indication of the fast convergence of the perturbation series which is a situation similar to the chiral corrections of π 0 → γ ( * ) γ ( * ) decay (cf. [25,26]).
Using the most reliable QCD modeling of the F π 0 γ * γ * via the lowest-meson dominance approach [13] we agree with the estimate made in [3] of 2σ discrepancy between the theory (including all radiative corrections) and the experiment. Let us remind that this number is significantly smaller than usually quoted difference (3.3σ), however, let us stress that this bigger number was obtained from the rough estimates of the QED radiative corrections and it is thus an indication of the importance of the full two-loop calculation for this process.
On the other hand, still unsatisfactory situation in the first-principle modeling of the three-point vector-vectorpseudoscalar correlator leads to the possibility to use the precise measurement and the full radiative calculation of this process to set the hadronic form factor, represented for this process by the constant χ. The obtained value χ (r) (M ρ ) = 4.5±1.0 (see (30)) is slightly different from the usual estimations, however, represents the model independent prediction for this quantity, based on the KTeV experiment. It can be further used e.g. in the hadronic lightby-light contribution of the muon g − 2 (see e.g. [27,28] for details).