CP asymmetry from hadronic charm rescattering in $B^\pm \to \pi^-\pi^+\pi^\pm$ decays at the high mass region

A model for the $B^\pm \to \pi^-\pi^+\pi^\pm$ decay amplitude is proposed to study the large CP violation observed at the high mass region of the Dalitz plane. A short distance $ b \to u $ amplitude with the weak phase $\gamma$ is considered together with the contribution of a hadronic charm loop and a s-wave $D\bar{D}\to \pi\pi$ rescattering. In the model, the $\chi_c^0$ appears as a narrow resonant state of the $D\bar D$ system below threshold. It is introduced in an unitary two channel S-matrix model of the coupled $D\bar D$ and $\pi\pi$ channels, where the $\chi_c^0$ complex pole in $D\bar D$ channel shows its signature in the off-diagonal matrix element and in the associated $D\bar{D}\to \pi\pi$ transition amplitude. The strong phase of the resulting decay amplitude has a sharp sign change at the $D\bar D$ threshold, changing the sign of the CP asymmetry, as it is observed in the data. We conclude that the hadronic charm loop and rescattering mechanism are relevant to the broadening of the CP asymmetry around the $\chi_c^0$ resonance in the $\pi\pi$ channel. This novel mechanism provides a possible interpretation of the CP asymmetry defier experimental result presented by the LHCb collaboration for the $B^\pm \to \pi^-\pi^+\pi^\pm$ decay in the high mass region.

B ± → K ± π − π + with B ± → K − K + K ± and also B ± → K − K + π ± with B ± → π − π + π ± [1, 2, 3, 5, 6]. These two possible sources of CP violation were already confirmed for B ± → π − π + π ± [7,8] and B ± → K − K + π ± decay channels [9] through a recent amplitude analysis performed by the LHCb collaboration for the Run I data. However, there are also strong experimental evidences for CPV in the Dalitz phase-space along the high mass region in all of those charged charmless B three-body decays [3]. Although the source for this CPV is not yet identified, we can assume that the variation of the CP asymmetry in the Dalitz plane is originated by the running of the strong phase along the phase-space, since the dominant weak phase contributing to these decays is the CKM phase γ, which must be a constant.
Recently we proposed a new source of strong phase variation, associated with the possible DD → K + K − rescattering, which couples the B ± → DDK ± to B ± → K − K + K ± and the B ± c → DDπ ± to B + c → K − K + π + decay channels [10,11]. Where in the later we have also considered the contribution of the channel B ± c → DD s K ± through DD s → Kπ rescattering. In these studies, we concluded that the long distance hadronic loop originated by the double charm penguin contribution can produce a strong phase that changes along the Dalitz phase-space. The phase starts around -1 radian until the DD threshold, then it has a quick phase variation given by a sharp change from negative to positive values. This phase variation can be responsible to change the sign of the CP asymmetry, as observed in experimental data [3] 2 , if the associated amplitude is interfering with another one carrying the weak phase. This CP asymmetry sign change at high mass is much more apparent in the Mirandizing plot of the B ± → π − π + π ± decay [3] near the DD open channel.
More than two decades ago it was predicted a CP violation in the high mass region of the B ± → π − π + π ± decay phase-space due to the presence of χ 0 c resonance [12]. Produced from b → ccd transition at tree level, without weak phase, it can interfere with b → uūd tree diagram amplitude, with weak phase, leading to a strong CP asymmetry in the phase-space, with the possibility to extract the CKM weak phase γ [12,13]. One expects the χ 0 c would be finally observed soon with the Run II LHCb data. This conclusion is based on counting the number of events already seen in the Cabibbo allowed B ± → K ± χ 0 c decay, in the B ± → K − K + K ± and B ± → K ± π − π + decays. Amplitude analysis performed by the Babar collaboration found a fit fraction of 1% for these three-body final states [14] and [15], respectively. From that one can do a simple relation with these decays and the Cabibbo suppressed B ± → π ± χ 0 c ( sin 2 θ cos 2 θ ≈ 0.05) to estimate the number of events expected in LHCb Run II for the B ± → π − π + π ± decay, arriving up to a few hundred events involving this scalar charmonium resonance.
Although LHCb did not find yet contribution from the B ± → π ± χ 0 c amplitude in B ± → π − π + π ± decay [7,8], the Mirandizing distribution for Run I data [3] have shown already a clear and huge CP asymmetry around the χ 0 c 2 see LHCb-PAPER-2014-044 supplemental material at https://cds.cern.ch/record/1751517/files/. invariant mass. This asymmetry suggests the presence of this resonance interfering with the nonresonant amplitude placed in this region 3 . However, the distribution of CP asymmetry is much larger than the narrow width expected for this resonance, suggesting that part of the nonresonant background around the χ 0 c peak comes from the same physical process that produces this resonance. Also, it is observed a change of sign in the CP asymmetry around the DD threshold, that can be assigned to the amplitude proposed in [10].
The discussion about the importance of charm loops as a source of important contribution in heavy decay processes is not new [16,17,18,19,20]. In particular, Colangelo et al. [18,19] calculated the B − → K − χ c0 decay rate using a hadronic triangle loop in combination with QCD factorization (QCDF) approach and HMχT to describe the heavy-light mesons vertices, including the coupling between DD → cc resonances. They argue that only QCDF cannot predict correctly the experimental branching fractions of the B → K(cc) transition, and in this framework the B − → K − χ c0 process is not allowed. Indeed other models of the B ± → π − π + π ± decay amplitude were proposed in the literature using QCDF approach and none of them included χ c0 [21,22,23,24].
In this work, we explore the same mechanism used to describe the B ± → K − K + K ± decays [10] (also applied to the rare B + c → K + K − π + decay [11]), namely the hadronic charm loop and DD rescattering to the two light pseudoscalars, to investigate the B ± → π − π + π ± decay, in an attempt to extract the main qualitative features observed in the high mass region (M 2 ππ > 3 GeV 2 ) of the CPV Mirandizing distribution. The present study brings one important novelty to the S-matrix model with respect to the one used to describe the s-wave scattering in the coupled DD → K + K − channels and transition amplitude [11]. The S-matrix model is generalized to the s-wave scattering in the coupled channels DD − ππ and in addition it includes the χ 0 c resonance with mass (3414.7 ± 0.3) MeV and width (10.5 ± 0.8) MeV [25], suggested to be a tetraquark [26]. Furthermore, focusing on a mechanism that can generate CP violation in high mass regions, the hadronic charm loop with rescattering is added to a constant amplitude carrying the weak phase, as will be explained and fully explored in what follows.
Decay amplitude model. A CPV process has to be described by a decay amplitude that must have two interfering contributions carrying different strong and weak phases in order to be observed in charge conjugate channels. The standard mechanism at the quark level to produce CP asymmetry in charmless charged B decays, comes through the interference of tree and penguin amplitudes from the short distance physics proposed in BSS model [27]. In the case of B ± → π − π + π ± , we assume that the weak phase γ come from the tree level diagram with quark transition b → u. For simplicity we neglect the penguin contribution b → d to the direct B ± → π − π + π ± decay process. The hadronic decay channel having as source tree or loop diagrams at the partonic level can also contribute with a strong phase from the final state interaction, and also from low energy resonances. On the other side, the B decay in two charmed mesons have a hadronic penguin like topology, that together with the subsequent rescattering process is assumed to contribute with a strong phase.
Inspired by the isobar model description of three-body decays, the amplitude of B ± → π − π + π ± decay can be parametrised by two independent contributions as: where we assume that A DD amplitude is dominated by a charm hadronic loop, Fig. 1, and A ± tree which is the dominant topology, has weak (±γ) and strong phases. Furthermore, the χ 0 c will be introduced as a resonant state below threshold within the DD scattering amplitude. Note, however, that our model is a naive representation of the nature and does not meant to be complete. The determination of γ from the experimental data of three-body decay is much more complex than our naive model suggests, being far from trivial separate the b → u and b → d amplitudes in this case. Despite of that, we will exploit the model to find out the fingerprint of the DD rescattering, with χ 0 c as a resonant molecular-like state below the DD threshold, looking at the CP violation distribution of the B ± → π − π + π ± decay in the high mass region.
A remark on the implication of the CPT invariance to the CP asymmetry of the B ± → π − π + π ± decay to the present model is appropriate. In the framework developed by Wolfenstein [16] (see also [28]) where the hadronic final-state interaction and the CPT constraint were considered together, the CP asymmetry seen in channels coupled by the long distance QCD dynamics are related. The consequence of this framework is that the sum of the partial widths for channels coupled by the strong Hamiltonian should be identical to the sum in the charge conjugated channels. Such result is more restrictive than the general CPT condition that gives equal lifetime for a particle and its anti-particle. Therefore, the CP asymmetry induced by rescattering, i.e., the compound contribution [29], should be considered together with the usual CP-violating amplitude from the BSS mechanism at the partonic level [27], in a way that short-and long-range dynamics have their place in producing the observed asymmetries.
The Wolfenstein formalism was further elaborated in [30]. It was considered the hadronic transition matrix of different channels coupled by the strong interaction in leading order in the expansion of the CP violating heavy meson decay amplitude. This amplitude was shown, despite the approximation, to fulfill the CPT constraint. Restricted to two channels the leading order formalism was applied to study the CP asymmetries seen in the B ± → K − K + K ± and B ± → K ± π − π + in the mass region where the K + K − and π + π − channels are strongly coupled. It explained the remarkable opposite signs and the shape of the projected CP asymmetry as a function of the invariant mass of the pair. Furthermore, this shape notably resembled the corresponding one of the magnitude of the KK → ππ transition amplitude as a function of the invariant mass of the pair. This mechanism was confirmed by the LHCb collaboration amplitude analyses for B ± → K − K + π ± [9] which found 65% of asymmetry due to KK → ππ with a different sign of the one observed in B + → π + π + π − decays [8,7], although with less intensity.
We observe that the leading order formalism applies to the present model of the three-body B decay where the B ± → DDπ ± and B ± → π − π + π ± channels are coupled by the strong force and the associated DD and ππ S-matrix provides the final state interaction contribution to the decay amplitude. The CP asymmetry of the B ± → DDπ ± has to receive a corresponding contribution with opposite sign respecting CPT invariance if only this channel coupling is present. However, the DD channel can also coupled to KK as we already discussed in [11], suggesting that the CP asymmetry in B ± → DDπ ± would call for contributions tfrom final state interaction involving more hadronic channels, a discussion that is much beyond the scope of the present work.
Hadronic charm loop. The charm rescattering contribution to the B ± → π − π + π ± decay can be described by a triangle loop of D mesons as the source for DD, which makes a transition to π + π − , as one can see in Fig. 1, for two possible charge states as depicted in the diagram. In this case, because both possibilities are similar we consider only the neutral one, B + → D 0D0 π + , which is similar to our previous study of the B ± → K − K + K ± decay [10].
The technique to compute the triangle loop given in Fig. 1 was already developed in previous works within the context of hadronic three-body decays [10,11,31,32]. For the sake of clarity, we repeat some of the steps required to formulate and compute the loop integral.
We assume factorisation to built the B → DDπ vertex in the loop diagram of Fig. 1 , which is written as the product: where G F is the Fermi constant and V cb and V * cd the matrix elements (m.e.) of the CKM matrix. The currents m.e.'s are described by form factors with the single pole approximation, and for convenience we introduce the notation The former represents the vector current m.e. of the B + →D 0 transition and the latter takes into account the amplitude for the pair D 0 π produced from the W boson excitation from the vacuum through the vector resonance D * . From crossing one finds that V µ Dπ represents the vector current m.e. of the D → π transition obtained from the vector meson dominance (VMD), with a single D * contribution.
The vector current m.e. for the transition B + →D 0 is written as: where = P B − pD 0 is the momentum transfer and the vector and scalar form factors are, for simplicity, given by: where ∆ B * (k 2 ) = k 2 − m 2 B * + i , which follows from describing the form factor as suggested by VMD, with the coupling of the weak current to B * , namely the heavy meson vector ground state with mass m B * .
The vector current m.e. that represents the amplitude for the pair D 0 π produced from the vacuum through the resonance D * is parametrized as: where the form factor is F + pole . The product of both currents m.e.'s in Eq.
(2) is written in terms of where p π is the bachelor pion momentum and the contraction of the Lorentz structure of the m.e.'s is given by the invariant: where p D 0 = − p π , pD0 = P B − and s = (P B − p π ) 2 is the mass of the pion pair in the transition DD → ππ.
The full amplitude for the hadronic loop is obtained by integrating the momentum inside the triangle with off-shell propagators, taking into account the absorptive and dispersive part of the triangle. It has to include the vertex DD → ππ which will be discussed in the sequence. The integral expression is given by: with The scattering amplitude, T DD→ππ (s) acts on the s-wave and does not contribute to the invariant, Eq. (7), from the scalar product of the m.e.'s of the vector currents. Furthermore, we assume minimal unitarity to describe the Tmatrix in the coupled DD and ππ channels, with its m.e.'s depending only on the Mandelstam variable s that allows to factorize out T DD→ππ (s) from the loop integral. The loop integral is calculated using the Feynman technique, which gives: where The functions I xyz are Feynman integrals defined in Appendix A, which are computed numerically with meson masses and widths from Ref. [25].
S-matrix and DD → ππ transition amplitude. We modify our previous phenomenological amplitude for DD → KK [11] and generalize it for DD → ππ. Furthermore, χ 0 c is introduced as resonant state below the DD threshold. This is an improvement with respect to the previous approach and different from considering only the contribution of χ 0 c to the DD → ππ transition as a Breit-Wigner resonance. Generically, a unitary two channel S-matrix can be parametrized as where δ 1 and δ 2 are the phase-shifts of the ππ and DD elastic scattering. η is the inelasticity parameter, which accounts for the probability flux between the two channels. Unitarity demands that the off-diagonal S-matrix elements should have a magnitude lower than one, and its modulus square can be interpreted as the probability to occur the transition between the initial and final channels. We introduce a parametrization for the phase-shifts and inelasticity parameter based on the reasonings presented in [10,11,33], brought to the context of DD → ππ transition to estimate T DD→ππ (s), which is a key ingredient to the hadronic charm loop to form the pions in the final state. Of course one should, in principle, resort to the QCD theory to compute the S-matrix, which, is however, much beyond our work. Given the importance of the DD → ππ rescattering to understand the CPV violation pattern in the high mass region, we retrieve the considerations given in previous works [10,11,33] generalized to the channels of our present interest.
A proposal for the dependence of the transition probability with the twomeson invariant mass, s, in light-meson processes has been discussed in [33], in the context of P V → P X transitions, and here, these qualitative reasonings are brought in the light of the present case. In a naive description of the π + π − inelastic collision amplitude, the pions annihilate into a quark-antiquark pair that propagates before recombining to produce the heavy-meson pair. The intermediate virtual state propagation of the quark pair scales roughly with the inverse of Mandelstam invariant s. The breakup of the pion into a quarkantiquark pair brings another factor of s −1 , and similarly for the formation of the D meson for s >> m 2 c , with m c the charm quark mass. That provides a damping factor of the off-diagonal S-matrix element of ∼ s −3 , which combined with the threshold behaviour gives √ s − s th /s 2.5 , keeps the asymptotic form for s large. Therefore, we write: where the normalization factor N is chosen to keep the modulus of the S-matrix elements smaller than 1, as required by the unitarity constraint. If we chose N = Λ 6 = (1.24) 6 in Eq.(12) then the maximum value reaches ∼ 0.87, at √ s = 1.08 √ s th , which is close to example of the s−wave isospin zero ππ → KK, where the cross section drops fast and is relevant below √ s ∼ 1.6 GeV [34]. This qualitative formula is also consistent with one of the possible parametrizations for inelasticity parameter η(s) = 1 − |S ππ→KK (s)| 2 presented in Ref. [35].
The magnitude of the off-diagonal S-matrix element is then written as Eq. (12), which is valid for s > s thDD . However, the three-body phase-space for the B decay has two pions below the DD threshold, which makes necessary the analytic continuation for s < s thDD in the physical sheet of complex momentum, imposing that k 2 → iκ 2 with κ 2 = 1 2 √ s thDD − s. Furthermore, the amplitude has to be regulated at low values of s. One phenomenological way to introduce an infrared cutoff in Eq. (12) is: where s QCD is an infrared scale of QCD estimated to be s QCD ∼ Λ 2 QCD ∼ 0.2 GeV 2 . A factor s α in the non-physical region, expressing that the coupling between the open channel and the off-mass-shell DD pair is damped in the nonphysical region, but without changing the large momentum power-law of the amplitude.
Next, we discuss the parametrization of the elastic phase-shift in the π + π − channel that takes the form dictated by the effective range expansion: where k 1 = 1 2 √ s − s thππ with the respective threshold of s thππ = 4M 2 π , the parameters b = 1 GeV −1 and c = 0.2 GeV come from our previous study [11].
The new aspect of the S-matrix parametrization with respect to our previous work is the introduction of the χ 0 c as a resonant molecular-like state of the DD system below threshold, which contributes as a pole in the diagonal term and to the phase δ 2 : where k 2 = 1 2 √ s − s thDD , and the DD threshold is s thDD = (M D + MD) 2 . For DD channel, we choose a complex scattering length dominated parametrization. We define the scattering length such that the elastic DD amplitude presents a pole in the complex plane. The real part below threshold accounts for the χ c mass and the width moves the pole into the complex plane: For s < s thDD the transition amplitude becomes where with N a normalization factor. The above formula respects the unitarity of the S-matrix model.
For s ≥ s thDD the transition amplitude is written as: where s thDD 2s−s thDD β was introduced to modulate the shape of the amplitude bump above the DD threshold as we have already used in the study of the DD → KK [11]. We should observe that our naive power counting can have corrections, and indeed this is the case as it will be show in our numerical study.
In our naive modeling we left as free parameters in Eqs. (17) and (18), the exponents α and β which can be determined by a fit to the data. As a theoretical exercise, we compare the transition amplitude obtained for the same set of parameters found in the study of B + c → K − K + π + (model I: α = 7 and β = 2) [11], and vary α and β to find another set (model II: α = 4 and β = 0.5), which seems more suitable to provide a qualitative description of the experimental data for CPV in the B ± → π − π + π ± decay for the high mass region. As shown in Fig. 2, the amplitude from Eqs. (17) and (18) plotted as a function of the ππ invariant mass, can have quantitative different signatures depending of the choice of the two exponents, but keeping three common features: (i) the χ 0 c peak superposed to a wide bump below DD threshold; (ii) the zero at the threshold; (iii) a bump above the threshold; and (iv) a jump of the strong phase close to π when crossing the DD threshold. The parameters can only move the quantitative values of the transition amplitude magnitude, while keeping the qualitative features (i)-(iv). The phase is not affected by the particular choice of parameters α and β once it is connected to the dynamical choice of the amplitude. We just remind the reader that the B ± → π − π + π ± decay amplitude includes the loop integral. Results for B ± → π − π + π ± decay and CPV. The total amplitude model for the B ± → π − π + π ± decay, Eq. (1), is the sum of a tree amplitude and the hadronic charm loop with DD → ππ rescattering. The CP asymmetry of B ± decay in the three-body phase-space will be the result of the interference between the tree amplitude A ± tree and A DD . In what follows, we are only interested in the dynamics above 3 GeV 2 where the low mass resonances contributions come mainly from their tails. Therefore, the amplitude A ± tree can be approximated as a flat nonresonant (NR) amplitude with a constant weak phase, γ, and an arbitrary strong phase: where a 0 can be complex to accommodate a strong phase. The total amplitude was simulated using Laura++ software [36] with hundred thousands events. There are two main variables when two amplitudes interfere: the relative phase between them and the relative magnitude, in principle those quantities are fixed by a fit to data. In our toy model we have to chose a 0 and in order to have an insight on the typical results one get by changing this quantity, we present a systematic study with model II.
To start our simulations, it is interesting to check the signature of each amplitude A ± tree and A DD alone in the phase-space projected on the m ππ high invariant mass which is defined as the higher one from the pair of the two invariants with opposite charge. Thus, we integrate in the m ππ low invariant mass starting at m 2 ππ =3 GeV 2 to exclude the low energy interaction region. In Fig. 3, one can see the result from the flat NR amplitude deformed by the phasespace integral and the hadronic loop with model II, where the χ 0 c is included in the transition amplitude. Each of them alone does not lead to CP violation, as expected.
In Fig. 4, we present the study of how the amplitudes interfere by fixing the relative magnitude for the NR to be twice the charm loop with rescattering amplitude and change the relative global phase between them. As one can see, the different relative phases can result in completely different patterns, but with a clear mark at the resonance position. In the bottom left frame in Fig. 4, the phase difference of 180 o eliminates the χ 0 c peak and make it appears as a dip. Whereas with 0 o phase the peak is enhanced. In principle, we have the freedom to chose the relative phase and intensity of the decay amplitudes A ± tree and A DD besides the model parameters, which can be fitted to data. However, our goal in this study, is to check if the model is able to reproduce the main characteristics observed in the LHCb data [3]: a CP asymmetry (A CP ) positive above 3 GeV 2 until the region where the charm channel opens and A CP flips sign. We can retrieve such A CP pattern with model II and a weak phase of γ = 70 o [25] by chosing, guided by the study presented in Fig. 4, the relative phase to be 45 o with magnitude of the NR amplitude twice the one for the hadronic charm loop with rescattering.
In the left frame of Fig. 5, we show that we can indeed produce the desired characteristics for A CP described above for the projection in the three-body phase-space. We also checked the CP violation signature produced by the interference of the same flat NR amplitude with a simple Breit-Wigner representation of the χ 0 c in an isobar model with the same relative phase and magnitude as above. We have found that the CP asymmetry is localized in a much smaller region around χ 0 c compared what we have observed with the rescattering model.  [4] applied to expose the CP violation in different regions of the B ± → π − π + π ± phase-space for model II.
In order to study the CPV signature between the B + and B − decays in the three-body phase-space in the high mass region, we use the Miranda technique [4] and present the CPV distribution in the right frame of Fig. 5, which can be compare to the LHCb data for the CP asymmetry in B ± → π − π + π ± decays in the Dalitz plane [3]. It is clear the signature of the χ 0 c as DD resonant molecular-type state below threshold with the peak widened by the DD → ππ rescattering. One can see in the right frame of Fig. 5 the red band for positive CP asymmetry in the χ 0 c region followed by a blue band pointing for a change in sign around the region of the DD threshold. A similar pattern can be identified in the experimental data [3]. We recall that there are other contributions that could spread the CP asymmetry of the B ± → π − π + π ± decay in high mass region, which were not considered here, like the tails of the low mass resonances, the excited states of the D's system, still coupled to ππ channels, and/or three-body rescattering in the DDπ channel.
Summary. We developed a model for the B ± → π − π + π ± decay amplitude, which has contribution from a hadronic charm loop with a s-wave DD → ππ rescattering, where the χ 0 c resonance is introduced as a resonant state of the DD system below threshold with the narrow experimental width. The χ 0 c pole of the elastic DD scattering amplitude modifies the DD → ππ transition amplitude due to the assumed S-matrix unitarity of the two-channel model. The narrow resonance appears in the ππ channel together with a wide distribution that spreads out the B ± decay intensity in a region of about 1 GeV 2 around the resonance. With this simple model for B ± → π − π + π ± decay amplitude we were able to mimic qualitatively the CP asymmetry distribution reported by LHCb Run I data in the high mass region [3], giving a possible interpretation of the mechanism behind these challenging experimental results. Therefore, we strongly encourage the experimentalists to incorporate the present model in their amplitude analyses for the next data generation in order to improve our understanding of the nature of CP violation in charmless three-body B decays in the high mass region.
The loop integral can be done using the standard Feynman technique: where the denominator is given by: For the specific case of B ± → π − π + π ± the four independent functions in Eq. (9), I D0D0 B * , I D0D0D * , ID 0D * B * and I D0D * B * , are obtained from the numerical integration of Eq. (A.2), with the denominators written explicitly as: and for the numerical integration we use a finite value of = 0.01 GeV.