Status of the BabaYaga event generator

The Monte Carlo event generator BabaYaga has been developed for high precision simulation of QED processes (e+e− → e+e−, e+e− → μ+μ− and e+e− → γγ) at flavour factories, chiefly for luminometry purposes, with an estimated theoretical accuracy at the 0.1% level or better. QED radiative corrections are included by means of a Parton Shower in QED matched with exact next-to-leading order corrections to reach the required accuracy. The latter is assessed by means of consistent comparisons to independent calculations and an estimate of the size of missing higher-order corrections. The main theoretical framework is overviewed and the status of the generator is summarized.


Introduction
The knowledge of the luminosity L is an important ingredient for any measurement at e + e − machines. The common strategy is to calculate it through the relation L = N obs /σ th , where σ th is the theoretical cross section of a QED process, namely e + e − → e + e − (Bhabha), e + e − → µ + µ − or e + e − → γγ, and N obs is the number of observed events. QED processes are the best choice because of their clean signal, low background and the possibility to push the theoretical accuracy up to the 0.1% level or better. The latter requires the inclusion of the relevant radiative corrections (RCs) in the cross sections calculation and their implementation into Monte Carlo (MC) event generators (EGs) in order to easily account for realistic event selection criteria.
Modern EGs used for luminometry simulate QED processes by including the exact next-to-leading order (NLO) QED corrections and/or a leading-logarithmic (LL) approximation of higher-order (h.o.) effects [1][2][3][4][5][6]. The consistent inclusion and matching of NLO and h.o. LL contributions is mandatory in view of the required theoretical accuracy. In the following, it is discussed how this is achieved and implemented in the BabaYaga EG.
In its first version [1], the generator relied upon a QED Parton Shower (PS) to account for the LL photonic corrections, resummed up to all orders in perturbation theory. The PS is a MC algorithm which gives an exact iterative solution of the Dokshitzer-Gribov-Lipatov-Altarelli-Parisi (DGLAP) equation in QED for the non-singlet QED structure function (SF) D(x, Q 2 ), which reads where P + (x) is the regularized Altarelli-Parisi vertex, x is the fraction of energy lost because of radiation and Q 2 is the scale of the process. The QED SFs account for photon radiation emitted by both initial-state and final-state fermions and allow to include the universal virtual and real photonic corrections, resummed up to all orders of perturbation theory. The advantage of the PS solution is that the kinematics of the emitted photons can be recovered (within some approximation) and hence an exclusive event generation can be performed, i.e. all the momenta of the final state particles (fermions and an indefinite number of photons) can be reconstructed. Within the SF approach the corrected cross section can be written as Despite its advantages, the PS is intrinsically accurate at the LL level and a precision better than 0.5-1% can not be expected in the calculation of the cross section (2). In order to improve the accuracy, a matching with the exact NLO RCs is mandatory, in such a way that the features of the PS are preserved (i.e. exclusive event generation and resummation of LL corrections up to all orders) while avoiding the double counting of the O(α) LL corrections, present both in the PS approach and in the NLO calculation.
An original matching algorithm has been implemented in the latest version of BabaYaga [2] (BabaYaga@NLO), which includes exact NLO corrections into a PS framework and achieves an accuracy at the 0.1% level in the calculation of the cross sections.
Without spelling out the details, the fully differential PS cross section implicit in Eq. (2) can be recast in the form where Π(Q 2 , ε) is the Sudakov form factor accounting for virtual and soft (up to x = ε) radiation, M n,PS is the amplitude for the emission of n real photons in the PS approximation and dΦ n is the exact phase space for the emission of n real photons (with x ≥ ε), divided by the incoming flux factor. Eq. (3) can be improved to include the missing NLO contributions by defining the following correction factors where S V stands for soft and virtual photon corrections, H for non-soft real photon corrections, [α, PS ] stands for the O(α) expansion of the PS contribution, i runs over the emitted photons and dσ 0 is the lowest-order differential cross section. With these definitions, the matched differential cross section can be written in the form which is the master formula according to which event generation and cross section calculation are performed in BabaYaga@NLO.
A few comment are in order with respect to the master formula (5): 1. by construction, the factors in Eq. (4) are infrared and collinear safe quantities and they let the O(α) expansion of dσ ∞ matched in Eq. (5) (4) tend to be larger in those phase space regions where the PS is more unreliable, typically in phase space regions where the photon is hard and/or not collinear to one of the charged particles; 4. Eq. (5) is cast in a completely differential form, so that events can be generated exclusively as in the pure PS approach; 5. the theoretical error is shifted to corrections of order α 2 , i.e. at the next-to-next-to-leading-order (NNLO) level.
Eq. (5) is used in BabaYaga@NLO to generate events for Bhabha, e + e − → µ + µ − and e + e − → γγ [3] processes, with an indefinite number of extra photons to account for NLO and h.o. QED radiative corrections. Some numerical results and the estimate of the theoretical accuracy of the approach are sketched in the next section.

Results and estimate of the theoretical error
Table (1) shows the impact of different classes of QED RCs in the determination of large-angle Bhabha cross section within typical event selection criteria used for luminometry at flavour factories. Setups a-b correspond to √ s = 1.02 GeV and c-d to √ s = 10 GeV, for large (a-c) and narrower (b-d) angular acceptance cuts. In the second paper of Ref. [2] the detailed description of the setups, which tend to select elastic events, is reported. The first row of table (1) is the lowest order cross section, the second includes the effects of vacuum polarization (VP), the third one is the exact NLO result, the fourth is the O(α) expansion of   (5). In order to assess the theoretical accuracy of the approach implemented in BabaYaga@NLO, an important step is performing tuned comparisons to independent MC EGs.
As an example, in Figs. (1) and (2), the acollinearity and final-state invariant mass distributions obtained with BabaYaga@NLO and the independent EG Bhwide [4] are compared for setup a. While the difference on the integrated cross sections is below the 0.01% level, some differences at the 1% level appear (see the insets in the figures) in the differential distributions and are due to the different theoretical frameworks on which the two EGs are based. However, it is worth noticing that such differences appear where the differential cross section drops by several order of magnitude with respect to the elastic region and where the different treatment of hard radiation beyond NLO can have an impact.
A further step to estimate the theoretical accuracy is to compare with existing calculations of the NNLO corrections to Bhabha scattering, which have been published over the years [7] (see also Ref. [8] and references therein). The dominant part of such corrections are already included in the master formula (5), which can be expanded up to NNLO and whose NNLO terms can be unambigu-ously compared to analytical exact results: any difference should be considered as a theoretical error in the formulation of BabaYaga@NLO. A detailed and non-trivial comparison has been carried out in the second paper of Ref. [2], finding that the NNLO corrections not included in Eq. (5) impact the cross section at the level of few units in 10 −4 , when typical event selection criteria are taken into account.
In the Luminosity section of Ref. [8], also other sources of theoretical uncertainties were considered, such as the uncertainty on the hadronic contribution to VP and extra light-pair emission (for which see also Refs. [9]). The main conclusion of that section is that a sound estimate of the present theoretical error on the luminosity determination via Bhabha scattering at flavour factories lies in the region of 0.1% or slightly below, which is enough for the present experimental requirements. Any further improvement would require the inclusion of the full exact NNLO RCs into MC EGs for Bhabha scattering, matched with LL h.o. contributions. This should be considered a feasible, although not trivial, task.
In a similar way, the matching algorithm of Eq. (5) has been also applied to the e + e − → γγ process and a phenomenological study at flavour factories has been presented in Ref. [3]. As an example, in table (2) the impact of different classes of radiative corrections for typical event selection criteria at flavour factories is shown, in analogy with table (1) for Bhabha scattering. In the e + e − → γγ case, the impact of NLO RC is in the 5-10% range while h.o. corrections change the cross sections by 0.1-0.5%: both effects must be accounted for to reach a theoretical accuracy at the 0.1% level. To the best of our knowledge, BabaYaga@NLO is the only MC EG implementing a matching of NLO and h.o. corrections in a PS approach for the γγ final state, being able to reach such a theoretical accuracy.
Before concluding, it is worth mentioning that an extension of BabaYaga has been developed for dark photon searches at low energy via the radiative return method [10] and has been extensively used by the KLOE-2 collaboration to set limits on the dark photon couplings and mass [11].
Furthermore, it is planned to exploit the same theoretical framework described above for the development of a high-precision EG for the QED process µe → µe, which is needed for the measurement of the hadronic contribution to the running of α QED in the space-like region and, in turn, for an independent determination of the leading order hadronic corrections to the muon anomalous magnetic moment. The proposed experiment (µONe) and the underlying ideas are discussed in Ref. [12].

Conclusions
The theoretical framework and the status of the BabaYaga and BabaYaga@NLO event generators has been overviewed and summarized. BabaYaga@NLO is a Monte Carlo event generator developed for high-precision simulation of QED processes (Bhabha scattering, e + e − → µ + µ − and e + e − → γγ) at flavour factories, mainly for luminometry purposes, and it is based on an original algorithm to match exact NLO with higher-order QED radiative corrections in a Parton Shower approach. Within typical event selection criteria for luminometry, the cross sections can be calculated with an overall theoretical error at the level of 0.1% or better. The assessment and the scrutiny of the theoretical error are based on consistent and detailed comparisons with independent generators and existing exact NNLO calculations of the relevant radiative corrections.