Model-independent observation of exotic contributions to $B^0\to J/\psi K^+\pi^-$ decays

An angular analysis of $B^0\to J/\psi K^+\pi^-$ decays is performed, using proton-proton collision data corresponding to an integrated luminosity of 3~${\rm fb}^{-1}$ collected with the LHCb detector. The $m(K^+\pi^-)$ spectrum is divided into fine bins. In each $m(K^+\pi^-)$ bin, the hypothesis that the three-dimensional angular distribution can be described by structures induced only by $K^\ast$ resonances is examined, making minimal assumptions about the $K^+\pi^-$ system. The data reject the $K^\ast$-only hypothesis with a large significance, implying the observation of exotic contributions in a model-independent fashion. Inspection of the $m(J/\psi\pi^-)$ versus $m(K^+\pi^-)$ plane suggests structures near $m(J/\psi\pi^-)=4200$~MeV and 4600~MeV.

In the Standard Model, the quark model allows for hadrons comprising any number of valence quarks, as long as they are colour-singlet states.Yet, after decades of searches, the reason why the vast majority of hadrons are built out of only quark-antiquark (meson) or three-quark (baryon) combinations remains a mystery.The best known exception is the Z(4430) − resonance with spin-parity 1 − and width Γ = 172 ± 13 MeV [1,2], 1 which has minimal quark content ccū d, and is therefore manifestly exotic, i.e., has components that are neither quark-antiquark or three-quark combinations.The only confirmed decay of the Z(4430) − state is via Z(4430) − → ψ(2S)π − , as seen in B 0 → ψ(2S)K + π − decays [1]. 2he corresponding Z(4430) − → J/ψ π − decay rate is suppressed by at least a factor of ten [3].The authors of Ref. [4] surmise that in a dynamical diquark picture, this owes to a better overlap of the Z(4430) − radial wavefunction with the excited state ψ(2S) than with the ground state J/ψ .For the B 0 → J/ψ K + π − channel, the Belle collaboration [3] has reported the observation of a new exotic Z(4200) − resonance decaying to J/ψ π − , that might correspond to the structure in m(ψ(2S)π − ) seen in Ref. [1] at around the same mass.
A generic concern in searches for broad exotic states like the Z(4430) − resonance is disentangling contributions from non-exotic components.For B 0 → ψ ( ) K + π − decays, 3the latter comprise different K * J resonances with spin J, that decay to K + π − .Figure 1 shows the K * J spectrum, which has multiple, overlapping, and poorly measured states.The bulk of the measurements come from the LASS K + π − scattering experiment [5].In particular, the decay B 0 → J/ψ K + π − is known to be dominated by K * J resonances, with an exotic fit fraction of only 2.4% [3], compared to a 10.3% contribution from the Z(4430) − for B 0 → ψ(2S)K + π − [6].This smaller exotic fit fraction for the J/ψ case makes it pertinent to study the evidence of exotic contributions in a manner independent of the dominant but poorly understood K * J spectrum.The BaBar collaboration [8] has performed a model-independent analysis of B 0 → ψ ( ) K + π − decays making minimal assumptions about the K * J spectrum, using two-dimensional (2D) moments in the variables m(K + π − ) and the K + helicity angle.The key feature of this approach is that no information on the exact content of the K * J states, including their masses, widths and m(K + π − )-dependent lineshapes, is required.An amplitude analysis would require the accurate description of the K * J lineshapes which depend on the underlying production dynamics.The model-independent procedure bypasses these problems, requiring only knowledge of the highest spin, J max , among all the contributing K * J states, for a given m(K + π − ) bin.Within uncertainties, the m(J/ψ π − ) spectrum in the BaBar data was found to be adequately described using just K * J states, without the need for exotic contributions.
In this Letter, a four-dimensional (4D) angular analysis of B 0 → J/ψ K + π − decays with J/ψ → µ + µ − is reported, employing the Run 1 LHCb dataset.The data sample corresponds to a signal yield approximately 40 and 20 times larger than those of the corresponding BaBar [8] and Belle [6] analyses, respectively.The larger sample size allows analysis of the differential rate as a function of the four variables, m(K + π − ), θ V , θ l and χ, that fully describe the decay topology.The lepton helicity angle, θ l , and the azimuthal angle, χ, between the (µ + µ − ) and (K + π − ) decay planes, were integrated over in the BaBar 2D analysis [8].The present 4D analysis therefore benefits from a significantly better sensitivity to exotic components than the previous 2D analysis.
The LHCb detector is a single-arm forward spectrometer covering the pseudorapidity range 2 < η < 5 and is described in detail in Ref. [9].Samples of simulated events are used to obtain the detector efficiency and optimise the selection.The pp collisions are generated using Pythia [10] with a specific LHCb configuration [11].Decays of hadronic particles are described by EvtGen [12], in which final-state radiation is generated using Photos [13].Dedicated control samples are employed to calibrate the simulation for agreement with the data.
The selection procedure is the same as in Refs.[14,15] for the rare decay B 0 → µ + µ − K + π − , with the additional requirement that the m(µ + µ − ) mass is constrained to the known J/ψ mass via a kinematic fit [16].The data sample is divided into 35 fine bins in m(K + π − ) such that the m(K + π − )-dependence can be neglected inside a given bin, and each subsample is processed independently.The bin-widths vary depending on the data sample size in a given m(K + π − ) region.Backgrounds from B + → J/ψ K + , B 0 s → J/ψ K + K − and Λ 0 b → J/ψ pK − decays are reduced to a level below 1% of the signal yield at the selection stage using the excellent tracking and particle-identification capabilities of the LHCb detector, and are subsequently removed by a background subtraction procedure.The B 0 (s) → J/ψ K + π − signal lineshape in the m(J/ψ K + π − ) spectrum is described by a bifurcated Gaussian core and exponential tails on both sides.A sum of two such lineshapes is used for the signal template for the mass fit, while the background lineshape is a falling exponential.The exponential tails in the signal lineshape are fixed from the simulation and all other parameters are allowed to vary in the fit, performed as a binned χ 2 minimisation.An example mass fit result is given in the Appendix.The cumulative signal yield in the m(K + π − ) ∈ [745, 1545] MeV region is 554,500 ± 800.
The strategy in this analysis is to examine the hypothesis that non-exotic K * J contributions alone can explain all features of the data.Under the approximation that the muon mass can be neglected and within a narrow m(K + π − ) bin, the CP -averaged transition matrix element squared is [17,18] where H η,J λ are the K * J helicity amplitudes and d j m ,m are Wigner rotation matrix elements.The helicities of the outgoing lepton and K * J are η = ±1 and λ ∈ {0, ±1}, respectively.Parity conservation in the electromagnetic J/ψ → µ − µ + decay leads to the relation The differential decay rate of B 0 → J/ψ (→ µ + µ − )K + π − with the K + π − system including spin-J partial waves with J ≤ J k max can be written as where the angular part in Eq. 1 has been expanded in an orthonormal basis of angular functions, f i (Ω).Here, k enumerates the m(K + π − ) bin under consideration and dΩ = dcos θ dcos θ V dχ is the angular phase space differential element.The angular basis functions, f i (Ω), are constructed from spherical harmonics, Y m l ≡ Y m l (θ l , χ), and reduced spherical harmonics, P m l ≡ √ 2πY m l (θ V , 0), and are given in are given in the Appendix.The Γ k i moments are observables that have an overall m(K + π − ) dependence, but within a narrow m(K + π − ) bin, this dependence can be neglected.The number of moments for the k th bin, n k max , depends on the allowed spin of the highest partial wave, J k max , and is given by [18] n Thus, for spin 3 onward, each additional higher spin component leads to 12 additional moments.In contrast to previous analyses, dcos θ dχ is not integrated over, which would have resulted in integrating over 10 out of these 12 moments, for each additional spin.Due to the orthonormality of the f i (Ω) basis functions, the angular observables, Γ k i , can be determined from the data in an unbiased fashion using a simple counting measurement [17].For the k th m(K + π − ) bin, the background-subtracted raw moments are estimated as where Ω p refers to the set of angles for a given event in this m(K + π − ) bin.The corresponding covariance matrix is Here, n k sig and n k bkg correspond to the number of candidates in the signal and background regions, respectively.The signal region is defined within ±15 MeV of the known B 0 mass, and the background region spans the range m(J/ψ K + π − ) ∈ [5450, 5560] MeV.The scale factor, x k , is the ratio of the estimated number of background candidates in the signal region divided by the number of candidates in the background region and is used to normalise the background subtraction.
To unfold effects from the detector efficiency including event reconstruction and selection, an efficiency matrix, E k ij , is used.It is obtained from simulated signal events generated according to a phase space distribution, uniform in Ω, as The w k p weight factors correct for differences between data and simulation, and the summation is over simulated and reconstructed events.They are derived using the B 0 → J/ψ K * (892) 0 control mode, as described in Refs.[14,15].The efficiency-corrected moments and covariance matrices are estimated as The first moment, Γ k 1 , corresponds to the overall rate.The remaining moments and the covariance matrix are normalised to this overall rate as Γ for i, j ∈ {2, . . ., n k max }.The normalisation with respect to the total rate renders the analysis insensitive to any overall systematic effect not correlated with dΩ in a given m(K + π − ) bin.The uncertainty from limited knowledge of the background is included in the second term in Eq. 5.The effect on the normalised moments, Γ k i , due to the uncertainty in the x k scale factors from the mass fit, is found to be negligible.The effect due to the limited simulation sample size compared to the data is small and accounted for using pseudoexperiments.The last source of systematic uncertainty is the effect of finite resolution in the reconstructed angles.The estimated biases in the measured Γ k i moments are added as additional uncertainties.The dominant contributions to B 0 → J/ψ K + π − are from the K * (892) 0 and K * 2 (1430) 0 states.To maximise the sensitivity to any exotic component, the dominant K * (892) 0 region that serves as a background for any non-K * J component, the analysis is performed on the m(K + π − ) ∈ [1085, 1445] MeV region, as marked by the dot-dashed lines in Fig. 1.The value of J k max depends on m(K + π − ), with higher spin states suppressed at lower m(K + π − ) values, due to the orbital angular momentum barrier factor [19].As seen from Fig. 1, only states with spin J = {0, 1} contribute below m(K + π − ) ∼ 1300 MeV and spin J = {0, 1, 2} below m(K + π − ) ∼ 1600 MeV.As a conservative choice, J k max is taken to be one unit larger than these expectations, Any exotic component in the J/ψ π − or J/ψ K + system will reflect onto the entire basis of K * J partial waves and give rise to nonzero contributions from P l (cos θ V ) components for l larger than those needed to account for K * J resonances.From the completeness of the f i (Ω) basis, a model with large enough J k max also describes any exotic component in the data.For a given value of m(K + π − ), there is a one-to-one correspondence between cos θ V and the variables m(J/ψ π − ) or m(J/ψ K + ).Therefore a complete basis of P l (cos θ V ) partial waves also describes any arbitrary shape in m(J/ψ π − ) or m(J/ψ K + ), for a given m(K + π − ) bin.The series is truncated at a value large enough to describe the relevant features of the distribution in data, but not so large that it follows bin-by-bin statistical fluctuations.A value of J k max = 15 is found to be suitable.For the k th m(K + π − ) bin, the probability density function (pdf) for the J k max model is Simulated events generated uniformly in Ω, after incorporating detector efficiency effects and weighting by the pdf in Eq. 11, are expected to match the background-subtracted data.The background subtraction is performed using the sP lot technique [20], where the weights are determined from fits to the invariant m(J/ψ K + π − ) distributions described previously.Figure 2 shows this comparison between the background-subtracted data and weighted simulated events in the m(K + π − ) ∈ [1085, 1265] MeV region.The J k max = 2 model clearly misses the peaking structures in the data around m(J/ψ π − ) = 4200 MeV and 4600 MeV.This inability of the J k max = 2 model to describe the data, even though the first spin 2 state, K * 2 (1430) 0 , lies beyond this mass region, strongly points toward the presence of exotic components.These could be four-quark bound states, meson molecules, or possibly dynamically generated features such as cusps.
To obtain a numerical estimate of the significance of exotic states, the likelihood ratio test is employed between the null hypothesis (K * J -only, from Eq. 10) and the exotic hypothesis (J k max = 15) pdfs, denoted P k K * J and P k exotic , respectively.The test statistic used in the likelihood ratio test is defined as for the k th m(K + π − ) bin, where (Ω) denotes the 3-dimensional angular detector efficiency in this bin, derived from the simulation weighted to match the data in the B 0 production kinematics.The last term in Eq. 12 ensures normalization of the relevant pdf and is calculated from simulated events that pass the reconstruction and selection criteria From Eq. 3 the number of degrees-of-freedom (ndf) increases by 12 for each additional spin-J wave in each m(K + π − ) bin.From Eq. 10, for the J k max = 2 and 3 choices, ∆ndf = 12 × (15 − 2) = 156 and 12 × (15 − 3) = 144, respectively, between the exotic and K * J -only pdf's for each m(K + π − ) bin.Each additional degree-of-freedom between the exotic and K * J -only pdf adds approximately one unit to the computed ∆(−2 log L) in the data due to increased sensitivity to the statistical fluctuations, and ∆(−2 log L) is therefore not expected to be zero even if there is no exotic contribution in the data.The expected ∆(−2 log L) distribution in the absence of exotic activity is evaluated using a large number of pseudoexperiments.For each m(K + π − ) bin, 11,000 pseudoexperiments are generated according to the K * J -only model with the moments varied according to the covariance matrix.The number of signal and background events for each pseudoexperiment are taken to be those measured in the data.The detector efficiency obtained from simulation is parameterised in 4D.Each pseudoexperiment is analyzed in exactly the same way as the data, where an independent efficiency matrix is generated for each pseudoexperiment.This accounts for the limited sample size of the simulation for the efficiency unfolding.The pseudoexperiments therefore represent the data faithfully at every step of the processing.
Figure 3 shows the distribution of ∆(−2 log L) from the pseudoexperiments in the m(K + π − ) ∈ [1085, 1445] MeV region comprising six m(K + π − ) bins each with the J k max = 2  3, shows a deviation of more than 10σ from the null hypothesis, corresponding to the distribution of the pseudoexperiments.The uncertainty due to the quality of the Gaussian profile fit in Fig. 3 is found to be negligible.The choice of large J k max for P k exotic , as well as the detector efficiency and calibration of the simulation, are systematically varied in pseudoexperiments, with significance for exotic components in excess of 6σ observed in each case.
In summary, employing the Run 1 LHCb dataset, non-K * J contributions in B 0 → J/ψ K + π − are observed with overwhelming significance.Compared to the previous BaBar analysis [8] of the same channel, the current study benefits from a 40-fold increase in signal yield and a full angular analysis of the decay topology.The method relies on a novel orthonormal angular moments expansion and, aside from a conservative limit on the highest allowed K * J spin for a given m(K + π − ) invariant mass, makes no other assumption about the K + π − system.Figure 4 shows a scatter plot of m(J/ψ π − ) against m(K + π − ) in the background-subtracted data.While the model-independent analysis performed here does not identify the origin of the non-K * J contributions, structures are visible at m(J/ψ π − ) ≈ 4200 MeV, close to the exotic state reported previously by Belle [3], and at m(J/ψ π − ) ≈ 4600 MeV.To interpret these structures as exotic tetraquark resonances and measure their properties will require a future model-dependent amplitude analysis of  The four kinematic variables for the process B 0 → J/ψ (→ µ + µ − )K + π − are the invariant mass m(K + π − ), and the three angles {θ l , θ V , χ}.The angle conventions for B 0 and B 0 are depicted in Fig. 5. Assuming negligible direct CP violation and production asymmetry, the rate expression remains the same between the charge-conjugate modes.C. Further comparison between the J k max = 2 and J k max = 15 models

D. Angular moments definitions
The transversity basis amplitudes, H J { ,⊥} , are defined as and the amplitudes for spin J ∈ {0, 1, 2} are denoted as S, H {0, ,⊥} and D {0, ,⊥} , respectively.For K * J contributions up to J = 2, there are 28 angular moments from the expansion Eq. 1, as explicitly listed in Table 1 in terms of the transversity amplitudes.The addition of K * J states from spin-3 onward results in 12 moments for each additional spin.The form of the moments are listed in Table 2, leading to the expression appearing in Eq. 3 of the main text.Further details can be obtained from Refs.[17,18].Table 1: The transversity-basis moments of the 28 orthonormal angular functions f i (Ω) in Eq. 2 till spin-2 in the K + π − system.

Figure 1 :
Figure 1: Spectrum of K *J resonances from Ref.[7], with the vertical span of the boxes indicating ±Γ 0 , where Γ 0 is the width of each resonance.The horizontal dashed lines mark the m(K + π − ) physical region for B 0 → J/ψ K + π − decays, while the dot-dashed lines mark the specific region, m(K + π − ) ∈ [1085, 1445] MeV, employed for determining the significance of exotic contributions.

Figure 2 :
Figure 2: Comparison of m(Jψπ − ) in the m(K + π − ) ∈ [1085, 1265] MeV region between the background-subtracted data and simulated events weighted by moments models with J k max = 2 and J k max = 15.

Figure 3 :
Figure 3: Likelihood-ratio test for exotic significance.The data shows a 10σ deviation from the pseudoexperiments generated according to the null hypothesis (K * J -only contributions).

Figure 7
Figure7shows a comparison between the J k max = 2 and J k max = 15 moments models in the m(K + π − ) = 895 ± 15 MeV bin.Since the spin-1 K * (892) 0 resonance strongly dominates here, the two models are compatible.

Figure 8
Figure 8 shows a comparison between the J k max = 2 and J k max = 15 moments models in the m(K + π − ) ∈ [1265, 1445] MeV.

Figure 7 :
Figure 7: Comparison of m(J/ψπ − ) in the m(K + π − ) ∈ [880, 910] MeV region between the background-subtracted data and simulation data weighted by moments models with J k max = 2 and J k max = 15.

Figure 8 :
Figure 8: Comparison of m(J/ψπ − ) in the m(K + π − ) ∈ [1265, 1445] MeV region between the background-subtracted data and simulation data weighted by moments models with J k max = 2 and J k max = 15.

Table 2 :
The 12 angular terms for each additional spin-J wave in the K + π − system, for J ≥ 3. Università di Roma La Sapienza, Roma, Italy l AGH -University of Science and Technology, Faculty of Computer Science, Electronics and Telecommunications, Kraków, Poland m LIFAELS, La Salle, Universitat Ramon Llull, Barcelona, Spain n Hanoi University of Science, Hanoi, Vietnam o Università di Padova, Padova, Italy p Università di Pisa, Pisa, Italy q Università degli Studi di Milano, Milano, Italy r Università di Urbino, Urbino, Italy s Università della Basilicata, Potenza, Italy t Scuola Normale Superiore, Pisa, Italy u Università di Modena e Reggio Emilia, Modena, Italy v H.H. Wills Physics Laboratory, University of Bristol, Bristol, United Kingdom w MSU -Iligan Institute of Technology (MSU-IIT), Iligan, Philippines x Novosibirsk State University, Novosibirsk, Russia y Sezione INFN di Trieste, Trieste, Italy z School of Physics and Information Technology, Shaanxi Normal University (SNNU), Xi'an, China aa Physics and Micro Electronic College, Hunan University, Changsha City, China ab Lanzhou University, Lanzhou, China † Deceased a Universidade Federal do Triângulo Mineiro (UFTM), Uberaba-MG, Brazil b Laboratoire Leprince-Ringuet, Palaiseau, France c P.N.Lebedev Physical Institute, Russian Academy of Science (LPI RAS), Moscow, Russia d Università di Bari, Bari, Italy e Università di Bologna, Bologna, Italy f Università di Cagliari, Cagliari, Italy g Università di Ferrara, Ferrara, Italy h Università di Genova, Genova, Italy i Università di Milano Bicocca, Milano, Italy j Università di Roma Tor Vergata, Roma, Italy k