1 Introduction

The Standard Model (SM) assumes lepton universality (LU) implying that once the mass difference is taken into account, all SM interactions treat the three charged leptons identically. The mass difference results in a different phase space between decays involving \(\tau ^-\) and the lighter \(e ^-\) and \(\mu ^-\) leptons.Footnote 1 LU can be tested by measuring the ratio of decay rates, ensuring that the Cabibbo–Kobayashi–Maskawa matrix elements, as well as most of the form factors, cancel in the ratio. This results in more accurate theoretical predictions and in the cancellation of many experimental systematic uncertainties. One type of these LU measurements is performed using semileptonic \(B \) decays of the form \({b} \rightarrow {c} {\ell ^-} {{\overline{\nu }} _\ell } \), commonly known as measurements of \(\mathcal {R}({H_c})\), defined as

$$\begin{aligned} {\mathcal {R}({H_c})} = \frac{{\mathcal {B}} ({H_b} \rightarrow {H_c} {\tau ^-} {{\overline{\nu }} _\tau })}{{\mathcal {B}} ({H_b} \rightarrow {H_c} \ell ^-{{\overline{\nu }} _\ell })} \, , \end{aligned}$$
(1)

where \(H_b\) and \(H_c\) are a \(b \) and \(c \) hadron, respectively, and \(\ell \) is either an electron or muon.

Several measurements of \(\mathcal {R}({H_c})\) have been performed by the LHCb, Belle and BaBar experiments. For \(\mathcal {R}({D})\), on which this paper is focused, the predicted value [1,2,3,4] is

$$\begin{aligned} {\mathcal {R}({D})} = {\mathcal {R}({{D} ^+})} = {\mathcal {R}({{D} ^0})} = 0.299 \pm 0.003, \end{aligned}$$
(2)

which assumes isospin symmetry. The average of the measured value of \(\mathcal {R}({D})\) is \(0.349\pm 0.027\pm 0.015\) [5,6,7], where the first uncertainty is statistical and the second systematic. Even though \(\mathcal {R}({D})\) differs from the SM prediction by only \(1.4\sigma \), it is remarkable that the deviation from the SM of the combined \(\mathcal {R}({D})\) and \(\mathcal {R}({{D} ^*})\) observables is \(3.1\sigma \) [8].

Radiative corrections were long thought to be negligible at the level of precision of measurements and predictions of \(\mathcal {R}({D})\). Recently, however, de Boer et al. [9] presented a new evaluation of the long-distance electromagnetic (QED) contributions to \({{\overline{B}{}} {}^0} \rightarrow {{D} ^+} {\ell ^-} {{\overline{\nu }} _\ell } \) and \({{{B} ^-}} \rightarrow {{D} ^0} {\ell ^-} {{\overline{\nu }} _\ell } \) decays, where \({\ell ^-} = {\mu ^-}, {\tau ^-} \). They point out that these soft-photon corrections are different for \(\mu ^-\) and \(\tau ^-\) decays, such that they do not cancel in the ratios \(\mathcal {R}({{D} ^+})\) and \(\mathcal {R}({{D} ^0})\). According to the authors of Ref. [9], the proper evaluation of the radiative corrections alters the SM predictions of the \(\mathcal {R}({D})\) and \(\mathcal {R}({{D} ^*})\) values and increases their uncertainty. The current tension between the SM and measurements could be weakened or strengthened if radiative corrections are not properly taken into account.

All experiments testing these types of LU are dependent on the simulation of QED radiative corrections in decays of particles and resonances. The widely used package to simulate these corrections is Photos  [10, 11], which is used by all three experiments measuring \(\mathcal {R}({D})\) and \(\mathcal {R}({{D} ^*})\).

This paper starts by comparing the radiative corrections on \(\mathcal {R}({{D} ^+})\) and \(\mathcal {R}({{D} ^0})\) from Ref. [9] with those simulated by Photos in Sect. 2. The sensitivity of measurements of \(\mathcal {R}({{D} ^+})\) and \(\mathcal {R}({{D} ^0})\) to radiative corrections in the \(\mu ^-\) and \(\tau ^-\) decay modes is studied with pseudo-experiments in an LHCb-like environment, with different assumptions on the shape of the total energy of the radiated photons. The method and the results of this study are reported in Sect. 3. Conclusions and recommendations are summarised in Sect. 4.

2 Radiative corrections in PHOTOS

Photos  [11, 12] is a universal Monte Carlo algorithm that simulates the effects of QED corrections in decays of particles and resonances. It exploits the factorisation property of QED coming from the exponentiation method used to improve the convergence of the perturbative expansion. Any particle-decay process accompanied by bremsstrahlung photons can be factorised into a tree term and bremsstrahlung factor. The latter depends only on the four-momenta of those particles taking part in the decay, and not on the underlying process. This approximation, which takes into account both real and virtual corrections, converges to an exact expression in the soft-photon region of phase space. It is worth noting that Photos does not incorporate the emission of photons depending on the hadronic structure. These so called structure-dependent (SD) photons impact the spin of the decay particle, and may also interfere with bremsstrahlung photons. The effect of SD photons depends on the specific decay under study and, as was the case for kaon decays [13], may not be negligible.

The latest versions of Photos include multi-photon emissions, and interference between final-state photons. The validity of Photos has been tested successfully by comparing its results to full calculations available in various processes involving \(W \), \(Z \) and hadronic \(B \) decays into scalar mesons [11, 14]. Because of the universal treatment of photon emission in Photos, its performances in specific processes should always be checked, especially when high precision is desired or when signal extraction is sensitive to detailed simulation of a phase space corner of the decay.

The calculation by de Boer et al. in Ref. [9] is the first that studies the impact of soft-photon corrections on \(\mathcal {R}({{D} ^+})\) and \(\mathcal {R}({{D} ^0})\). It is valid in the regime in which the maximum energy of the radiated photons is smaller than the lepton mass, which is the muon mass in this case. This calculation includes more effects than Photos does, in particular the interference between initial- and final-state photons, and the Coulomb correction. The latter increases the decay rate of decays with charged particles in the final state. It should be noted that the contribution of the Coulomb correction is singular for null relative velocity between final-state charged particles.

To compare QED corrections between Photos and Ref. [9], four samples (\({{\overline{B}{}} {}^0} \rightarrow {{D} ^+} {\ell ^-} {{\overline{\nu }} _\ell } \) and \({{{B} ^-}} \rightarrow {{D} ^0} {\ell ^-} {{\overline{\nu }} _\ell } \), where \({\ell ^-} = {\mu ^-}, {\tau ^-} \)) with three million \(B \)-meson decays are generated by Pythia 8 [15, 16]. The decays are simulated by EvtGen  [17], and the radiative corrections by Photos v.3.56, with the “option with interference” switched on. QED corrections are applied by Photos by modifying the charged track’s four-momentum in the event record filled by EvtGen every time a photon is added.

The four-momentum of the total radiated photons, \(p _{\gamma }\), is defined as

$$\begin{aligned} p _{\gamma } = p _{B} - \left( p _{{D}} + p _{{\ell ^-}} + p _{{{\overline{\nu }} _\ell }} \right) , \end{aligned}$$
(3)

where \(p _{B}\), \(p _{{D}}\), \(p _{{{\overline{\nu }} _\ell }}\), and \(p _{{\ell ^-}}\) are the four-momenta of the \(B \), \(D \), \(\ell ^-\) and \({\overline{\nu }} _\ell \) particles, respectively, taken from the event record updated by Photos. This means that, in agreement with Ref. [9], the radiation of the \(D \) decay products is not taken into account. The total energy of the radiated photons, \(E_{\gamma }\), is computed in the \(B \) rest frame. As in Ref. [9], the variable \(E_{\mathrm{max}}\) is defined as the maximum value that \(E_{\gamma }\) is allowed to have to consider \({B} \rightarrow {D} \ell {\bar{\nu }}_\ell (\gamma )\) as signal.

Fig. 1
figure 1

Distribution of the total energy of the radiated photons, \(E_{\gamma }\), up to 100 MeV for \({\overline{B}{}} {}^0\) \(\rightarrow \) \({D} ^+\) \(\tau ^-\) \({\overline{\nu }} _\tau \) and \({\overline{B}{}} {}^0\) \(\rightarrow \) \({D} ^+\) \(\mu ^-\) \({\overline{\nu }} _\mu \) decays as simulated by Photos. Once a value of \(E_{\mathrm{max}}\) is chosen, e.g. 40 MeV as in the plot, all events with higher \(E_{\gamma }\) values are discarded

The QED correction, \(\delta _{\mathrm{\mathrm QED}}\), is given by the relative variation of the branching ratio when events with total radiated energy greater than \(E_{\mathrm{max}}\) are discarded. This can be calculated as follows:

$$\begin{aligned} {\delta _{\mathrm{\mathrm QED}}} = \frac{\int _{0}^{{E_{\mathrm{max}}}} N(E_{\gamma } ) dE_{\gamma } }{\int _{0}^{\infty } N(E_{\gamma }) dE_{\gamma }} - 1, \end{aligned}$$
(4)

where \(N(E_{\gamma } )\) is the distribution of events with \(E_{\gamma }\). This distribution is shown for \({\overline{B}{}} {}^0\) \(\rightarrow \) \({D} ^+\) \(\tau ^-\) \({\overline{\nu }} _\tau \) and \({\overline{B}{}} {}^0\) \(\rightarrow \) \({D} ^+\) \(\mu ^-\) \({\overline{\nu }} _\mu \) decays in Fig. 1. The considered energy range is up to 100 MeV, which covers the majority of radiative photons, namely 98% of the \(\mu ^-\) decays and 99.7% for the \(\tau ^-\) decays generated by Photos.

Fig. 2
figure 2

Radiative corrections to the branching ratios of \(\bar{B}^0\rightarrow D^+\ell ^-\bar{\nu _\ell }\) (left) and \(B^-\rightarrow D^0\ell ^-\bar{\nu _\ell }\) (middle) decays, as a function of \(E_{\mathrm{max}}\). The long-distance QED corrections to \(\mathcal {R}({{D} ^+})\) (orange) and \(\mathcal {R}({{D} ^0})\) (violet) as a function of \(E_{\mathrm{max}}\) (right). The plots are obtained from simulated data (solid lines, which include the statistical uncertainty) and from Ref. [9] (dashed lines, filled with transparent colours when the uncertainties are significant)

Comparisons between radiative corrections from Photos and Ref. [9] are shown in Fig. 2 for the \({{\overline{B}{}} {}^0} \rightarrow {{D} ^+} {\ell ^-} {{\overline{\nu }} _\ell } \) (left panel) and \({{{B} ^-}} \rightarrow {{D} ^0} {\ell ^-} {{\overline{\nu }} _\ell } \) (middle panel) branching fractions. These plots show differences of up to 2% for \({\overline{B}{}} {}^0\) decays, and 0.5–1% for \({{B} ^-} \) decays. Unfortunately, the effect does not cancel in the ratios of branching fractions. This is clearly visible from the right panel of Fig. 2 where radiative corrections on \(\mathcal {R}({D})\), \(\delta _{QED}(\mathcal {R})\), are shown as a function of \(E_{\mathrm{max}}\). Photos predicts a QED correction that is 0.5% lower than the one in Ref. [9] for \(\mathcal {R}({{D} ^+})\), while it is 0.5% higher than the one in Ref. [9] for \(\mathcal {R}({{D} ^0})\).

2.1 Coulomb correction

A significant part of the radiative corrections in Ref. [9] originates from Coulomb interactions, which are not included in Photos. Note that the Coulomb correction is relevant for the \({D} ^+\) mode, but not for the \({D} ^0\) mode. For a fermion-scalar (and fermion–fermion) pair, this correction is given by

$$\begin{aligned} {\varOmega _{\mathrm{C}}} = \frac{2\pi \alpha }{\beta _{D\ell }} \frac{1}{1-e^{-\frac{2\pi \alpha }{\beta _{D\ell }}}}, \end{aligned}$$
(5)

where \(\alpha = 1/137\) and \(\beta _{D\ell }\) is the relative velocity between the \(D \) meson and the lepton, defined as

$$\begin{aligned} \beta _{D\ell } = \left[ 1 - \frac{4m_D^2 m_\ell ^2}{(s_{D\ell }-m_D^2 - m_\ell ^2)^2} \right] ^{1/2}, \end{aligned}$$
(6)

where \(s_{D\ell } = (p_D+p_\ell )^2\). A well-known approximation of the Coulomb correction by Atwood and Marciano [18], yields \({\varOmega _{\mathrm{C}}} = (1+\pi \alpha ) \, \approx 1.023\) which occurs when \(\beta _{D\ell } \approx 1\). This is accurate for decays with light leptons, but not for those with \(\tau ^-\) leptons. For the semitauonic mode, the typical relative velocity is 0.5–0.9, resulting in a Coulomb correction between 2.5 and 5.0%.

QED corrections from Photos for the \({D} ^+\) mode are also compared with predictions not including the Coulomb correction from Ref. [9]. This reduces the difference of the corrections to the branching ratios between Photos and the theoretical calculations to about 1% and brings the corrections on \(\mathcal {R}({{D} ^+})\) in close agreement, as shown in Fig. 3 (left and middle, respectively).

Fig. 3
figure 3

Radiative corrections to the branching ratios of \(\bar{B}^0\rightarrow D^+\ell ^-\bar{\nu _\ell }\) (left) and \(\mathcal {R}({{D} ^+})\) (middle) in the case that no Coulomb correction is applied. The plot on the right shows the ratio \({\delta _{\mathrm{\mathrm QED}}} ({\mathcal {R}({{D} ^+})})/{\delta _{\mathrm{\mathrm QED}}} ({\mathcal {R}({{D} ^0})})\). The plots are obtained from simulated data (solid lines, which include the statistical uncertainty) and from Ref. [9] (dashed lines, filled with transparent colours when the uncertainties are significant)

Figure 3 (right) shows the ratio of QED corrections on \(\mathcal {R}({{D} ^+})\) over those on \(\mathcal {R}({{D} ^0})\). It is worth noting that both Photos and the calculation in Ref. [9] without Coulomb correction conserve isospin symmetry (\(\delta _{\mathrm{\mathrm QED}}\) values for \(\mathcal {R}({{D} ^+})\) and \(\mathcal {R}({{D} ^0})\) agree within the errors), while the Coulomb correction introduces an isospin-breaking term.

3 Effects on LHCb-like analysis

The comparison in the previous section only holds for values of \(E_{\gamma }\) up to 100 MeV. For higher energies, no calculations relevant to \(\mathcal {R}({D})\) are available.Footnote 2 Nevertheless, Photos generates photons with energies larger than those treated in Ref. [9], which is the range where the effect of SD photons might be relevant.

To study the effect of under- or overestimating radiative corrections in simulations used for measurements of \(\mathcal {R}({D})\), a simplified analysis is performed in an LHCb-like environment. Also here, the radiation emitted by the decay products of the \(D \) mesons is neglected because their contributions largely cancel out in the ratio.

The strategy of this study consists of fitting a data sample with templates describing the \({B} \rightarrow {D} {\mu ^-} {{\overline{\nu }} _\mu } \) and \({B} \rightarrow {D} {\tau ^-} {{\overline{\nu }} _\tau } \) components. The fits are performed with templates built under the hypothesis that no radiation with \(E_\gamma \) above a certain value \(E_{\mathrm{max}}\) is emitted. In particular, five \(E_{\mathrm{max}}\) values were chosen to cut on \(E_\gamma \): 100, 300, 500, 800 and 1500 MeV. The bias on \(\mathcal {R}({D})\), determined from these fits, is an indication of the importance of the simulation of the \(E_\gamma \) distribution in the high-energy region.

This analysis follows a strategy similar to the one used in Ref. [20], where \(\mathcal {R}({{D} ^*})\) is measured using a three-dimensional templated fit. The data samples, referred to as pseudo-experiments in the following text, are generated from a mixture of \({B} \rightarrow {D} {\mu ^-} {{\overline{\nu }} _\mu } \) and \({B} \rightarrow {D} {\tau ^-} {{\overline{\nu }} _\tau } \) decays, with radiative corrections generated by Photos. Here \(\mathcal {R}({D})\) is assumed to be 0.3 as predicted by the SM.

The variables used in the templated fit performed to extract the \({B} \rightarrow {D} {\mu ^-} {{\overline{\nu }} _\mu } \) and \({B} \rightarrow {D} {\tau ^-} {{\overline{\nu }} _\tau } \) yields from the pseudo-experiments are: the muon energy computed in the \(B \) meson rest frame, \(E_{\mu }\); the missing mass squared, \({m_{\mathrm{\mathrm miss}}^2} = (p _{{B}}-p _{{D}}-p _{\mu })^2\); and the squared four-momentum transferred to the lepton system, \({q^2} =(p _{{B}}-p _{{D}})^2\). The variables are binned as follows: four bins in \(q^2\) in the range \(-0.4<{q^2} <12.6\) GeV, 40 bins in \(m_{\mathrm{\mathrm miss}}^2\) between \(-2<{m_{\mathrm{\mathrm miss}}^2} <10\) GeV\(^2\), and 30 bins in muon energy in the range \(100<{E_{\mu }} <2500\) MeV, consistent with Ref. [20]. In this case study, only the signal (\({B} \rightarrow {D} {\tau ^-} {{\overline{\nu }} _\tau } \)) and normalisation (\({B} \rightarrow {D} {\mu ^-} {{\overline{\nu }} _\mu } \)) components are considered, while all backgrounds are ignored.

Basic selection requirements are applied to mimic the acceptance of the LHCb detector and its trigger following Ref. [21]. Both production and \(B \)-decay vertex positions are smeared to simulate the resolution of the LHCb detector. The resolution on the production vertices is 13\({\,\upmu \mathrm {m}}\) in x and y, and 70\({\,\upmu \mathrm {m}}\) in z direction. For the \(B \)-decay vertices a resolution of 20\({\,\upmu \mathrm {m}}\) in x and y, and 200\({\,\upmu \mathrm {m}}\) in z direction is used, after which the \(B \) direction is computed. The \(D \) mesons decay as \({{D} ^0} \rightarrow {{K} ^-} {{\pi } ^+} \) and \({{D} ^+} \rightarrow {{K} ^-} {{\pi } ^+} {{\pi } ^+} \), the \(\tau ^-\) lepton as \({\tau ^-} \rightarrow {\mu ^-} {{\overline{\nu }} _\mu } {{\nu } _\tau } \). The muons and all decay products from the \(D \) mesons are required to be in the pseudorapidity range between 1.9 and 4.9. In addition, the momentum of each of these particles is required to be larger than 5 GeV, and its component transverse to the beam direction must be larger than 250 MeV. The distance between the production and \(B \)-decay vertex should be at least 3 mm, similar to the requirements applied in a typical trigger selection.

Due to the missing neutrino and unknown effective centre-of-mass energy of the collision, the \(B \)-meson momentum cannot be reconstructed in an \(\mathcal {R}({H_c})\) analysis at LHCb. Therefore the momentum of the \(B \) meson in the z direction, \((p _{{B}})_z\), is approximated as \((p _{{B}})_z=(m_{{B}}/m_\mathrm{vis})(p _\mathrm{vis})_z\), where \(m_{{B}}\) is the \(B \) mass, and \(m_\mathrm{vis}\) and \((p _\mathrm{vis})_z\) are the momentum in the z direction and the mass of the visible decay products of the \(B \) meson, respectively. This directly follows the approach from Ref. [20]. After computing the \(B \) momentum with the above approximation and applying the selection criteria described in this section, \(q^2\), \(m_{\mathrm{\mathrm miss}}^2\) and \(E_{\mu }\) are calculated. The distributions for the signal and control samples are shown in Fig. 4. Even using this simplified detector description, these distributions show the same key features as those in Ref. [20].

Fig. 4
figure 4

Shapes of the \(q^2\), \(m_{\mathrm{\mathrm miss}}^2\), and \(E_{\mu }\) distributions after applying basic selection requirements and the rest frame approximation for the \({{\overline{B}{}} {}^0} \rightarrow {{D} ^+} {\ell ^-} {{\overline{\nu }} _\ell } \) decays where \({\ell ^-} = {\mu ^-}, {\tau ^-} \)

When applying cuts on \(E_{\gamma }\), the templates shapes change. This is most clearly seen in the distributions of \(m_{\mathrm{\mathrm miss}}^2\), shown for the \({\overline{B}{}} {}^0\) decay in Fig. 5. Especially in the \(\mu ^-\) decay mode the effect is large, altering the shape at high values of \(m_{\mathrm{\mathrm miss}}^2\). Since this feature is not present in the \(\tau ^-\) mode, this does not cancel when measuring the ratio \(\mathcal {R}({{D} ^+})\). For completeness, changes in the shape of \(q^2\) and \(E_{\mu }\) for the \(\mu ^-\) mode are shown in Appendix A.

Fig. 5
figure 5

Ratios of the cut \(m_{\mathrm{\mathrm miss}}^2\) distribution over the default \(m_{\mathrm{\mathrm miss}}^2\) distributions for the \(\bar{B}^0\rightarrow D^+\ell ^-\bar{\nu _\ell }\) decays, for various \(E_{\mathrm{max}}\) cuts. On the left for the \(\mu ^-\) and on the right for the \(\tau ^-\) mode

The number of events generated to simulate data is determined from the estimated number of events that LHCb gathered during their Run II data-taking period. The estimate takes into account the \(B \)-production cross-section at 13 TeV, branching fractions, and assumes the average reconstruction efficiency is the same as in Ref. [20]. This results in data samples of \(1.0\times 10^6\) and \(0.5\times 10^5\) for the \({{\overline{B}{}} {}^0} \rightarrow {{D} ^+} {\ell ^-} {{\overline{\nu }} _\ell } \) decays, and \(4.4\times 10^5\) and \(2.3\times 10^4\) for the \({{{B} ^-}} \rightarrow {{D} ^0} {\ell ^-} {{\overline{\nu }} _\ell } \) decays, where the first yield represents the \(\mu ^-\) sample, and the second the \(\tau ^-\) sample. In an actual analysis, the efficiencies for \({B} \rightarrow {D} \ell {{\overline{\nu }} _\ell } \) decays are likely higher than those for \({B} \rightarrow {D} ^{*}\ell {{\overline{\nu }} _\ell } \) decays, where \({D} ^{*}\) is reconstructed in the \({D} \pi \) decay mode.

The measured value of \(\mathcal {R}({D})\) is determined from two components. The first is the ratio of reconstruction efficiencies \(\varepsilon _{\mu }\) and \(\varepsilon _{\tau }\) for the \(\mu ^-\) and \(\tau ^-\) samples, respectively, which takes into account the selection requirements described earlier in this section. The second component is the fraction of semitauonic decays in the sample, \(f_{\tau }\), determined from the three-dimensional template fit (the absence of background events in the simulated samples implies that the fraction of \(\mu ^-\) and \(\tau ^-\) components add up to one). These are combined to measure \(\mathcal {R}({D})\) as

$$\begin{aligned} {\mathcal {R}({D})} = \frac{f_{\tau }}{1-f_{\tau }} \frac{\varepsilon _{\mu }}{\varepsilon _{\tau }}. \end{aligned}$$
(7)

The exercise of generating pseudo-experiments is repeated 10,000 times after which the spread of the measured values of \(\mathcal {R}({D})\) is taken as the statistical uncertainty.

The resulting values of \(\mathcal {R}({{D} ^+})\) and \(\mathcal {R}({{D} ^0})\) as a function of \(E_{\mathrm{max}}\) are shown in Fig. 6. From here it is clear that there is a significant effect in underestimating the QED radiative corrections which could be up to 0.02 for both \(\mathcal {R}({{D} ^+})\) and \(\mathcal {R}({{D} ^0})\) values, corresponding to a relative bias of 7.5%. The largest contribution to the observed bias is due to the fit fraction \(f_{\tau }\), which is strongly affected by the shapes of the \(\mu ^-\) templates. Instead, the ratio of efficiencies \({\varepsilon _{\mu }}/{\varepsilon _{\tau }}\) is only marginally dependent on \(E_{\mathrm{max}}\). However, this last statement holds only for this specific case study. Different sets of selection cuts or different experimental environments could indeed introduce a significant bias also in the ratio of efficiencies. The observed bias can be understood when looking at the \(m_{\mathrm{\mathrm miss}}^2\) distribution in Figs. 4 and 5. When cutting on \(E_{\gamma }\), part of the tail of the \(\mu ^-\) distribution is removed, which is compensated by a higher \(\tau ^-\) fraction in the fit.

Fig. 6
figure 6

Values of \(\mathcal {R}({D})\) as a function of \(E_{\mathrm{max}}\) measured in a simplified LHCb-like analysis of \(\mathcal {R}({{D} ^+})\) (left) and \(\mathcal {R}({{D} ^0})\) (right). The error bars reflect the statistical uncertainties from the generated MC samples. The blue bands correspond to fit results obtained with no \(E_{\mathrm{max}}\) cuts

Fig. 7
figure 7

Coulomb correction as a function of \(q^2\), \(m_{\mathrm{\mathrm miss}}^2\), and \(E_{\mu }\) for \({{\overline{B}{}} {}^0} \rightarrow {{D} ^+} {\ell ^-} {{\overline{\nu }} _\ell } \) decays, where \({\ell ^-} = {\mu ^-}, {\tau ^-} \)

In an actual analysis there are radiative corrections in MC and radiated photons in data. Therefore, it is useful to check the above approach using an alternative strategy. In this case, the templates include all QED corrections predicted by Photos while an \(E_{\mathrm{max}}\) cut is applied on the pseudo-experiments. This approach leads to an overestimate on the QED corrections, resulting in a negative bias on the \(\mathcal {R}({D})\) values. The results for \(\mathcal {R}({{D} ^+})\) and \(\mathcal {R}({{D} ^0})\) as a function of \(E_{\mathrm{max}}\) are reported in Appendix A. The corrections are of the exact same size as those in the baseline approach.

It is worth to note that, despite the fact that LHCb does not cut explicitly on \(E_{\gamma }\) in its analyses, indirect cuts on the total radiated energy are applied through e.g. requirements on isolation variables or inefficient reconstruction algorithms for low momentum particles. This could alter the observed bias if \(E_{\gamma }\) is not simulated correctly.

These studies show that radiative corrections play a crucial role in \(\mathcal {R}({D})\) measurements. Since part of these corrections are already simulated in Photos, the above exercise shows the effect of a worst-case scenario. Nevertheless, additional effects such as the Coulomb correction, as detailed in next section, or the calculation for energies greater than 100 MeV are becoming fundamental in view of the increased experimental precision expected in the coming years. Also, these quantitative effects strongly depend on explicit or implicit cuts on radiative photons and must be carefully evaluated for each analysis measuring \(\mathcal {R}({D})\).

3.1 Coulomb correction

Beyond affecting the SM prediction of \(\mathcal {R}({{D} ^+})\), the Coulomb correction impacts the experimental results by changing the shape of the fit templates. This is evaluated in the \({{\overline{B}{}} {}^0} \rightarrow {{D} ^+} {\ell ^-} {{\overline{\nu }} _\ell } \) decay by weighting each event by the term \(\varOmega _{\mathrm{C}}\).Footnote 3 The changes in the shape of the \(q^2\), \(m_{\mathrm{\mathrm miss}}^2\) and \(E_{\mu }\) distributions are shown in Fig. 7. While for the \(\mu ^-\) mode \(\varOmega _{\mathrm{C}}\) is mostly constant, for the \(\tau ^-\) mode there is a dependence on each of the three variables due to the smaller relative velocity. To quantify the effect of the Coulomb correction, the above analysis is repeated, without including any \(E_{\mathrm{max}}\) cuts. The Coulomb correction is applied to the pseudo-experiments, but not to the fit templates, resulting in a relative shift of about – 1.0% on \(\mathcal {R}({{D} ^+})\). This effect can even be amplified by selecting certain regions of phase space.

4 Conclusions and recommendations

The work in Ref. [9] describes QED corrections which are not fully included in Photos. These corrections affect the semimuonic and semitauonic modes differently at the level of a few percent. Ignoring the Coulomb correction, there is more radiated energy in the calculation in Ref. [9] than in Photos for the \({\overline{B}{}} {}^0\) decays, while this is the other way around for the \({{B} ^-} \) decays. In the ratio \(\mathcal {R}({D})\), this small discrepancy mostly cancels out. However, the main difference between the QED corrections on \(\mathcal {R}({{D} ^+})\) and \(\mathcal {R}({{D} ^0})\), which is up to 1%, is due to the Coulomb correction that only affects \(\mathcal {R}({{D} ^+})\).

Coulomb interactions are not simulated by Photos and mainly affect the kinematics of semitauonic decays, which in turn influence the shape of distributions used to determine the signal yields in LHCb, BaBar, and Belle analyses. These effects can alter values of \(\mathcal {R}({D})\) up to 1% in an LHCb-like analysis, and should be evaluated precisely for each measurement.

Using a simplified LHCb-like analysis, it is shown that over- or underestimating radiative corrections could bias measurements of \(\mathcal {R}({D})\) up to 7% in an extreme case. This results in a bias of 0.02 on the value of \(\mathcal {R}({D})\), and should be studied further when performing these types of measurements, including a realistic evaluation of cuts on \(E_{\gamma }\). These effects could potentially be enhanced in measurements from Belle  II [22] where the resolution on the kinematic variables is better than at LHCb.

When measuring values of \(\mathcal {R}({D})\) with higher precision, additional calculations of QED corrections for \({B} \rightarrow {D} \ell {{\nu } _\ell } \) decays are necessary. Especially calculations involving high-energy and structure-dependent photons are currently mostly missing.