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

TL;DR

This study performs a model-independent angular analysis of B0 decays to J/psi, K+, and pi-, revealing significant exotic contributions beyond known K* resonances, with potential new structures near 4200 and 4600 MeV.

Contribution

It introduces a minimal-assumption, model-independent method to detect exotic contributions in B0 decay angular distributions, challenging the K* resonance-only hypothesis.

Findings

Rejection of the K* resonance-only hypothesis with high significance.

Evidence of exotic contributions in B0 decay channels.

Identification of potential new structures near 4200 and 4600 MeV.

Abstract

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.

Figures

Tables

Table 1. Table 1: The transversity-basis moments of the 28 orthonormal angular functions f i ​ ( Ω ) subscript 𝑓 𝑖 Ω f_{i}(\Omega) in Eq.  2 till spin-2 in the K + ​ π − superscript 𝐾 superscript 𝜋 {{K}^{+}}{{\pi}^{-}} system.

$i$	$f_i\left(\Omega \right)$	${\Gamma }_i^{\mathrm{tr}}\left(q^2\right)$
1	$P_0^0{Y}_0^0$	$\left[{\left|H_0\right|}^2+{\left|H_{\parallel }\right|}^2+{\left|H_{⟂}\right|}^2+{\left|S\right|}^2+{\left|D_0\right|}^2+{\left|D_{\parallel }\right|}^2+{\left|D_{⟂}\right|}^2\right]$
2	$P_1^0{Y}_0^0$	$2\left[\frac{2}{\sqrt{5}}Re\left(H_0{D}_0^{\ast }\right)+Re\left(S{H}_0^{\ast }\right)+\sqrt{\frac{3}{5}}Re\left(H_{\parallel }{D}_{\parallel }^{\ast }+H_{⟂}{D}_{⟂}^{\ast }\right)\right]$
3	$P_2^0{Y}_0^0$	$\frac{\sqrt{5}}{7}$ (${\left|D_{\parallel }\right|}^2$ + ${\left|D_{⟂}\right|}^2$) - $\frac{1}{\sqrt{5}}$ (${\left|H_{\parallel }\right|}^2$ + ${\left|H_{⟂}\right|}^2$) + $\frac{2}{\sqrt{5}}$ ${\left|H_0\right|}^2$ + $\frac{10}{7\sqrt{5}}$ ${\left|D_0\right|}^2$ + $2$ $Re\left(S{D}_0^{\ast }\right)$
4	$P_3^0{Y}_0^0$	$\frac{6}{\sqrt{35}}\left[-Re\left(H_{\parallel }{D}_{\parallel }^{\ast }+H_{⟂}{D}_{⟂}^{\ast }\right)+\sqrt{3}Re\left(H_0{D}_0^{\ast }\right)\right]$
5	$P_4^0{Y}_0^0$	$\frac{2}{7}\left[-2\left({\left|D_{\parallel }\right|}^2+{\left|D_{⟂}\right|}^2\right)+3{\left|D_0\right|}^2\right]$
6	$P_0^0{Y}_2^0$	$\frac{1}{2\sqrt{5}}\left[\left({\left|D_{\parallel }\right|}^2+{\left|D_{⟂}\right|}^2\right)+\left({\left|H_{\parallel }\right|}^2+{\left|H_{⟂}\right|}^2\right)-2{\left|S\right|}^2-2{\left|D_0\right|}^2-2{\left|H_0\right|}^2\right]$
7	$P_1^0{Y}_2^0$	$\left[\frac{\sqrt{3}}{5}Re\left(H_{\parallel }{D}_{\parallel }^{\ast }+H_{⟂}{D}_{⟂}^{\ast }\right)-\frac{2}{\sqrt{5}}Re\left(S{H}_0^{\ast }\right)-\frac{4}{5}Re\left(H_0{D}_0^{\ast }\right)\right]$
8	$P_2^0{Y}_2^0$	$\left[\frac{1}{14}\left({\left|D_{\parallel }\right|}^2+{\left|D_{⟂}\right|}^2\right)-\frac{2}{7}{\left|D_0\right|}^2-\frac{1}{10}\left({\left|H_{\parallel }\right|}^2+{\left|H_{⟂}\right|}^2\right)-\frac{2}{5}{\left|H_0\right|}^2-\frac{2}{\sqrt{5}}Re\left(S{D}_0^{\ast }\right)\right]$
9	$P_3^0{Y}_2^0$	$-\frac{3}{5\sqrt{7}}\left[Re\left(H_{\parallel }{D}_{\parallel }^{\ast }+H_{⟂}{D}_{⟂}^{\ast }\right)+2\sqrt{3}Re\left(H_0{D}_0^{\ast }\right)\right]$
10	$P_4^0{Y}_2^0$	$-\frac{2}{7\sqrt{5}}\left[{\left|D_{\parallel }\right|}^2+{\left|D_{⟂}\right|}^2+3{\left|D_0\right|}^2\right]$
11	$P_1^1\sqrt{2}Re\left(Y_2^1\right)$	$-\frac{3}{\sqrt{10}}\left[\sqrt{\frac{2}{3}}Re\left(H_{\parallel }{S}^{\ast }\right)-\sqrt{\frac{2}{15}}Re\left(H_{\parallel }{D}_0^{\ast }\right)+\sqrt{\frac{2}{5}}Re\left(D_{\parallel }{H}_0^{\ast }\right)\right]$
12	$P_2^1\sqrt{2}Re\left(Y_2^1\right)$	$-\frac{3}{5}\left[Re\left(H_{\parallel }{H}_0^{\ast }\right)+\sqrt{\frac{5}{3}}Re\left(D_{\parallel }{S}^{\ast }\right)+\frac{5}{7\sqrt{3}}Re\left(D_{\parallel }{D}_0^{\ast }\right)\right]$
13	$P_3^1\sqrt{2}Re\left(Y_2^1\right)$	$-\frac{6}{5\sqrt{14}}\left[2Re\left(D_{\parallel }{H}_0^{\ast }\right)+\sqrt{3}Re\left(H_{\parallel }{D}_0^{\ast }\right)\right]$
14	$P_4^1\sqrt{2}Re\left(Y_2^1\right)$	$-\frac{6}{7\sqrt{2}}Re\left(D_{\parallel }{D}_0^{\ast }\right)$
15	$P_1^1\sqrt{2}Im\left(Y_2^1\right)$	$3\left[\frac{1}{\sqrt{15}}Im\left(H_{⟂}{S}^{\ast }\right)+\frac{1}{5}Im\left(D_{⟂}{H}_0^{\ast }\right)-\frac{1}{5\sqrt{3}}Im\left(H_{⟂}{D}_0^{\ast }\right)\right]$
16	$P_2^1\sqrt{2}Im\left(Y_2^1\right)$	$3\left[\frac{1}{7\sqrt{3}}Im\left(D_{⟂}{D}_0^{\ast }\right)+\frac{1}{5}Im\left(H_{⟂}{H}_0^{\ast }\right)+\frac{1}{\sqrt{15}}Im\left(D_{⟂}{S}^{\ast }\right)\right]$
17	$P_3^1\sqrt{2}Im\left(Y_2^1\right)$	$\frac{6}{5\sqrt{14}}\left[2Im\left(D_{⟂}{H}_0^{\ast }\right)+\sqrt{3}Im\left(H_{⟂}{D}_0^{\ast }\right)\right]$
18	$P_4^1\sqrt{2}Im\left(Y_2^1\right)$	$\frac{6}{7\sqrt{2}}Im\left(D_{⟂}{D}_0^{\ast }\right)$
19	$P_0^0\sqrt{2}Re\left(Y_2^2\right)$	$-\frac{3}{2\sqrt{15}}\left[\left({\left|H_{\parallel }\right|}^2-{\left|H_{⟂}\right|}^2\right)+\left({\left|D_{\parallel }\right|}^2-{\left|D_{⟂}\right|}^2\right)\right]$
20	$P_1^0\sqrt{2}Re\left(Y_2^2\right)$	$-\frac{3}{5}\left[Re\left(H_{\parallel }{D}_{\parallel }^{\ast }\right)-Re\left(D_{⟂}{H}_{⟂}^{\ast }\right)\right]$
21	$P_2^0\sqrt{2}Re\left(Y_2^2\right)$	$\frac{\sqrt{3}}{2}\left[-\frac{1}{7}\left({\left|D_{\parallel }\right|}^2-{\left|D_{⟂}\right|}^2\right)+\frac{1}{5}\left({\left|H_{\parallel }\right|}^2-{\left|H_{⟂}\right|}^2\right)\right]$
22	$P_3^0\sqrt{2}Re\left(Y_2^2\right)$	$\frac{3}{5}\sqrt{\frac{3}{7}}\left[Re\left(H_{\parallel }{D}_{\parallel }^{\ast }\right)-Re\left(D_{⟂}{H}_{⟂}^{\ast }\right)\right]$
23	$P_4^0\sqrt{2}Re\left(Y_2^2\right)$	$\frac{2}{7}\sqrt{\frac{3}{5}}\left({\left|D_{\parallel }\right|}^2-{\left|D_{⟂}\right|}^2\right)$
24	$P_0^0\sqrt{2}Im\left(Y_2^2\right)$	$\sqrt{\frac{3}{5}}\left[Im\left(H_{⟂}{H}_{\parallel }^{\ast }\right)+Im\left(D_{⟂}{D}_{\parallel }^{\ast }\right)\right]$
25	$P_1^0\sqrt{2}Im\left(Y_2^2\right)$	$\frac{3}{5}Im\left(H_{⟂}{D}_{\parallel }^{\ast }+D_{⟂}{H}_{\parallel }^{\ast }\right)$
26	$P_2^0\sqrt{2}Im\left(Y_2^2\right)$	$\sqrt{3}\left[\frac{1}{7}Im\left(D_{⟂}{D}_{\parallel }^{\ast }\right)-\frac{1}{5}Im\left(H_{⟂}{H}_{\parallel }^{\ast }\right)\right]$
27	$P_3^0\sqrt{2}Im\left(Y_2^2\right)$	$-\frac{3}{5}\sqrt{\frac{3}{7}}Im\left(D_{⟂}{H}_{\parallel }^{\ast }+H_{⟂}{D}_{\parallel }^{\ast }\right)$
28	$P_4^0\sqrt{2}Im\left(Y_2^2\right)$	$-\frac{4}{7}\sqrt{\frac{3}{5}}Im\left(D_{⟂}{D}_{\parallel }^{\ast }\right)$

Table 2. Table 2: The 12 angular terms for each additional spin- J 𝐽 J wave in the K + ​ π − superscript 𝐾 superscript 𝜋 {{K}^{+}}{{\pi}^{-}} system, for J ≥ 3 𝐽 3 J\geq 3 .

$i$	$f_i\left(\Omega \right)$
$1$	$P_{2J-1}^0{Y}_0^0$
$2$	$P_{2J}^0{Y}_0^0$
$3$	$P_{2J-1}^0{Y}_2^0$
$4$	$P_{2J}^0{Y}_2^0$
$5$	$P_{2J-1}^1\sqrt{2}Re\left(Y_2^1\right)$
$6$	$P_{2J}^1\sqrt{2}Re\left(Y_2^1\right)$
$7$	$P_{2J-1}^1\sqrt{2}Im\left(Y_2^1\right)$
$8$	$P_{2J}^1\sqrt{2}Im\left(Y_2^1\right)$
$9$	$P_{2J-1}^0\sqrt{2}Re\left(Y_2^2\right)$
$10$	$P_{2J}^0\sqrt{2}Re\left(Y_2^2\right)$
$11$	$P_{2J-1}^0\sqrt{2}Im\left(Y_2^2\right)$
$12$	$P_{2J}^0\sqrt{2}Im\left(Y_2^2\right)$

Equations

$$|\mathcal{M}{|}^2$$

$${\left[\frac{\mathrm{d}{\Gamma }^k}{\mathrm{d}{\Omega }^{}}\right]}_{J_{\max }^k}\propto \sum _{i=1}^{n_{\max }^k}{f}_i(\Omega ){\Gamma }_i^k,$$

$$n_{\max }^k=28+12\times (J_{\max }^k-2),\text{ for }{J}_{\max }^k>2.$$

$${\Gamma }_{i,\mathrm{raw}}^k=\sum _{p=1}^{n_{\mathrm{sig}}^k}{f}_i({\Omega }_p)-x^k\sum _{p=1}^{n_{\mathrm{bkg}}^k}{f}_i({\Omega }_p),$$

$$C_{ij,\mathrm{raw}}^k=\sum _{p=1}^{n_{\mathrm{sig}}^k}{f}_i({\Omega }_p)f_j({\Omega }_p)+(x^k{)}^2\sum _{p=1}^{n_{\mathrm{bkg}}^k}{f}_i({\Omega }_p)f_j({\Omega }_p).$$

$$E_{ij}^k=\sum _{p=1}^{n_{\mathrm{sim}}^k}{w}_p^k{f}_i({\Omega }_p)f_j({\Omega }_p).$$

$${\Gamma }_i^k$$

$$C_{ij}^k$$

$${\overline{C}}_{ij,\mathrm{stat}}^k$$

$$J_{\max }^k=\{\begin{array}{l}2\mathrm{for}1085\le m(K^{+}{\pi }^{-})<1265\text{MeV},\\ 3\mathrm{for}1265\le m(K^{+}{\pi }^{-})<1445\text{MeV}.\end{array}$$

$${\mathcal{P}}_{J_{\max }^k}(\Omega )$$

$$\displaystyle\Delta(-2\log{\cal L})\Big{|}_{k}\equiv-\displaystyle\phantom{x^{k}}\sum_{p=1}^{n^{k}_{\rm sig}}2\left[\log\left(\frac{{\cal P}^{k}_{K^{\ast}_{J}}(\Omega_{p})}{{\cal P}^{k}_{\rm exotic}(\Omega_{p})}\right)\right]+\displaystyle x^{k}\sum_{p=1}^{n^{k}_{\rm bkg}}2\left[\log\left(\frac{{\cal P}^{k}_{K^{\ast}_{J}}(\Omega_{p})}{{\cal P}^{k}_{\rm exotic}(\Omega_{p})}\right)\right]+$$

$$2\times (n_{\mathrm{sig}}^k-x^k{n}_{\mathrm{bkg}}^k)\times \log \left(\frac{\int {\mathcal{P}}_{K_J^{\ast }}^k(\Omega )ϵ(\Omega )\mathrm{d}\Omega }{\int {\mathcal{P}}_{\mathrm{exotic}}^k(\Omega )ϵ(\Omega )\mathrm{d}\Omega }\right),$$

$$E_i^k$$

$$\int {\mathcal{P}}_{J_{\max }^k}(\Omega )ϵ(\Omega )\mathrm{d}\Omega$$

$${\mathcal{H}}_{\pm }^J$$

Full text

EUROPEAN ORGANIZATION FOR NUCLEAR RESEARCH (CERN)

​​​

CERN-EP-2018-330

LHCb-PAPER-2018-043

**Model-independent observation of exotic contributions to **${{B}^{0}}\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}{{K}^{+}}{{\pi}^{-}}$decays

LHCb collaboration†††Authors are listed at the end of this paper.

An angular analysis of ${{B}^{0}}\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}{{K}^{+}}{{\pi}^{-}}$ decays is performed, using proton-proton collision data corresponding to an integrated luminosity of 3$\text{\,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\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}{{\pi}^{-}})$ versus $m({{K}^{+}}{{\pi}^{-}})$ plane suggests structures near $m({{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}{{\pi}^{-}})=4200$  MeV and 4600  MeV.

Published in Phys. Rev. Lett. 122, 152002 (2019)

© 2024 CERN for the benefit of the LHCb collaboration. CC-BY-4.0 licence.

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 $\Gamma=172\pm 13$ MeV [1, 2],111Natural units with $\hslash=c=1$ are used throughout the document. which has minimal quark content $c\bar{c}\bar{u}\bar{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)^{-}\rightarrow{\psi{(2S)}}\pi^{-}}$, as seen in ${B^{0}\rightarrow{\psi{(2S)}}K^{+}\pi^{-}}$ decays [1].222The inclusion of charge-conjugate decay modes is implied throughout. The corresponding ${Z(4430)^{-}\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}\pi^{-}}$ 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 is because of a larger overlap of the $Z(4430)^{-}$ radial wavefunction with the excited state $\psi{(2S)}$ than with the ground state ${J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}$. For the $B^{0}\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}K^{+}\pi^{-}$ channel, the Belle collaboration [3] has reported the observation of a new exotic $Z(4200)^{-}$ resonance decaying to ${{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}\pi^{-}$, that might correspond to the structure in $m({\psi{(2S)}}{{\pi}^{-}})$ 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}\rightarrow\psi^{(^{\prime})}K^{+}\pi^{-}$ decays,333Here $\psi$ denotes the ground state ${{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}$ and $\psi^{\prime}$ denotes the excited state $\psi{(2S)}$. the latter comprise different $K^{\ast}_{J}$ resonances with spin $J$, that decay to $K^{+}\pi^{-}$. Figure 1 shows the $K^{\ast}_{J}$ spectrum, which has multiple, overlapping, and poorly measured states. The bulk of the measurements come from the LASS $K^{+}\pi^{-}$ scattering experiment [5]. In particular, the decay $B^{0}\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}K^{+}{{\pi}^{-}}$ is known to be dominated by $K^{\ast}_{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}\rightarrow{\psi{(2S)}}K^{+}{{\pi}^{-}}$ [6]. This smaller exotic fit fraction for the ${J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}$ case makes it pertinent to study the evidence of exotic contributions in a manner independent of the dominant but poorly understood $K^{\ast}_{J}$ spectrum.

The BaBar collaboration [8] has performed a model-independent analysis of ${B^{0}\rightarrow\psi^{(^{\prime})}K^{+}\pi^{-}}$ decays making minimal assumptions about the $K^{\ast}_{J}$ spectrum, using two-dimensional (2D) moments in the variables $m({{K}^{+}}{{\pi}^{-}})$ and the $K^{+}$ helicity angle, $\theta_{V}$. The key feature of this approach is that no information on the exact content of the $K^{\ast}_{J}$ states, including their masses, widths and $m({{K}^{+}}{{\pi}^{-}})$-dependent lineshapes, is required. An amplitude analysis would require the accurate description of the $K^{\ast}_{J}$ lineshapes which depend on the underlying production dynamics. The model-independent procedure bypasses these problems, requiring only knowledge of the highest spin, $J_{\rm max}$, among all the contributing $K^{\ast}_{J}$ states, for a given $m({{K}^{+}}{{\pi}^{-}})$ bin. Within uncertainties, the $m({{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}\pi^{-})$ spectrum in the BaBar data was found to be adequately described using just $K^{\ast}_{J}$ states, without the need for exotic contributions.

In this Letter, a four-dimensional (4D) angular analysis of $B^{0}\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}K^{+}\pi^{-}$ decays with ${{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}\rightarrow{\mu^{+}}{\mu^{-}}$ 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}^{+}}{{\pi}^{-}})$, $\theta_{V}$, $\theta_{l}$ and $\chi$, that fully describe the decay topology. The lepton helicity angle, $\theta_{l}$, and the azimuthal angle, $\chi$, between the $({\mu^{+}}{\mu^{-}})$ and $(K^{+}\pi^{-})$ 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<\eta<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, *Sjostrand:2007gs] with a specific LHCb configuration [12]. Decays of hadronic particles are described by EvtGen [13], in which final-state radiation is generated using Photos [14]. Dedicated control samples are employed to calibrate the simulation for agreement with the data.

The selection procedure is the same as in Refs. [15, 16] for the rare decay ${B^{0}\rightarrow\mu^{+}\mu^{-}K^{+}\pi^{-}}$, with the additional requirement that the $m(\mu^{+}\mu^{-})$ mass is constrained to the known ${J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}$ mass via a kinematic fit [17]. The data sample is divided into 35 fine bins in $m({{K}^{+}}{{\pi}^{-}})$ such that the $m({{K}^{+}}{{\pi}^{-}})$-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}^{+}}{{\pi}^{-}})$ region. Backgrounds from $B^{+}\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}K^{+}$, $B^{0}_{s}\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}K^{+}K^{-}$ and ${{\mathchar 28931\relax}^{0}_{b}}\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}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)}\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}K^{+}\pi^{-}$ signal lineshape in the $m({{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}K^{+}\pi^{-})$ 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 $\chi^{2}$ minimisation. An example mass fit result is given in the Appendix. The cumulative signal yield in the $m({{K}^{+}}{{\pi}^{-}})\in[745,1545]$  MeV region is $554,\!500\pm 800$.

The strategy in this analysis is to examine the hypothesis that non-exotic $K^{\ast}_{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}^{+}}{{\pi}^{-}})$ bin, the $C\!P$-averaged transition matrix element squared is [18, 19]

[TABLE]

where $\mathcal{H}^{\eta,J}_{\lambda}$ are the $K^{\ast}_{J}$ helicity amplitudes and $d^{j}_{m^{\prime},m}$ are Wigner rotation matrix elements. The helicities of the outgoing lepton and $K^{\ast}_{J}$ are ${\eta=\pm 1}$ and ${\lambda\in\{0,\pm 1\}}$, respectively. Parity conservation in the electromagnetic ${{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}\rightarrow{\mu^{-}}{\mu^{+}}$ decay leads to the relation ${\mathcal{H}^{+,J}_{\lambda}=\mathcal{H}^{-,J}_{\lambda}\equiv\mathcal{H}^{J}_{\lambda}}$. The differential decay rate of ${{B}^{0}}\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}(\rightarrow{\mu^{+}}{\mu^{-}}){{K}^{+}}{{\pi}^{-}}$ with the ${K}^{+}$ ${\pi}^{-}$ system including spin-$J$ partial waves with $J\leq J^{k}_{\rm max}$ can be written as

[TABLE]

where the angular part in Eq. 1 has been expanded in an orthonormal basis of angular functions, $f_{i}(\Omega)$. Here, $k$ enumerates the $m({{K}^{+}}{{\pi}^{-}})$ bin under consideration and ${\mathrm{d}\Omega=\mathrm{d}\!\cos\theta_{\ell}\,\mathrm{d}\!\cos\theta_{V}\,\mathrm{d}\chi}$ is the angular phase space differential element. The angular basis functions, $f_{i}(\Omega)$, are constructed from spherical harmonics, $Y^{m}_{l}\equiv Y^{m}_{l}(\theta_{l},\chi)$, and reduced spherical harmonics, ${P^{m}_{l}\equiv\sqrt{2\pi}Y^{m}_{l}(\theta_{V},0)}$, and are given in are given in the Appendix.

The $\Gamma^{k}_{i}$ moments are observables that have an overall $m({{K}^{+}}{{\pi}^{-}})$ dependence, but within a narrow $m({{K}^{+}}{{\pi}^{-}})$ bin, this dependence can be neglected. The number of moments for the $k^{\rm th}$ bin, $n_{\rm max}^{k}$, depends on the allowed spin of the highest partial wave, $J_{\rm max}^{k}$, and is given by [19]

[TABLE]

Thus, for spin 3 onward, each additional higher spin component leads to 12 additional moments. In contrast to previous analyses, $\mathrm{d}\!\cos\theta_{\ell}\,\mathrm{d}\chi$ 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}(\Omega)$ basis functions, the angular observables, $\Gamma^{k}_{i}$, can be determined from the data in an unbiased fashion using a simple counting measurement [18]. For the $k^{\rm th}$ $m({{K}^{+}}{{\pi}^{-}})$ bin, the background-subtracted raw moments are estimated as

[TABLE]

where $\Omega_{p}$ refers to the set of angles for a given event in this $m({{K}^{+}}{{\pi}^{-}})$ bin. The corresponding covariance matrix is

[TABLE]

Here, $n^{k}_{\rm sig}$ and $n^{k}_{\rm bkg}$ correspond to the number of candidates in the signal and background regions, respectively. The signal region is defined within $\pm 15$  MeV of the known ${B}^{0}$ mass, and the background region spans the range $m({{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}{{K}^{+}}{{\pi}^{-}})\in[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 $\Omega$, as

[TABLE]

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}}\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}K^{\ast}(892)^{0}$ control mode, as described in Refs. [15, 16]. The efficiency-corrected moments and covariance matrices are estimated as

[TABLE]

The first moment, $\Gamma^{k}_{1}$, corresponds to the overall rate. The remaining moments and the covariance matrix are normalised to this overall rate as $\overline{\Gamma}^{k}_{i}\equiv\Gamma^{k}_{i}/\Gamma^{k}_{1}$ and

[TABLE]

for $i,j\in\{2,\ldots,n^{k}_{\rm max}\}$.

The normalisation with respect to the total rate renders the analysis insensitive to any overall systematic effect not correlated with $\mathrm{d}\Omega$ in a given $m({{K}^{+}}{{\pi}^{-}})$ bin. The uncertainty from limited knowledge of the background is included in the second term in Eq. 5. The effect on the normalised moments, $\overline{\Gamma}^{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 $\overline{\Gamma}^{k}_{i}$ moments are added as additional uncertainties.

The dominant contributions to ${{B}^{0}}\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}{{K}^{+}}{{\pi}^{-}}$ are from the $K^{\ast}(892)^{0}$ and $K^{\ast}_{2}(1430)^{0}$ states. To maximise the sensitivity to any exotic component, the dominant $K^{\ast}(892)^{0}$ region that serves as a background for any non-$K^{\ast}_{J}$ component, the analysis is performed on the ${m({{K}^{+}}{{\pi}^{-}})\in[1085,1445]}$  MeV region, as marked by the dot-dashed lines in Fig. 1. The value of $J^{k}_{\rm max}$ depends on $m({{K}^{+}}{{\pi}^{-}})$, with higher spin states suppressed at lower $m({{K}^{+}}{{\pi}^{-}})$ values, due to the orbital angular momentum barrier factor [20]. As seen from Fig. 1, only states with spin $J=\{0,1\}$ contribute below $m({{K}^{+}}{{\pi}^{-}})\sim 1300$  MeV and spin $J=\{0,1,2\}$ below $m({{K}^{+}}{{\pi}^{-}})\sim 1600$  MeV. As a conservative choice, $J^{k}_{\rm max}$ is taken to be one unit larger than these expectations,

[TABLE]

Any exotic component in the ${{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}{{\pi}^{-}}$ or ${{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}{{K}^{+}}$ system will reflect onto the entire basis of $K^{\ast}_{J}$ partial waves and give rise to nonzero contributions from $P_{l}(\cos\theta_{V})$ components for $l$ larger than those needed to account for $K^{\ast}_{J}$ resonances. From the completeness of the $f_{i}(\Omega)$ basis, a model with large enough $J^{k}_{\rm max}$ also describes any exotic component in the data. For a given value of $m({{K}^{+}}{{\pi}^{-}})$, there is a one-to-one correspondence between $\cos\theta_{V}$ and the variables $m({{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}{{\pi}^{-}})$ or $m({{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}{{K}^{+}})$. Therefore a complete basis of $P_{l}(\cos\theta_{V})$ partial waves also describes any arbitrary shape in $m({{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}{{\pi}^{-}})$ or $m({{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}{{K}^{+}})$, for a given $m({{K}^{+}}{{\pi}^{-}})$ 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}_{\rm max}=15$ is found to be suitable.

For the $k^{\rm th}$ $m({{K}^{+}}{{\pi}^{-}})$ bin, the probability density function (pdf) for the $J^{k}_{\rm max}$ model is

[TABLE]

Simulated events generated uniformly in $\Omega$, 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 $sPlot$ technique [21], where the weights are determined from fits to the invariant $m({{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}K^{+}\pi^{-})$ distributions described previously. Figure 2 shows this comparison between the background-subtracted data and weighted simulated events in the $m({{K}^{+}}{{\pi}^{-}})\in[1085,1265]$  MeV region. The $J^{k}_{\rm max}=2$ model clearly misses the peaking structures in the data around $m({{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}\pi^{-})=4200$  MeV and $4600$  MeV. This inability of the $J^{k}_{\rm max}=2$ model to describe the data, even though the first spin 2 state, $K^{\ast}_{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^{\ast}_{J}$-only, from Eq. 10) and the exotic hypothesis ($J^{k}_{\rm max}=15$) pdfs, denoted ${\cal P}^{k}_{K^{\ast}_{J}}$ and ${\cal P}^{k}_{\rm exotic}$, respectively. The test statistic used in the likelihood ratio test is defined as

[TABLE]

for the $k^{\rm th}$ $m({{K}^{+}}{{\pi}^{-}})$ bin, where $\epsilon(\Omega)$ 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. S0.Ex1 ensures normalization of the relevant pdf and is calculated from simulated events that pass the reconstruction and selection criteria

[TABLE]

Results from individual $m({{K}^{+}}{{\pi}^{-}})$ bins are combined to give the final test statistic {\Delta(-2\log{\cal L})=\displaystyle\sum_{k}\Delta(-2\log{\cal L})\Big{|}_{k}}.

From Eq. 3 the number of degrees-of-freedom (ndf) increases by 12 for each additional spin-$J$ wave in each $m({{K}^{+}}{{\pi}^{-}})$ bin. From Eq. 10, for the $J^{k}_{\rm max}=2$ and $3$ choices, $\Delta{\rm ndf}=12\times(15-2)=156$ and $12\times(15-3)=144$, respectively, between the exotic and $K^{\ast}_{J}$-only pdf’s for each $m({{K}^{+}}{{\pi}^{-}})$ bin. Each additional degree-of-freedom between the exotic and $K^{\ast}_{J}$-only pdf adds approximately one unit to the computed $\Delta(-2\log{\cal L})$ in the data due to increased sensitivity to the statistical fluctuations, and $\Delta(-2\log{\cal L})$ is therefore not expected to be zero even if there is no exotic contribution in the data. The expected $\Delta(-2\log{\cal L})$ distribution in the absence of exotic activity is evaluated using a large number of pseudoexperiments. For each $m({{K}^{+}}{{\pi}^{-}})$ bin, 11,000 pseudoexperiments are generated according to the $K^{\ast}_{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 $\Delta(-2\log{\cal L})$ from the pseudoexperiments in the $m({{K}^{+}}{{\pi}^{-}})\in[1085,1445]$  MeV region comprising six $m({{K}^{+}}{{\pi}^{-}})$ bins each with the $J^{k}_{\rm max}=2$ or $3$ choice. A fit to a Gaussian profile gives $\Delta(-2\log{\cal L})\approx 2051$ between the null and exotic hypothesis, even in the absence of any exotic contributions. This value is consistent with the naïve expectation $\Delta({\rm ndf})=1800$ from the counting discussed earlier. The value of $\Delta(-2\log{\cal L})$ for the data, as marked by the vertical line in Fig. 3, shows a deviation of more than $10\sigma$ 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}_{\rm max}$ for ${\cal P}^{k}_{\rm 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\sigma$ observed in each case.

In summary, employing the Run 1 LHCb dataset, non-$K^{\ast}_{J}$ contributions in ${B^{0}\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}K^{+}\pi^{-}}$ 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^{\ast}_{J}$ spin for a given $m({{K}^{+}}{{\pi}^{-}})$ invariant mass, makes no other assumption about the $K^{+}\pi^{-}$ system. Figure 4 shows a scatter plot of $m({{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}\pi^{-})$ against $m({{K}^{+}}{{\pi}^{-}})$ in the background-subtracted data. While the model-independent analysis performed here does not identify the origin of the non-$K^{\ast}_{J}$ contributions, structures are visible at ${m({{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}\pi^{-})\approx 4200}$  MeV, close to the exotic state reported previously by Belle [3], and at $m({{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}\pi^{-})\approx 4600$  MeV. To interpret these structures as exotic tetraquark resonances and measure their properties will require a future model-dependent amplitude analysis of the data.

Acknowledgements

We express our gratitude to our colleagues in the CERN accelerator departments for the excellent performance of the LHC. We thank the technical and administrative staff at the LHCb institutes. We acknowledge support from CERN and from the national agencies: CAPES, CNPq, FAPERJ and FINEP (Brazil); MOST and NSFC (China); CNRS/IN2P3 (France); BMBF, DFG and MPG (Germany); INFN (Italy); NWO (Netherlands); MNiSW and NCN (Poland); MEN/IFA (Romania); MSHE (Russia); MinECo (Spain); SNSF and SER (Switzerland); NASU (Ukraine); STFC (United Kingdom); NSF (USA). We acknowledge the computing resources that are provided by CERN, IN2P3 (France), KIT and DESY (Germany), INFN (Italy), SURF (Netherlands), PIC (Spain), GridPP (United Kingdom), RRCKI and Yandex LLC (Russia), CSCS (Switzerland), IFIN-HH (Romania), CBPF (Brazil), PL-GRID (Poland) and OSC (USA). We are indebted to the communities behind the multiple open-source software packages on which we depend. Individual groups or members have received support from AvH Foundation (Germany); EPLANET, Marie Skłodowska-Curie Actions and ERC (European Union); ANR, Labex P2IO and OCEVU, and Région Auvergne-Rhône-Alpes (France); Key Research Program of Frontier Sciences of CAS, CAS PIFI, and the Thousand Talents Program (China); RFBR, RSF and Yandex LLC (Russia); GVA, XuntaGal and GENCAT (Spain); the Royal Society and the Leverhulme Trust (United Kingdom); Laboratory Directed Research and Development program of LANL (USA).

**Appendix

**

**A. Angle conventions

**

The four kinematic variables for the process ${{B}^{0}}\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}(\rightarrow{\mu^{+}}{\mu^{-}}){{K}^{+}}{{\pi}^{-}}$ are the invariant mass $m({{K}^{+}}{{\pi}^{-}})$, and the three angles $\{\theta_{l},\theta_{V},\chi\}$. The angle conventions for ${{B}^{0}}$ and ${{\kern 1.79993pt\overline{\kern-1.79993ptB}}{}^{0}}$ are depicted in Fig. 5. Assuming negligible direct $C\!P$ violation and production asymmetry, the rate expression remains the same between the charge-conjugate modes.

**B. Example mass fit result in a particular bin

**

Figure 6 shows an example mass fit result for the $m({{K}^{+}}{{\pi}^{-}})\in[1235,1265]$  MeV bin.

**C. Further comparison between the $\bm{J^{k}_{\rm max}=2}$ and $\bm{J^{k}_{\rm max}=15}$ models

**

Figure 7 shows a comparison between the $J^{k}_{\rm max}=2$ and $J^{k}_{\rm max}=15$ moments models in the $m(K^{+}\pi^{-})=895\pm 15$ MeV bin. Since the spin-1 $K^{\ast}(892)^{0}$ resonance strongly dominates here, the two models are compatible.

Figure 8 shows a comparison between the $J^{k}_{\rm max}=2$ and $J^{k}_{\rm max}=15$ moments models in the $m({{K}^{+}}{{\pi}^{-}})\in[1265,1445]$  MeV.

**D. Angular moments definitions

**

The transversity basis amplitudes, $\mathcal{H}^{J}_{\{\parallel,\perp\}}$, are defined as

[TABLE]

and the amplitudes for spin $J\in\{0,1,2\}$ are denoted as $S$, $H_{\{0,\parallel,\perp\}}$ and $D_{\{0,\parallel,\perp\}}$, respectively. For $K^{\ast}_{J}$ contributions up to $J=2$, there are 28 angular moments from the expansion of Eq. 1, as explicitly listed in Table 1 in terms of the transversity amplitudes. The addition of $K^{\ast}_{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. [18, 19].