CP Violation in Multibody $B$ Decays from QCD Factorization

We test a data-driven approach based on QCD factorization for charmless three-body $B$-decays by confronting it to measurements of CP violation in $B^- \to \pi^- \pi^+ \pi^-$. While some of the needed non-perturbative objects can be directly extracted from data, some others can, so far, only be modelled. Although this approach is currently model dependent, we comment on the perspectives to reduce this model dependence. While our model naturally accommodates the gross features of the Dalitz distribution, it cannot quantitatively explain the details seen in the current experimental data on local CP asymmetries. We comment on possible refinements of our simple model and conclude by briefly discussing a possible extension of the model to large invariant masses, where large local CP asymmetries have been measured.


Introduction
Hadronic multi-body decays with more than two final-state hadrons constitute a large part of the branching fraction for heavy hadron non-leptonic decays. In principle, three-and more body decays have non-trivial kinematics and the phase space distributions contain far more information than the two-body decays.
For D decays there exists a large amount of data on multi-body decays. However, the charm quark mass is not large enough for heavy quark methods, since the typical invariant masses m ij of final state hadron pairs are roughly m c / √ N , where N is the final state multiplicity. Already for three-body decays of charmed hadrons this is outside the perturbative region and hence there is no chance to discuss the resulting amplitudes on the basis of some factorization theorem.
For B decays the situation may be slightly better, as we pointed out in a recent publication [1]. For three-body decays such as B → πππ the bottom-quark mass turns out to be still too small to allow for a complete factorization in the central region of the Dalitz plot. However, it seems that one can make use of a "partial factorization" at the edges of the Dalitz distribution, where the invariant mass of two of the pions is small.
It has been discussed in [1,2] that for this part of the phase space the same proof of factorization as for the two-body decays [3][4][5][6] is valid. However, the non-perturbative input given by the matrix elements of the factorized operators is different: in the case of three-body decays the light-cone distribution of two collinearly moving pions and the soft B → ππ form factor are needed. At least the edges of the Dalitz plot can be described in terms of these quantities; however, as has been argued in [1] this formalism may be extrapolated to the more central parts of the Dalitz plot.

JHEP10(2017)117
Hadronic multi-body decays are also interesting for studies of CP violation. Although the integrated (direct) CP asymmetries are small, local CP asymmetries (i.e. the CP asymmetry for fixed values of the final-state invariant masses) are measured to be large in some regions of the phase space and exhibit a rich structure [7][8][9][10][11][12]. Assuming the well known CKM mechanism for CP violation, the requirement for its appearance is an interference between at least two amplitudes with different weak and strong phases. Since the weak phases are independent of the kinematics, any dependence on the kinematical variables of the local CP asymmetries reflects kinematics-dependent strong-phase differences.
In the present paper we discuss an approach for three-body decay amplitudes based on QCD factorization. We will take into account the leading term only, which is equivalent to adopting naive factorization for the hadronic matrix elements. The main purpose is to study to which extent such a framework can properly describe the observed Dalitz distribution and local CP asymmetries in B − → π − π + π − .
In the next section we summarize the QCDF formula for three-body decays, then we discuss the non-perturbative input needed in the factorization formula. In section 4 we compute the Dalitz distributions and the local CP asymmetries in our framework, with a fit to experimental data. We conclude with a discussion of the results.
2 QCD-factorization for B − → π − π + π − In the following we will discuss charmless hadronic three-body decays and as a concrete example we consider B − → π − π + π − . We define the external momenta where p B = k 1 + k 2 + k 3 and, for massless pions, such that s 12 + s 13 + s 23 = 1. For B − → π − π + π − the Dalitz distribution is symmetric in s 12 and s 23 . Experimentally these variables cannot be distinguished, and we define k 1 and k 3 by s low ± ≡ s 12 and s high ± ≡ s 23 , with s low ± < s high ± . The application to other combinations of charges as well as to final states with kaons is obvious. As we have discussed in our previous paper [1], the structure of the amplitude in the region s low ± 1 is very similar to the two body case within the QCD-factorization framework. The only difference is that the matrix elements of the operators will eventually induce new non-perturbative quantities.
In this region, the B − → π − π + π − amplitude at leading order in α s and at leading twist is given by [

JHEP10(2017)117
The quantities λ p ≡ V pb V * pd encode the CKM factors, a 1,2 and a u,c 4 are constructed from Wilson coefficients, loop functions and convolutions with light-cone distributions (see section 4 and ref. [6]), and the objects f π , f + , F em π and F I t are non-perturbative quantities to be discussed in section 3. The amplitude in eq. (2.3) is the key formula in this paper.
At this order and twist, this formula coincides with the result obtained by applying the "naive-factorization" ansatz (see e.g. [13]), and we find it convenient to use some of its notation here. To this end, we can simply take the QCD-factorized effective Hamiltonian from [3], which reads: (2.6) The notation of the operators means that the matrix element is to be evaluated in the factorized form as a product of two matrix elements. To be specific, for the case at hand we have the two cases The relevant non-perturbative objects in the leading-order amplitude are the pion decay constant f π , the B → π form factor f + , the time-like helicity B → ππ form factors F I=0 t and F I=1 t , and the pion form factor in the time-like region F π . In the following section we give proper definitions for these objects and specify how they will be fixed in our approach.

Non-perturbative input
The strength of our QCD-factorization based model is that non-perturbative inputs may be obtained from data. The pion decay constant and the B → π form factors can both be taken as real, but the new B → ππ and pion form factors contain non-perturbative strong phases which will be driving the CP asymmetry distribution.  Figure 1. Pion vector form factor F em π (k 2 ) = |F em π |e iδ in the time-like region.

The pion decay constant and the timelike pion form factors
We define the pion decay constant in the usual way with the numerical value f π ∼ 130 MeV. The pion form factor is defined by and can be obtained from electromagnetic probes. Note that in the time-like region this form factor picks up a non-trivial strong phase. Here we use the parametrization of ref. [14] fitted to the measurements of e + e − → π + π − (γ) [15] (see also ref. [16]). The absolute value and the phase of this form factor are shown in figure 1. Unfortunately, while the absolute value is very precisely measured up to k 2 ∼ 3.5 GeV 2 , its phase is not so well constrained. This will add to the level of model dependence of our approach. We will also need the corresponding form factor for the scalar current: where the mass factors are chosen such that a proper chiral limit exists [17]. This form factor can be obtained using a coupled channel analysis. We use the results of ref. [18], which are valid up to around k 2 3 GeV 2 , as shown in figure 2. Similar results have been obtained in ref. [17] in connection with a study on B → J/ψππ. We note that the shape of F S π (k 2 ) around low-lying scalar resonances such as the f 0 (500) does not even remotely resemble the shape of a Breit-Wigner function.

The B → π form factor
We use the following definitions for the vector form factors [19]:  Figure 2. Pion scalar form factor F S π (k 2 ) = |F S π |e iδ S in the time-like region.
When applying the factorization formula (2.8), this expression is contracted with the time-like vector form factor for the two other pions. Using the fact that the current of these two pions is conserved, we get for (2.8) wherek ≡ k 1 − k 2 . For the form factor f + we use the LCSR calculation in ref. [20].

The B → ππ form factors
The form factors appearing in the B → ππ transitions have been studied in [21,22] for B → ππ ν and we use the definitions from these papers. However, when applying (2.7) we only need the contraction with the matrix element (3.1), and hence only a single form factor appears where we used that k 2 3 = m 2 π , and defines the polar angle θ π of the π − in the rest frame of the dipion, where β 2 is the Källén function. The two pions can have isospin I = 0 or I = 1, such that (3.8) The isovector form factor F I=1 t has been studied using QCD light-cone sum rules in [23,24], (see also [25] for a similar study of the other P -wave form factors). Analogous studies of the isoscalar form factor F I=0 t have not been performed. Here we model the form factor F I=1 t by assuming that the decay B → ππ proceeds only resonantly through JHEP10(2017)117 B → ρ → ππ. In this approximation we need the form factors for the B → ρ transition via the left-handed current. In general, this requires four form factors, but when applying (2.7) this reduces to a single form factor for the axial vector current where is the polarization vector of the ρ meson with momentum k, and q is the momentum transfer.
Treating the ρ as an intermediate resonance we obtain where we sum over the ρ polarizations and introduce the Breit-Wigner function where Γ P is the total decay width of the particle P . The decay matrix element for the ρ → ππ transition is defined as and g ρπ − π + can be obtained from the decay width of the ρ resonance.
Combining the various ingredients and contracting with q = k 3 gives Replacing the outgoing two pion state by a ρ resonance described by a simple Breit-Wigner shape is clearly a crude approximation for both the absolute value and the phase. We refine this approximation in the following way: we use the same model for the time-like form factor, and we determine a replacement for the Breit-Wigner function in terms of the measured pion form factor in figure 1. First, we have: The ρ-meson decay constant is then defined by

JHEP10(2017)117
which allows us to write the pion form factor as We can now solve for g ρππ B P and insert this into (3.13). Finally, using (3.6) yields where f ρ = 0.209 GeV and A 0 (m 2 π ) A 0 (0) = 0.36 ± 0.04 [26,27]. A similar procedure can be applied to the I = 0 channel, assuming dominance of a scalar resonance, B → S 0 → ππ. We write where the relevant part of the form factor for the B → S 0 transition is defined as (see e.g. [28]) with q = p B − k. Similarly, we can write F S π in the same way where the decay constant and the strong coupling constant are defined by Finally, we substitute again g Sπ − π + B S for F S π . However, F BS 0 and f S are unknown for the lightest scalar resonances. We thus model the isoscalar form factor through where the model parameters β and φ can be obtained from a fit to data. The approximations made to obtain the B → ππ form factors in terms of the pion form factors are currently unavoidable. In the future this modelling might be circumvented using QCD sum rules that employ the pion distribution amplitudes [24,25]. In addition, these models can be fitted separately to the light-cone sum rules with B distribution amplitudes, as done in ref. [23].

CP violation in
We start from the amplitude in eq. (2.3), use F I=1 t in (3.17) and we model F I=0 t using (3.21). Our model thus contains only two free parameters: β and the phase φ.
The CP asymmetry is given by

JHEP10(2017)117
whereĀ is equal to A with all weak phases conjugated. The required weak phase difference is given through the different structures with λ u = V ub V * ud = |V ud V ub |e −iγ , where γ is the corresponding weak phase of the Unitarity Triangle and λ c = V cb V cd = |V cb V cd | is real within our convention. We will use the values quoted in [29].
At tree level a u 4 = a c 4 and the coefficients a i are given by the Wilson coefficients C i [3, 5] where N C = 3 denotes the number of colors. At O(α s ) the coefficients a i also acquire perturbative strong phases [3,30,31]. These can be included using the partial QCD factorization formalism discussed in ref. [1], which requires taking into account the convolutions of the hard kernels with the generalized 2π distribution amplitude (DA) (see also [32][33][34]). At leading order, the pion DA and the generalized 2π DA reduce to their local limits, corresponding to eqs. (3.1) and (3.2), respectively. Since these O(α) correction cannot generate large CP asymmetries, we work at leading order, where the coefficients a i are real, leaving higher-order effects for future studies. The required strong phase difference to generate CP violation should thus come from the interference between the form factors F em π and F t . In our model, the phases of F em π and F I=1 t are identical. For elastic scattering (below the threshold of the first inelasticity in ππ scattering) this is a general statement following from Watson's theorem. This condition has been emphasized within the framework of QCD sum rules in ref. [23].
Inserting the amplitude in eq. (2.3) into (4.1), one finds that the CP asymmetry is proportional to A CP (s low ± , cos θ π ) = β sin γ sin(δ S (s low ± ) + φ − δ(s low ± )) cos θ π |F S π (s low ± )| |F em π (s low ± )| g(s low ± ) , (4.3) where g(s low ± ) is a real function that can be computed from eqs. (2.3) and (4.1). We have replaced the s high ± variable with cos θ π following eq. (3.7). We see that only the interference between F em π and F I=0 t terms contribute to the CP asymmetry. Therefore, the specific parametrizations for F em π and F I=0 t are of crucial importance. Here we use the parametrizations discussed in section 3 and depicted in figures 1 and 2, which allow us to perform a first analysis of our QCD-based model. The model dependence of our approach could be reduced in the future when more data for the form factors is available. We do not take into account uncertainties for the pion form factors.
Using the data from the LHCb Collaboration [7], we may fit our model parameters β and φ directly. Unfortunately, the full efficiency-and background-corrected Dalitz distribution is not provided by the LHCb analysis. Therefore, we use the projections of the data given for B + and B − decays, separated for cos θ π < 0 and cos θ π > 0. We show these two regions in the B − → π − π + π − Dalitz distribution in figure 3, as given in ref. [7],  . Dalitz distribution for B − → π − π + π − (a) as measured by the LHCb Collaboration [7] where the region below the line corresponds to cos θ π > 0 line as discussed in the text (b) Dalitz distribution of our model including (5.1).
where cos θ π > 0 corresponds to s high ± < 1 2 (1 − s low ± ), i.e. the lower part of the distribution. Figure 4 shows the projections of the LHCb data for cos θ π < 0 and cos θ π > 0 in bins of 0.05 GeV for the variable m 12 = k 1 + k 2 .
We now perform a fit to the data to determine the most likely values for the model parameters. The fit is performed by a standard χ 2 minimization. Our model predicts the decay rate for each bin; since the measurement of the absolute branching ratio is not available we have to scale our results to match the arbitrary units used in figure 4. Fitting this scaling parameter together with our parameters β and φ gives: The yield predictions with these best fit parameters are also included in figure 4. These figures show that our fit represents the data for B + at cos θ > 0 best, although in general our fit describes the data very poorly. We therefore refrain from giving an error to our fit parameters. These results call for refinements in the modelling of the form factors, which at this stage has been relatively simplistic. A number of possibilities will be mentioned later. A more clear picture of the situation is obtained by scrutinizing the CP asymmetry in more detail. In figure 6 we show the complete CP distribution as provided by LHCb [7] in a specific binning that ensures that each bin has the same number of events. The projections for the B − -B + yield differences are also given by LHCb [7]. We show these in figure 5 together with the outcome of our fit. The resulting CP violation in our model is much smaller than that seen in the data, nevertheless it reproduces the gross structures except for the region around 1.3 GeV.
In the region around m ρ we expect our model to most accurate. The differences as seen in the CP asymmetries might be due to the simplistic model used for F I=1 t in (3.17).
To study the effect of relaxing this assumption, we added another fit parameter to F I=I t .
Performing then the χ 2 analysis, leads to a slightly better agreement around the ρ peak,  but the total fit still remains a poor description of the data. The neglected higher-order terms might also give small modifications in this region. However, our model clearly fails to describe the interesting behavior of the CP asymmetry around 1.3 GeV. Here there is a positive CP asymmetry for both of the cos θ regions. In our model, the small CP violation in this region switches sign as does the CP asymmetry in the ρ region. This is because A CP in eq. (4.3) is only generated by a vector-scalar interference, which always comes with a cos θ π term, and hence the CP asymmetry always switches sign when comparing cos θ > 0 and cos θ < 0 (see also [35][36][37] for an elaborate discussion on these issues). However, if the CP asymmetry were dominated by two S or S-D wave interferences, the CP asymmetry would not switch sign, which could be an explanation for the behaviour in this region. Additional S or D-wave terms might still arise in our approach when including higher-order (twist) corrections, we leave the study of these corrections for future work. In addition, we note  Figure 5. Difference between the B − and B + yield in our best fit compared to the LHCb data [7]. that this region is also at the boundary of where the scalar form factor depicted in figure 2 can be trusted. Therefore, inelasticities may also play an important role.
For this first study, we have compared to the available projections of the LHCb data. Therefore, important information about the CP asymmetry in the s high +− variable is essentially washed out. Our simple model does not give a good quantitative description of the CP asymmetries, however, several refinements are possible. For future studies it would be beneficial and desirable to have the full information on the Dalitz and CP distributions.

Conclusion and comment on charm resonances
We have discussed a data-driven model based on QCD factorization to study CP violation in B → πππ, which depends on the model parameters β and φ (a strong phase). The form factors for the B → ππ transition as well as the time-like pion form factors have nonperturbative strong phases that lead to a complicated phase structure of the amplitudes T u and T c . Although we have shown that our simple model can describe some of the features of the decay rates and CP asymmetries, it cannot capture all the physics which is relevant for the local CP asymmetries. We have discussed some possible refinements of our model to accommodate these features. Nonetheless, beyond the particular model-dependent choices adopted in this analysis, the aim of being able to fix the amplitude in eq. (2.3) completely from data on the time-like pion form factors is an important one. In this way one can avoid the use of isobar assumptions [38,39] and Breit-Wigner-shaped resonance models (e.g. [40][41][42][43]). We are confident that progress will be made in this direction.
Since we use only the parametrizations of the scalar and vector form factors of the pion, the modelled non-perturbative phases are only the ones related to the final-state interactions of the two opposite-sign pions. This means that we can only expect this simple model to work within the regions where this is the dominant effect, i.e. at the corresponding edges of the phase space. Nevertheless, one might consider a possible extension of our  model, especially when considering the measured CP asymmetry in figure 6. These measurements find large local CP asymmetries at high s low +− and in regions of the phase space where there seem to be not many events when comparing with the Dalitz distribution in figure 3. Unfortunately, projections of this high momentum region are not (yet) available.
It is possible to extrapolate the pion form factors up to larger invariant mass. However, since there is no extra "structure" in this region this would fix the phases to around 180 • everywhere, suppressing the local CP asymmetry at high s low +− in contrast to the observation. The observed CP asymmetry might be created by subleading effects that were thus far assumed to be suppressed, but that might give significant effects at such high momenta. We note that the amplitudes T u and T c differ by the fact that T c contains penguins with charm, and it is thus sensitive to the heavy charm-quark mass. Subleading terms in QCD factorization for two-body decays generate perturbatively calculable strong phases for the coefficient a 4 which generates the CP violation in the B → ρπ decay. In the three-body decay one might expect a similar effect from the charm quarks, which would modify the shape of the local CP asymmetry. The details of this contribution will depend on the non-perturbative interaction of the two charm quarks in T c .
Clearly we do not have a way to actually compute this, so we have to make use of some modelling to get a qualitative picture. We first regard the region close to the charm threshold (2m c ) as the relevant region where sharp charmed resonances may affect the CP asymmetry. The simplest way to introduce non-trivial phases is to consider a resonancelike structure in T c described by a Breit-Wigner shape. Thus we modify our model by a simple addition: where T (0) q is the leading order amplitude given in eq. (2.3) by the term proportional to λ q . As an example, we fix the constant g to be 0.02 and Γ = 0.15 GeV, and take m c = 1.6 GeV for definiteness.
In figure 3 we show the resulting logarithmic Dalitz distribution, compared to the measured Dalitz distribution. Clearly the Dalitz plot is dominated by the resonance structure given by the time-like pion form factor, and the subleading term modelled by eq. (5.1) yields indeed small contributions as expected (and can be tuned with the constant g). However, as mentioned earlier it is not possible to qualitatively compare our Dalitz distribution with the measured one [7] because the latter is not background subtracted. A more thorough comparison in the line of that in section 4 would require the data projections for the high momentum part as well.
In figure 6 we also show the corresponding local CP asymmetry distribution. The subleading term (5.1) now generates a non-perturbative, phase-space dependent phase difference between T c and T u which induces sizeable local CP asymmetries in the region around the charm threshold. We note that the actual CP distribution depends on the values of g, Γ, β and φ. It might be interesting to include this charm-resonance model into an amplitude analysis to obtain further insights on the behaviour on the CP asymmetry at high momenta.
Obviously this is only a crude model. However, we note that the qualitative structure is in agreement with the observations by LHCb [7]. The data are not yet very precise, but the sizeable CP asymmetries observed by LHCb are compatible with structures originating from charm-threshold effects as we model them in eq. (5.1). However, it is extremely difficult to achieve a quantitative understanding of these effects from QCD. The issue of the charm contributions has also been a hot topic of discussion in two-body decays, but it is in three-body decays where there might be a chance to measure their effect and to interpret it. We believe this qualitative discussion may provide some motivation.