Long-distance effects in $B\to K^*\ell\ell$ from Analyticity

We discuss a novel approach to systematically determine the long-distance contribution to $B\to K^*\ell\ell$ decays in the kinematic region where the dilepton invariant mass is below the open charm threshold. This approach provides the most consistent and reliable determination to date and can be used to compute Standard Model predictions for all observables of interest, including the kinematic region where the dilepton invariant mass lies between the $J/\psi$ and the $\psi(2S)$ resonances. We illustrate the power of our results by performing a New Physics fit to the Wilson coefficient $C_9$. This approach is systematically improvable from theoretical and experimental sides, and applies to other decay modes of the type $B\to V\ell\ell$, $B\to P\ell\ell$ and $B\to V\gamma$.


I. INTRODUCTION
B → K * decays are sensitive to modified shortdistance physics from sources beyond the Standard Model (SM), and a great deal of experimental and theoretical work has been devoted to extract short-distance information from them. However, long-distance physics within the SM also contributes significantly to the decay, and its effects are very difficult to assess reliably from first principles. On the other hand, tighter experimental constraints from increasingly precise measurements of b → s processes have significantly limited the size of allowed New Physics (NP) effects in B → K * , which are now comparable to current SM uncertainties. Thus, our inability to reliably constrain these long-distance contributions to acceptable levels stands in the way of obtaining unambiguous information on physics beyond the SM.
The B → K * decay is conveniently described by the K * transversity amplitudes (λ =⊥, , 0) where C 7,9,10 are short-distance Wilson coefficients, and N λ are normalization factors. The non-trivial matter from the theory point of view is the determination of the "local" and "non-local" long-distance effects encoded in the functions F (T ) λ (q 2 ) and H λ (q 2 ), respectively, which depend on the dilepton invariant mass squared q 2 .
Here we focus on the functions H λ (q 2 ), which are related to the contribution from 4-quark and chromomagnetic operators in the Weak Effective Hamiltonian, and emerge from the "non-local" matrix element where p = q + k, η is the polarization vector of the K * , and K(x, y) is a bi-local operator. The most relevant contribution to this matrix element in the SM arises from the current-current operators O 1,2 , since they come with large Wilson coefficients. In this letter we consider only this contribution -the so-called "charm-loop effect" -for which the object K µ (x, y) is given by: with j µ em (x) = q Q qq (x)γ µ q(x) the electromagnetic current. The scalar functions H λ (q 2 ) are given by the Lorentz decomposition: where S αµ λ are a set of structures given in the appendix. In the heavy b-quark limit and for very small q 2 , the functions H λ (q 2 ) factorize into non-perturbative form factors and light-cone distribution amplitudes, up to perturbatively calculable "hard" functions [12]. However this perturbative expansion breaks down when q 2 approaches 4m 2 c , leading to questionable predictions for q 2 6 GeV 2 . The integral in Eq.
(2) is in fact dominated by the region x 2 (2m c − q 2 ) −2 [13], so for q 2 4m 2 c one may expand the operator K µ (x, 0) around x 2 = 0 (a light-cone operator-product expansion, or LCOPE). This leads to an expansion of Eq. (2) in powers of (2m c − q 2 ) −1 , with matrix elements of operators that are non-local only along the light cone. This theory framework has been worked out up to NLO in α s [12,14] and including subleading terms in the LCOPE [13], and can be safely applied for q 2 4m 2 c (preferably at q 2 < 0). However, reliable predictions for larger values of q 2 remain a challenge.
In this letter we consider a consistent, modelindependent and systematically-improvable approach to determine the long-distance contributions H λ (q 2 ) to B → K * in the region q 2 14 GeV 2 . It provides genuine SM predictions even in the presence of NP in semileptonic operators. In addition, this approach provides access to the inter-resonance region 10 GeV 2 q 2 13 GeV 2 . The idea is the following: We determine the analytic properties of the functions H λ (q 2 ) in the complex plane, and use this information to write down general and modelindependent parametrizations. We then use two pieces of information to constrain the parametrized functions: data on B → K * J/ψ and B → K * ψ(2S), which is independent of NP in semileptonic operators; and theory at q 2 < 0, where it is reliable. This method, which builds upon Refs. [13,15], gives the most reliable and consistent a-priori determination of the functions H λ (q 2 ) to date. We use these results to compute SM predictions (assuming no NP in O 1,2 ), and to perform a NP fit to C 9 . All our numerical computations are performed with the help of EOS [16], which has been modified for this purpose [17].

II. ANALYTIC STRUCTURE AND PARAMETRIZATION
It is a standard assumption in quantum field theory that the only analytic singularities of a correlation function -as a complex function of all its complexified kinematic invariants -are those required by unitarity [18]. Unitarity, in turn, relates analytic singularities with onshell intermediate states: poles for one-particle states, and branch cuts for multi-particle states. Thus, the analytic structure of a correlation function can be learned by analysing its on-shell cuts.
In the case at hand, inspection of the correlation function (2) reveals the following analytic properties of the scalar functions H λ (q 2 ): On-shell cuts in the variable q 2 include: two poles at q 2 = M 2 J/ψ 9 GeV 2 and q 2 = M 2 ψ(2S) 14 GeV 2 corresponding to one-particle intermediate states through B → K * ψ n (→ + − ), with ψ 1 = J/ψ and ψ 2 = ψ(2S); a branch cut starting at q 2 = t + ≡ 4M 2 D corresponding to two-particle intermediate states through B → K * [DD](→ + − ), plus other "cc" cuts with higher thresholds; and "light-hadron" branch cuts starting at q 2 0 from intermediate states such as B → K * [3π](→ + − ), which include finite-width effects of J/ψ and ψ(2S). The effects of these "light-hadron" cuts are OZI suppressed [19][20][21]. Given the limited precision of current data, we will neglect these OZI suppressed contributions, as well as the effects of other light hadron resonances.
On-shell cuts in the variable (q + k) 2 (the "forward" or "decay" channel) include branch cuts from intermediate states such as B →DD s → K * + − . The physical point (q + k) 2 = M 2 B lies on these cuts, which implies that the functions H λ (q 2 ) are complex-valued for all values of q 2 . But this imaginary part is not associated with any singularity in the variable q 2 . Thus, one can write H λ (q 2 ) = H (re) satisfying the analytic properties of the previous point as functions of q 2 , and obeying the same dispersion relation.
These properties can be exploited to write down a general parametrization for the correlator consistent with unitarity. A convenient way to do so is to re-express the functions H λ (q 2 ) in terms of the "conformal" variable z: ). This transformation maps the cc branch cut in the q 2 plane to the unit circumference |z| = 1, and the entire first Riemann sheet in the q 2 plane to the interior of the unit circle |z| < 1. Our choice for t 0 implies that within the relevant interval −7 GeV 2 ≤ q 2 ≤ M 2 ψ(2S) , |z| < 0.52. The approach now resembles and is inspired by the zparametrization used for the form factors [22,23]. The functions H λ (z) ≡ H λ (q 2 (z)) are meromorphic in |z| < 1, with two simple poles at z J/ψ ≡ z(M 2 J/ψ ) 0.18 and z ψ(2S) ≡ z(M 2 ψ(2S) ) −0.44. Therefore, multiplying by the corresponding zeroes will give an analytic function in |z| < 1 that can be Taylor-expanded around z = 0. This expansion should converge reasonably well in the region of interest, where |z| < 0.52. This is the basis of our proposed parametrization.
In order to assure that the leading terms in the expansion will capture most of the features of the function (thus improving convergence), we use two more pieces of information: First, the correlator inherits all the singularities of the form factor (e.g. the M B * s pole), and the leading OPE contribution to the correlator is indeed proportional to the form factor. Therefore it is better to parametrize the ratios H λ (q 2 )/F λ (q 2 ) instead. Second, the poles should not modify the asymptotic behaviour. This is achieved by introducing appropriate "Blaschke factors" [22]. All in all, we propose the following parametrization: withĤ where α (λ) k are complex coefficients, and the expansion is truncated after the term z K . This truncation unavoidably introduces some model dependence. The maximum value that can be chosen for K will depend on the available set of experimental measurements and theory inputs.
We use the data to produce two sets of five pseudoobservables (three magnitudes and two relative phases on each resonance): where r ψn λ ≡ Res .
The numerical values for these pseudo-observables are obtained from the posterior-predictive distributions of a Bayesian fit. The inputs for this fit and the results are provided for completeness in the appendix. These pseudo-observables will act as constraints on the parameters of the correlators at z = 0.18 and z = −0.44.

IV. THEORY CONSTRAINTS
At q 2 < 0 the functions H λ can be calculated with the current approaches for the large recoil region. We  use QCD-factorization at next-to-leading order in α s , including the form factor terms and hard-spectator contributions [12,31]. In addition, we include [32] the soft-gluon correction calculated via a LCSR in Ref. [13]. For the form factors we use the results from the LCSR with B-meson distribution amplitudes [2], in order to have a mutually consistent description of form factors and non-local contributions and benefit from theoretical correlations among both. In this way we compute the ratios H λ (q 2 )/F λ (q 2 ) at the points q 2 = {−7, −5, −3, −1} GeV 2 . These ratios are used as pseudoobservables to constrain the parameters in Eq. (6) at z = {0.52, 0.50, 0.48, 0.46}. Further details and results are presented for completeness in the appendix. We emphasize that no theory is used at q 2 ≥ 0 at all.

V. SM PREDICTIONS
We now perform a fit of Eq. (6) to the combined experimental and theoretical constraints described above in Sections III and IV. We find that Eq. (6) with K = 2 provides an excellent fit to all inputs, with a p-value of 0.91. All 1D-marginalised posteriors are reasonably symmetric around their modes. The result of this fit is a set of correlated values for the complex parameters α (λ) k , which are summarized in Table I. These values lead to a determination of the non-local correlator in Eq. (2) that is consistent with the B → K * ψ n measurements, the theory calculations at negative q 2 , and it is independent of new physics in semileptonic operators. Thus, unlike Ref. [33], this is a genuine SM determination.
The gray band in Fig. 1 shows the result of this "prior" fit for the case of the real part of H ⊥ (q 2 ). Similar plots for the other correlators are provided in the appendix for completeness.
With these results at hand, we can compute SM pre- dictions for all observables of interest within the range 0 ≤ q 2 14 GeV 2 . One of them is the angular observable P 5 [34], which is the visible face of the "B → K * µ + µ − anomaly" [35]. Our SM prediction for P 5 is represented by the gray band in Fig. 2. We find relatively small uncertainties and a clearly apparent tension with LHCb data (represented by purple boxes in Fig. 2).
Another interesting SM prediction that we obtain from our analysis is: in agreement with the world average [36]. The larger uncertainties as compared to Ref. [37] are due to our doubling of the form factor uncertainties. SM predictions for all other observables will be given elsewhere.

VI. NEW PHYSICS ANALYSIS
We now perform a fit to B → K * µ + µ − data using as prior information the SM predictions derived in Section V. We include the branching ratio and the angular observables S i [38] within the q 2 bins in the region 1 ≤ q 2 14 GeV 2 . We use the latest LHCb measurements [39,40], and perform different separate fits, using the results from the maximum-likelihood fit excluding (LLH) and including (LLH2) the inter-resonance bin, or using the results from the method of moments [41] (MOM and MOM2), and both including (NP fit) and not including (SM fit) a floating NP contribution to C 9 .
The fits provide posterior distributions for the correlator, for B → K * µ + µ − and B → K * γ observables, and for C 9 . We first discuss some illustrative results of the LLH2 fit. The posteriors for the real part of H ⊥ (q 2 ) are shown in Fig. 1, both for the SM and the NP fits. In this case it is reassuring that both are consistent within errors with the result of the prior fit, indicating that modifying the long-distance contribution does not lead to improvement in the SM fit, and so the long-distance contribution is not likely to mimic a NP contribution.
The posterior NP prediction for P 5 (corresponding to the LLH2 fit) is shown in Fig. 2, exhibiting a much better agreement with the experimental measurements than the SM (prior) prediction.
The main conclusion of the fits is the following. The SM fits are relatively inefficient in comparison with the NP fits, with posterior odds [42] ranging from ∼ 2.7 to ∼ 10 (on the log scale) in favor of the NP hypothesis. The one-dimensional marginalized posteriors yield: The corresponding pulls with respect to the SM point C SM 9 (µ = 4.2 GeV) = 4.27 range from 3.4 to 6.1 standard deviations, and are illustrated in Fig. 3. These results, from a fit to B → K * µ + µ − data only, are in qualitative agreement with global fits [42][43][44][45][46][47][48], but rely on a more fundamented theory treatment.

VII. CONCLUSIONS
Analyticity provides strong constraints on the hadronic contribution to B → K * observables, and fixes the q 2 dependence up to a polynomial, which under some circumstances is an expansion in a small kinematical parameter. In this letter we have exploited this idea to propose a systematic approach to determine the non-local contributions, which at this time are the main source of theory uncertainty. This approach is systematically improvable with more precise data on B → K * ψ n and/or more precise theory calculations at negative q 2 . In addition, this approach allows access to the inter-resonance region, which provides valuable information on shortdistance physics. We have focused on B → K * , but the approach applies to any other B → M modes such as B → {K, π, ρ} and B s → φ .
We have performed a numerical analysis implementing this idea, and conclude that significantly improved theory predictions can be obtained, leading to a more precise and robust interpretation of experimental data and an improved sensitivity to short-distance physics. We thus believe that this approach will become very useful in fu- The effective Lagrangian that governs b → s transitions contains the following terms relevant for our analysis: with the current-current operators defined as such that C 2 (µ = M W ) = 1 + O(α s ), with P L(R) = (1∓γ 5 )/2 and T A the generators of SU(3), and the dipole and semileptonic operators given by  TABLE II. Uncorrelated priors for the CKM parameters in our analysis, taken from the tree-level-only fit in Ref. [49].
In the SM, the values of the Wilson coefficients are: The Lorentz decomposition for the form factors that we use in this analysis is given by denoting Γ ≡ K * (k, η)|Γ|B(q + k) . The non-local correlator H αµ is decomposed analogously according to Eq. (4): The Lorentz structures S αµ λ are given by: The experimental constraints in Section III are based on the experimental pseudo-observables in Eq. (9), which are obtained by fit to B → K * ψ n data. For B → K * J/ψ this data includes the branching ratio as measured by Belle [29], as well as the full set of angular observables F ⊥ , F , δ ⊥ and δ measured by BaBar [26] and LHCb [30]. For B → K * ψ(2S) the data includes the branching ratio  and the longitudinal polarization measured by Belle [28], and the full set of angular observables from BaBar [26]. For all measurements, correlations have been taken into account where available. More recent results for the full angular distributions, stemming from amplitude analyses that take into account tetra-quark contributions [28,29], are not used here. The ansatz involving tetra-quark amplitudes is incompatible with the basis of our analysis. Although we expect to be able to use these additional results in future studies, this requires further dedicated work.
The relevant input to this fit are the CKM parameters, listed in Table II, and the form factors (see Eq. (10)). Since the experimental inputs are sensitive to the form factors at q 2 ≤ 10 GeV 2 , we use the combined fit to K *meson LCSR and Lattice results performed in Ref. [6]. However, we double the uncertainties quoted in [6] to ensure full agreement among the LCSRs and Lattice results. Note that we do not account for correlations among the r ψn λ due to correlations among the form factor parameters. The results for the pseudo-observables are given in Table III. The two sets of observables for J/ψ and ψ(2S) are correlated with correlation matrices: In both cases, the mean and standard deviations have been obtained from a fit to 10 6 samples of the posterior predictive distributions. On the other hand, the correlation coefficients have been obtained from the sample covariance of these 10 6 samples. It is noteworthy that none of the coefficients exceeds a level of 78% for the J/ψ and 68% for the ψ(2S), respectively. The theory constraints in Section IV are based on pseudo-observables at four different points at spacelike q 2 . The derived values including uncertainties and correlations are listed in Table IV.
Nominally, for K = 2 the fit would involve 18 realvalued parameters α (λ) k . Using the property that the longitudinal correlator must vanish at zero momentum transfer q 2 = 0,Ĥ 0 (z(q 2 = 0)) = 0, we can reduce the number of parameters by 2 through the replacement From the fit of the parameters α (λ) k to the theory constraints and the pseudo-observables r ψn λ we find all 1D-marginalized posteriors to be reasonably symmetric around their modes. As for the pseudo-observables, we obtain means and standard deviations from fits to 10 6 samples of the PDF, while the correlation coefficients are obtained through computation of the sample covariance. Our results are summarized in Table I. Finally, the set of predictive distributions for the ratios H λ /F λ are shown in Fig. 5 for all transversities, including real and imaginary parts.
A key result of our analysis is that modifications to the correlator parameters cannot bring the theory predictions and measurements into better agreement. However, modifications to the form factor parameters can reduce the tensions, albeit not completely. We show in Fig. 4(a) the impact on the form factor V (q 2 ), which is the one most affected. Moreover, since the form factor A 1 is almost unaffected, we find that the SM fit prefers a substantial violation of the large-recoil symmetry relation [50,51] involving V and A 1 , as shown in Fig. 4(b).
[  A-priori predictions of the ratios H λ /F λ as functions of z, and comparison with two a-posteriori results. Note that the SM fit posterior is multi-modal, and we therefore choose not to illustrate its best-fit curves.