Hunting for B + → K + τ + τ - imprints on the B + → K + μ + μ - dimuon spectrum.

We investigate the possibility of indirectly constraining the B + → K + τ + τ - decay rate using precise data on the B + → K + μ + μ - dimuon spectrum. To this end, we estimate the distortion of the spectrum induced by the B + → K + τ + τ - → K + μ + μ - re-scattering process, and propose a method to simultaneously constrain this (non-standard) contribution and the long-distance effects associated to hadronic intermediate states. The latter are constrained using the analytic properties of the amplitude combined with data and perturbative calculations. Finally, we estimate the sensitivity expected at the LHCb experiment with present and future datasets. We find that constraints on the branching fraction of O ( 10 - 3 ) , competitive with current direct bounds, can be achieved with the current dataset, while bounds of O ( 10 - 4 ) could be obtained with the LHCb upgrade-II luminosity.

Introduction

In recent years, discrepancies between the observed values and the Standard Model (SM) predictions of the lepton-flavour universality (LFU) ratios [1-5] and [6-9], characterizing the semileptonic transitions and , have sparked great interest. The pattern of anomalies seems to point to intriguing new-physics (NP) scenarios, with possible connections to the SM flavour puzzle. A large class of NP models proposed to explain these hints of physics beyond the SM, and in particular those aiming for a combined explanation of the and anomalies, imply dominant couplings to third-generation fermions, which should also enter other semileptonic b-quark decays.

A general expectation, confirmed by many explicit NP constructions, is that of a large enhancement of transitions (see e.g. [10-17]). While flavour-changing neutral-current (FCNC) decays with muon and electron pairs have been observed both at the exclusive and at the inclusive level, probing rare decays with a pair in the final state is experimentally very challenging. The current experimental limits for all processes mediated by the amplitude are still very far from the corresponding SM predictions [18, 19], leaving the NP expectation of possible large enhancements unchallenged.

In this work we investigate the possibility of indirectly constraining the amplitude via its imprint on the dimuon spectrum. In presence of a large NP enhancement, the amplitude would induce a distinctive distortion of the spectrum via the (QED-induced) re-scattering process [10]. The latter has a discontinuity at (), namely at the threshold where the tau leptons can be produced on-shell. This gives rise to a “cusp” in the dimuon-invariant mass spectrum, which could in principle be detected with sufficient experimental precision. More generally, the lightness of the -leptons implies a well-defined deformation of the spectrum, which is determined only by the analytic properties of the re-scattering amplitude.

It should be stressed that the phenomenon we are considering here is different from the QED mixing between dimension-six FCNC operators with different lepton species analysed in Ref. [20]. If NP is heavy and the amplitude is strongly enhanced, the operator mixing can give rise to sizable modifications of the Wilson coefficients of the dimension-six effective Hamiltonian relevant to decays (). However, this phenomenon cannot be distinguished in a model-independent way from other NP effects of short-distance origin (at least using low-energy data only). On the contrary, the non-local effect we are interested in can be unambiguously attributed to the re-scattering of light intermediate states characterised by the tau mass, hence it can be translated into a model-independent constraint on the amplitude.

The main difficulty in extracting such bound is obtaining a reliable description of the dimuon spectrum within the SM, or better in the limit where the re-scattering is negligible. This is non trivial, given that the spectrum is plagued by theoretical uncertainties originating from form factors and hadronic long-distance contributions. While the former are smooth functions in the region of interest and can be well described using lattice QCD [21, 22] and/or light-cone sum rules [23], long-distance effects induced by intermediate hadronic states, such as the charmonium resonances, are more problematic. They are genuine non-perturbative effects and introduce physical discontinuities below and above the threshold. Far from the resonance region, these effects can be estimated using perturbative constraints derived at , with , combined with a or expansion to incorporate the leading non-perturbative corrections [24, 25]. However, this approach is not suitable for our purpose, which requires a reliable description of the whole spectrum, and in particular of the resonance region. To achieve this goal, we adopt a data-driven approach which takes full advantage of the known analytic properties of the amplitude: knowing the precise location of all one- and two-particle hadronic thresholds, we use subtracted dispersion relations to describe the -dependence of the whole spectrum in terms of a series of (-independent) hadronic parameters, which are fitted from data. This method can be considered a generalisation of the approaches proposed in Ref. [26] and, to some extent, in Refs. [27-29], with a few key differences, the most notable ones being the use of subtracted dispersion relations and the explicit inclusion of two-particle thresholds. To reduce the number of independent free parameters, perturbative constraints derived from the low- region are also implemented. Proceeding this way we obtain a description of the spectrum that is flexible enough to extract the non-perturbative parameters characterising the various hadronic thresholds from data, but retains a significant predictive power in the smooth region within and below the two narrow charmonium states, allowing us to set useful constraints on the re-scattering.

The method we propose is particularly well suited for the LHCb experiment, which has already collected a large sample of events and has an excellent resolution in the dimuon spectrum [30]. In order to estimate the sensitivity of LHCb in view of the full run II dataset, we generate pseudo-experiments based on the yields and amplitudes obtained in Ref. [30], and calculate the expected limit under the background-only hypothesis using the CLs method [34].

The paper is organised as follows: in Sect. 2 we introduce the theoretical framework necessary to describe the dimuon spectrum within and beyond the SM, separating short-distance contributions (Sect. 2.1), long-distance contributions due to intermediate hadronic states (Sect. 2.3), and long-distance contributions due to the re-scattering (Sect. 2.4). The analysis of the LHCb sensitivity is presented in Sect. 3. The results are summarised in the Conclusions.

Theoretical framework

Effective Hamiltonian and differential decay rate

The dimension-six effective Langrangian describing transitions, renormalized at low energies [], can be decomposed aswhere the leading FCNC effective operators are defined asand the most relevant four-quark operators () asWithin the class of models we are considering, all relevant NP effects are encoded in the values of the Wilson coefficients . Given the normalisation in Eq. (1), and are real and O(1) within the SM, whereas (see Ref. [26] for the precise values of the Wilson coefficients and the complete basis of operators).

The matrix elements are non-vanishing at the tree level only in the case of the FCNC operators (with ). Considering only the contribution of the FCNC operators, which can be expressed in terms of the form factors, the decay rate can be written as:where is the kaon momentum in the B-meson rest frame, , and with are the vector, scalar and tensor form factors.

Non-local contributions: general considerations

The non-local contributions generated by the non-leptonic operators in and by the operator can be encoded in Eq. (4) by replacing with a -dependent function:where denotes the non-local contributions corresponding to the intermediate state , which can annihilate into a dimuon pair via a single-photon exchange.

The functions and encode non-perturbative hadronic contributions, which cannot be estimated reliably in perturbation theory, at least in a large fraction of the accessible spectrum. Adopting a notation similar to that of Ref. [26], we can express aswhere is defined by the gauge-invariant decomposition of the following non-local hadronic matrix elementwith . The function , containing the contribution of the subleading non-leptonic operators in , is defined in a similar way via the replacementOur main strategy is to write the non-perturbative functions and using hadronic dispersion relations. More precisely, for the leading charm contribution we consider one- (1P) and two-particle (2P) intermediate states (see Fig. 1), using dispersion relations subtracted at , while for the subleading function we consider only one-particle intermediate states and use unsubtracted dispersion relations. We stress that these dispersion relations, and the corresponding decomposition of , are not exhaustive of all the discontinuities of the four-point function in (7). However, our goal is not to determine completely via dispersion relations, but only to describe its functional dependence with respect to the variable .

Fig. 1

Diagrammatic representations of the long-distance contributions to . The left-hand side depicts the exchange of a single vector resonance. The graph on the right-hand side shows the contribution from two-particle intermediate states

Given these considerations, in Eq. (5) is finally decomposed according towith . In the next section we analyse the structure of , , and in detail. The expression of , which is the only term in Eq. (9) that can be fully evaluated in perturbation theory, is given in Sect. 2.4.

Long-distance hadronic contributions

The general structure of the subtracted dispersion relation used to determine isThe function is the spectral-density function describing the hadronic states , characterized by valence charm-quarks and invariant mass s, contributing as real intermediate states in the re-scattering . As noted before, we decompose into one- and two-particle intermediates states, ,neglecting the phase-space suppressed contribution with three or more particles.

Charmonium resonances

For the sake of simplicity, in Eq. (11) we have treated the single-particle states as infinitely narrow resonances. The effect of finite widths can be incorporated via Breit–Wigner functions, yieldingwhere the sum runs over all the charmonium vector resonances in the accessible kinematical range. Here and are real parameters which must be determined from data, similarly to what has been performed by the LHCb collaboration in [30].

We stress that our scope is not to compute the employing a facotrization hypothesis for the hadronic matrix element in (7), combined with experimental data on , as originally proposed in [31]. As clearly shown in [27], this approach leads to rather inconsistent results. We refrain from any attempt to compute the , which are genuine non-perturbative hadronic matrix elements: we treat them as free parameter which need to be determined by data.

In principle, both and are -dependent functions. However, if this dependence is smooth around the resonance poles (or the two-particle thresholds, for the analog parameters that we define in Sect. 2.3.2), we can treat them as constant terms for the purpose of our analysis. On general grounds, we expect the characteristic scale determining the -variation of the appearing in to be twice the charm mass. This -dependence is indeed related to the invariant mass of the hadronic intermediate states that can mix into the charmonia via re-scattering processes. Incidentally, we note that this is also what one would infer from a perturbative estimated of the using the factorization hypothesis. In view of this argument, we believe that it is a good approximation to treat the as constant terms in our analysis. As stated above, our goal is not to compute these parameters but only to fit them from data in order to have a sufficiently general description of the long-distance part the amplitude, able to reproduce all the known discontinuities related to hadronic intermediate states.

The fitted ’s can be put in one-to-one correspondence with the product of the and branching fractions viaThe expression (13) differs from the decomposition adopted in Ref. [30] by the term, which arises from the subtraction procedure in the dispersion relation. On the one hand, the use of subtracted dispersion relations for the charm contribution is necessary to ensure the convergence of the integral in the two-particle intermediate states (see Sect. 2.3.2). On the other hand, choosing the subtraction point at allows us to decouple the determination of the resonance parameters of the spectrum from the overall normalisation of the rate, and hence from the determination of from data. The price to pay is the appearance of the undetermined constant term in Eq. (9). This term plays no role in the description of the dimuon spectrum, but is relevant for the extraction of the value of . To this purpose, we note that the estimate presented in Ref. [26], which is based on a expansion and also takes next-to-leading corrections on the pure partonic result into account (see Sect. 2.3.4), yieldswhich is about of .

Two-particle intermediate states

Proceeding in a similar way, we can decompose the two-particle contributions aswhere are normalised spectral densities for the two-body intermediate states characterised by the threshold .

While we do not have a precise estimate of these spectral densities at generic kinematical points, an excellent description of their behaviour around the respective thresholds is obtained by approximating them with powers of the Källén function, with an exponent determined by the lowest partial wave allowed in the re-scattering. This is because higher-order partial waves, characterised by higher powers of the Källén function, are both phase-space suppressed and, most importantly, give rise to a less singular behaviour at the threshold. From angular momentum conservation we can then determine the leading partial wave and obtain the following estimates for the normalised spectral densities of the two-particle intermediate states of lowest mass:In the case of the intermediate state we have replaced the complete expression depending on both masses with a simplified one depending only on , which provides an excellent approximation. With these estimates in place we find:withand

It is worth noting that, while the lowest threshold is at , the contribution from the intermediate state is the only one which can occur in the S-wave, corresponding to a singular (square-root) behaviour at the threshold (see Fig. 2).

Fig. 2

Real (solid) and imaginary (dashed) parts of the normalised hadronic two-particle contributions to , as defined in Eq. (16)

Light resonances

The remaining hadronic contribution we need to estimate is , defined by Eqs. (6) and (7) via the replacement (8). The Wilson coefficients of the effective operators appearing in are either loop- or CKM-suppressed. As a result, we can limit ourselves to include only one-particle hadronic intermediate states. In principle, such operators describe transitions also to states with valence charm quarks; however, since we fit the hadronic coefficients from data, these terms are naturally absorbed in the appearing in . We are thus left only with vector resonances containing light valence quarks. Among them, we can further restrict the attention to the , , and resonances, since the leptonic decay rates of the heavier states are very small.

There is no clear advantage in using subtracted vs. unsubtracted dispersion relations in describing the contributions of the light vector resonances. The convergence of the dispersive integrals does not pose a problem, and the subtraction at is not particularly useful since the light-quark contributions are in a non-perturbative regime at . However, when fitting data, the subtraction at retains the advantage of decoupling the determination of the spectrum from that of the Wilson coefficient. As default option, we adopt unsubtracted dispersion relations. As discussed in Sect. 2.3.5, checking the stability of the result using subtracted vs. unsubtracted dispersion relations for the light vector resonances provides an estimate of the “model error” of the proposed approach.

Given these considerations, we decompose asin perfect analogy with the decomposition adopted in Ref. [30] for these light states.

Theoretical constraints on the hadronic parameters

The hadronic decompositions in Eqs. (13), (18) and (21) contain 12 free complex parameters: 6 in , 3 in , and 3 in . In principle, since they correspond to different functional forms, they could all be fitted from data. In practice however, an unconstrained fit would leave significant degeneracies in the parameter space. It is therefore useful to restrict the variability of such parameters using theoretical constraints. In the following we discuss three conservative conditions which can be imposed using perturbative arguments.

I. Constraint on the slope of at .

The lowest-order perturbative estimate of is obtained by factorising the matrix element in Eq. (7) and computing the charm-loop at :This expression is certainly not a good approximation of close to the resonance region; however, it is expected to provide a reasonable approximation at , up to corrections. We can thus use it to set bounds on the slope of in the vicinity of . The perturbative result implieswhere the numerical value has been obtained setting and  GeV. According to the analysis of Ref. [26], the inclusion of corrections (which involve new hadronic matrix elements) modifies the above prediction to . Given these considerations, in the numerical analysis we employ the following constraintswhere we slightly enlarged the error from (23), such that the range covers the difference between the central value in (23) and the one including corrections estimated in Ref. [26].

II. Upper bound on the in .

The comparison of the perturbative result with also allows us to define the natural range for the parameters, which are poorly constrained by data. Focusing the attention on the leading S-wave contribution, it turns out that the perturbative quark loop can be saturated by the meson loop, in the limit , setting On general grounds, each of the exclusive meson contributions should be significantly smaller than the inclusive quark contribution. As a result, in the following we set the upper limitIII. Upper bound on .

Using an unsubtracted dispersion relation and taking into account only one-particle intermediate states for the light-quark contributions implies for large , while . More precisely, one finds a power-like suppression of the type at large , whereas is not parametrically suppressed by any scale ratio. However, since the Wilson coefficients entering are strongly suppressed, either by loop factors or by subleading CKM factors, cannot be too large. Parametrically we expectTaking the size of the into account, we set the conservative bound1which should be interpreted as a constraint on the relative phases of the light resonances.

Estimate of the “model error”

Despite not being entirely dictated by first principles, the parameterisation of long-distance hadronic contributions discussed so far contains all the relevant one- and two-particle discontinuities of the amplitude, with free coefficients to be fixed by data. It should therefore provide a sufficiently general (and unbiased) description of the impact of hadronic contributions on the spectrum. Still, it may be worthwhile to assess whether the proposed parameterisation influences the extraction of information on and, correspondingly, the extraction of a bound on . An estimate of this “model error” can be obtained by examining the stability of the obtained bound on under small variations of the model assumptions. The latter include: (i) the use of subtracted vs. unsubtracted dispersion relations for the light resonances; (ii) the use of -dependent widths for both charmonium and/or light resonances; (iii) strengthening or relaxing the theoretical constraints in Eqs. (24), (25), and (27).

Tau-lepton contribution

The contribution from the intermediate -leptons can be computed in perturbation theory, yieldingwith the functions defined in Eq. (19). The functional form is identical to the one of the perturbative charm contribution and, to a large extent, to the one of the contribution, illustrated in Fig. 2. However, the cusp is located at , sufficiently well separated from the various hadronic thresholds.

In principle, the short-distance amplitude does not need to be controlled by the CKM matrix in a generic NP model. However, in most realistic scenarios the weak phases of all amplitudes are aligned to the SM one, implying . In the following, we adopt this (motivated) simplifying assumption.

An estimate of the maximal allowed size of can be derived from the experimental upper bound on at 90% CL by Babar [19], which is more than four orders of magnitude larger than [32]. Neglecting the contributions from operators other than and , we findIn the case () the Babar result then implies , to be compared to . As we discuss below, saturating this bound leads to a pronounced ditau cusp in the spectrum (see Fig. 3), opening the possibility of extracting a more stringent bound on from a precise measurement of the dimuon spectrum.

Fig. 3

Example pseudodata expected from the full run II dataset collected by the LHCb experiment assuming the SM. The distribution expected if the branching fraction were present at the current experimental limit of is overlaid

Analysis of the expected sensitivity at LHCb

In order to assess the sensitivity to the branching ratio at the LHCb experiment, we generate pseudo-experiments corresponding to the signal yields obtained in Ref. [30] and scaled to the full run II dataset, taking into account the collected luminosity and b-hadron cross-section increase at 13 TeV [33]. This leads to around 40,000 candidates (cutting the two narrow resonances). As the efficiency is reasonably flat as a function of dimuon mass and the background level is very low, we neglect these effects. Figure 3 shows the fit model with a dataset generated at the expected yield. This illustrates the visible sensitivity to a hypothetical signal component generated according to the current experimental limit [19].

The size and phase of the one-particle resonant contributions are determined from the branching fractions reported in Ref. [33], which are used to determine the initial values of and for the data to be generated. Due to the complicated experimental resolution effects near the and resonances, the regions and are excluded from the fit and the phase differences associated with these resonances are constrained to the uncertainties in Ref. [30]. Outside of this region, finite-resolution effects in are ignored as all the components are broad. In order to mimic the sensitivity one would have when fitting the data, Gaussian constraints are applied to the and resonant parameters according to the uncertainties reported in Ref. [33].

For the two-particle hadronic contribution, we conservatively allow the magnitude and phase of them to vary in the fit. As the shape of the and spectral densities are very similar, we combine them with an equal contribution to avoid large correlations in the fit.

The form factor uncertainties are taken from Ref. [22] and are implemented in the fit as a multivariate Gaussian constraint. The data slightly helps constrain the form factor parameters, but this affects the sensitivity on only in a mild way.

The expected sensitivity on the contribution is determined using the CLs method [34]. The sensitivity with the current dataset is reported in Table 1, along with two other potential future scenarios corresponding to the LHCb upgrade-II luminosity and a hypothetical improvement of the form factor uncertainties by a factor of three. The estimated sensitivity utilising the run I–II datset corresponds to a limit on the branching ratio which is slightly more stringent than the current constraints placed by the BaBar collaboration and is expected to compete with the projected sensitivity of the Belle-II experiment when more data is collected.

Table 1

Sensitivity to according to various LHCb scenarios []

Scenario	C9τ (90% CL)	B (C10τ=-C9τ)	B (C10τ=0)
Run I–II dataset	533	2.7×10^-3	0.8×10^-3
Run I–V dataset	139	1.8×10^-4	0.5×10^-4
Run I–II dataset, improved form factors	533	2.7×10^-3	0.8×10^-3
Run I–V dataset, improved form factors	127	1.5×10^-4	0.5×10^-4

Interplay with the hadronic contributions

The re-scattering leads to two main features in the dimuon spectrum: (i) the cusp in-between the and resonances, and (ii) a distortion in the shape of the spectrum before the two resonances. The effect after the peak is less relevant since in that region the spectrum is rather discontinuous due to the various one- and two-particle thresholds. In order to investigate the sensitivity to the cusp feature, we have performed a fit limited to the region between the and resonances: this leads to a sensitivity to diluted by a factor of four. We thus conclude that is the deformation of the spectrum, in particular before the , that generates the largest sensitivity to . This implies that neglecting the resolution is justified.

Since the deformation of the spectrum at low plays a relevant role, we deduce that the assumed shape of the charmonium contribution is an important ingredient in constraining the signal. The component which most closely resembles the signal is the contribution from two-particle hadronic intermediate states. This is reflected in a correlation coefficient between this amplitude and the signal of about 0.6. However, the two are clearly distinct given the different location of the thresholds. We also explicitly checked that the theoretical constraints described in Eqs. (24) and (25) do not affect our sensitivity estimate: the best fit value of the two-particle hadronic contribution lies far from these bounds. This leaves open the possibility of a further increase of sensitivity with more stringent constraints on the two-particle hadronic contribution, which could be derived using data. We finally note that the correlation between the hadronic contribution and the signal is fully taken into account in the sensitivity estimates using the CLs method.

Conclusions

If the branching ratio were significantly enhanced over its SM value, it would induce a peculiar distortion of the spectrum, characterised by a cusp at and by a distortion of the dimuon distribution. In this work we have proposed a method that uses this effect as a tool to extract a bound on from future precise measurements of .

A necessary ingredient to achieve this goal is a reliable description of the dimuon spectrum, within the SM, in the full kinematical range, especially in the region before and within the narrow charmonium states. As we have shown, this can be obtained by means of a data-driven approach which takes full advantage of the known analytic properties of the decay amplitude, supplemented by robust theoretical constraints. Our approach differs from previous attempts of including non-local hadronic contributions to the decay amplitude by three main points: (i) the use of dispersion relations subtracted at for the charmonium states; (ii) the inclusion of two-particle thresholds; (iii) the use of short-distance constraints at low to reduce the number of free parameters. In this way one separates the problem of the normalisation of the rate, and the corresponding extraction of short-distance Wilson coefficients, from the problem of obtaining a reliable description of the dimuon spectrum. While within our approach there is no significant progress on the first problem, there is a tangible advantage on the second one. The parameterisation of the amplitude we propose is flexible enough to allow the extraction of all the relevant parameters characterising hadronic thresholds in the dimuon spectrum from data, while retaining significant predictive power in the smooth region within and below the two narrow charmonium resonances. This fact is the key property which allows us to set useful constraints on the re-scattering from future precise measurements of .

The method we have proposed is particularly well suited for the LHCb experiment, which has already collected a large sample of events and has an excellent resolution in the dimuon spectrum [30]. As we have shown, the data already collected in run II should allow to set a bound on of , competitive with current direct bounds (see Table 1). Bounds of could be obtained with the LHCb upgrade-II luminosity.