Amplitude analysis of $B^+\to J/\psi \phi K^+$ decays

The first full amplitude analysis of $B^+\to J/\psi \phi K^+$ with $J/\psi\to\mu^+\mu^-$, $\phi\to K^+K^-$ decays is performed with a data sample of 3 fb$^{-1}$ of $pp$ collision data collected at $\sqrt{s}=7$ and $8$ TeV with the LHCb detector. The data cannot be described by a model that contains only excited kaon states decaying into $\phi K^+$, and four $J/\psi\phi$ structures are observed, each with significance over $5$ standard deviations. The quantum numbers of these structures are determined with significance of at least $4$ standard deviations. The lightest has mass consistent with, but width much larger than, previous measurements of the claimed $X(4140)$ state. The model includes significant contributions from a number of expected kaon excitations, including the first observation of the $K^{*}(1680)^+\to\phi K^+$ transition.

The LHCb collaboration did not see evidence for the narrow X(4140) peak in the analysis presented in Ref. [21], based on a data sample corresponding to 0.37 fb −1 of integrated luminosity, a fraction of that now available. Searches for the narrow X(4140) did not confirm its presence in analyses performed by the Belle [22,23] (unpublished) and BaBar [24] experiments. The X(4140) structure was observed however by the CMS (5σ) [25] collaboration. Evidence for it was also reported in B + → J/ψ φK + decays by the D0 (3σ) [26] collaboration. The D0 collaboration claimed in addition a significant signal for prompt X(4140) production in pp collisions [27]. The BES-III collaboration did not find evidence for X(4140) → J/ψ φ in e + e − → γX(4140) and set upper limits on its production cross-section at √ s = 4.23, 4.26 and 4.36 GeV [28]. Previous results related to the X(4140) structure are summarized in Table 1. Table 1: Previous results related to the X(4140) → J/ψ φ mass peak, first observed in B + → J/ψ φK + decays. The first (second) significance quoted for Ref. [27] is for the prompt (nonprompt) production components. The statistical and systematic errors are added in quadrature and then used in the weights to calculate the averages, excluding unpublished results (shown in italics). The last column gives a fraction of the total B + → J/ψ φK + rate attributed to the X(4140) structure. In an unpublished update to their B + → J/ψ φK + analysis [29], the CDF collaboration presented 3.1σ evidence for a second relatively narrow J/ψ φ mass peak near 4274 MeV. This observation has also received attention in the literature [30,31]. A second J/ψ φ mass peak was observed by the CMS collaboration at a mass which is higher by 3.2 standard deviations, but the statistical significance of this structure was not determined [25]. The Belle collaboration saw 3.2σ evidence for a narrow J/ψ φ peak at 4350.6 +4. 6 −5.1 ± 0.7 MeV in two-photon collisions, which implies J P C = 0 ++ or 2 ++ , and found no evidence for X(4140) in the same analysis [32]. The experimental results related to J/ψ φ mass peaks heavier than X(4140) are summarized in Table 2. Belle [32] γγ → J/ψ φ 4350.6 +4.6 −5.1 ±0.7 13 +18 − 9 ±4 3.2σ In view of the considerable theoretical interest in possible exotic hadronic states decaying to J/ψ φ, it is important to clarify the rather confusing experimental situation concerning J/ψ φ mass structures. The data sample used in this work corresponds to an integrated luminosity of 3 fb −1 collected with the LHCb detector in pp collisions at center-of-mass energies 7 and 8 TeV. Thanks to the larger signal yield, corresponding to 4289 ± 151 reconstructed B + → J/ψ φK + decays, the roughly uniform efficiency and the relatively low background across the entire J/ψ φ mass range, this data sample offers the best sensitivity to date, not only to probe for the X(4140), X(4274) and other previously claimed structures, but also to inspect the high mass region.
All previous analyses were based on naive J/ψ φ mass (m J/ψ φ ) fits, with Breit-Wigner signal peaks on top of incoherent background described by ad-hoc functional shapes (e.g. three-body phase space distribution in B + → J/ψ φK + decays). While the m φK distribution has been observed to be smooth, several resonant contributions from kaon excitations (hereafter denoted generically as K * ) are expected. It is important to prove that any m J/ψ φ peaks are not merely reflections of these conventional resonances. If genuine J/ψ φ states are present, it is crucial to determine their quantum numbers to aid their theoretical interpretation. Both of these tasks call for a proper amplitude analysis of B + → J/ψ φK + decays, in which the observed m φK and m J/ψ φ masses are analyzed simultaneously with the distributions of decay angles, without which the resolution of different resonant contributions is difficult, if not impossible. The analysis of J/ψ and φ polarizations via their decays to µ + µ − and K + K − , respectively, increases greatly the sensitivity of the analysis as compared with the Dalitz plot analysis alone. In addition to the search for exotic hadrons, which includes X → J/ψ φ and Z + → J/ψ K + contributions, the amplitude analysis of our data offers unique insight into the spectroscopy of the poorly experimentally understood higher excitations of the kaon system, in their decays to a φK + final state.
In this article, an amplitude analysis of the decay B + → J/ψ φK + is presented for the first time, with additional results for, and containing more detailed description of, the work presented in Ref. [33].

LHCb detector
The LHCb detector [34,35] is a single-arm forward spectrometer covering the pseudorapidity range 2 < η < 5, designed for the study of particles containing b or c quarks. The detector includes a high-precision tracking system consisting of a silicon-strip vertex detector surrounding the pp interaction region, a large-area silicon-strip detector located upstream of a dipole magnet with bending power of about 4 Tm, and three stations of siliconstrip detectors and straw drift tubes placed downstream of the magnet. The tracking system provides a measurement of momentum, p, of charged particles with relative uncertainty that varies from 0.5% at low momentum to 1.0% at 200 GeV. The minimum distance of a track to a primary vertex (PV), the impact parameter (IP), is measured with a resolution of (15 + 29/p T ) µm, where p T is the component of the momentum transverse to the beam, in GeV. Different types of charged hadrons are distinguished using information from two ring-imaging Cherenkov detectors. Photons, electrons and hadrons are identified by a calorimeter system consisting of scintillating-pad and preshower detectors, an electromagnetic calorimeter and a hadronic calorimeter. Muons are identified by a system composed of alternating layers of iron and multiwire proportional chambers. The online event selection is performed by a trigger, which consists of a hardware stage, based on information from the calorimeter and muon systems, followed by a software stage, which applies a full event reconstruction.

Data selection
Candidate events for this analysis are first required to pass the hardware trigger, which selects muons with transverse momentum p T > 1.48 GeV in the 7 TeV data or p T > 1.76 GeV in the 8 TeV data. In the subsequent software trigger, at least one of the final-state particles is required to have p T > 1.7 GeV in the 7 TeV data or p T > 1.6 GeV in the 8 TeV data, unless the particle is identified as a muon in which case p T > 1.0 GeV is required. The final-state particles that satisfy these transverse momentum criteria are also required to have an impact parameter larger than 100 µm with respect to all of the primary pp interaction vertices (PVs) in the event. Finally, the tracks of two or more of the final-state particles are required to form a vertex that is significantly displaced from the PVs. In the subsequent offline selection, trigger signals are required to be associated with reconstructed particles in the signal decay chain.
[MeV]  Figure 1: Distribution of m K + K − near the φ peak before the φ candidate selection. Non-B + backgrounds have been subtracted using sPlot weights [36] obtained from a fit to the m J/ψ K + K − K + distribution. The default φ selection window is indicated with vertical red lines. The fit (solid blue line) of a Breit-Wigner φ signal shape plus two-body phase space function (dashed red line), convolved with a Gaussian resolution function, is superimposed.
The offline data selection is very similar to that described in Ref. [21], with J/ψ → µ + µ − candidates required to satisfy the following criteria: p T (µ) > 0.55 GeV, p T (J/ψ ) > 1.5 GeV, χ 2 per degree of freedom for the two muons to form a common vertex, χ 2 vtx (µ + µ − )/ndf < 9, and mass consistent with the J/ψ meson. Every charged track with p T > 0. 25 GeV, missing all PVs by at least 3 standard deviations (χ 2 IP (K) > 9) and classified as more likely to be a kaon than a pion according to the particle identification system, is considered a kaon candidate. The quantity χ 2 IP (K) is defined as the difference between the χ 2 of the PV reconstructed with and without the considered particle. Combinations of K + K − K + candidates that are consistent with originating from a common vertex with χ 2 vtx (K + K − K + )/ndf < 9 are selected. We combine J/ψ candidates with K + K − K + candidates to form B + candidates, which must satisfy χ 2 vtx (J/ψ K + K − K + )/ndf < 9, p T (B + ) > 2 GeV and have decay time greater than 0.25 ps. The J/ψ K + K − K + mass is calculated using the known J/ψ mass [37] and the B + vertex as constraints [38].
Four discriminating variables (x i ) are used in a likelihood ratio to improve the back-ground suppression: the minimal χ 2 IP (K), χ 2 vtx (J/ψ K + K − K + )/ndf, χ 2 IP (B + ), and the cosine of the largest opening angle between the J/ψ and the kaon transverse momenta. The latter peaks at positive values for the signal as the B + meson has high transverse momentum. Background events in which particles are combined from two different B decays peak at negative values, whilst those due to random combinations of particles are more uniformly distributed. The four signal probability density functions (PDFs), P sig (x i ), are obtained from simulated B + → J/ψ φK + decays. The background PDFs, P bkg (x i ), are obtained from candidates in data with J/ψ K + K − K + invariant mass between 5.6 and 6.4 GeV. We require −2 4 i=1 ln[P sig (x i )/P bkg (x i )] < 5, which retains about 90% of the signal events.
Relative to the data selection described in Ref. [21], the requirements on transverse momentum for µ and B + candidates have been lowered and the requirement on the multivariate signal-to-background log-likelihood difference was loosened. As a result, the B + signal yield per unit luminosity has increased by about 50% at the expense of somewhat higher background.
[MeV] The distribution of m K + K − for the selected B + → J/ψ K + K − K + candidates is shown in Fig. 1 (two entries per candidate). A fit with a P-wave relativistic Breit-Wigner shape on top of a two-body phase space distribution representing non-φ background, both convolved with a Gaussian resolution function with width of 1.2 MeV is superimposed. Integration of the fit components gives (5.3 ± 0.5)% of nonresonant background in the |m K + K − − 1020 MeV| < 15 MeV region used to define a φ candidate. To avoid reconstruction ambiguities, we require that there is exactly one φ candidate per J/ψ K + K − K + combination, which reduces the B + yield by 3.2%. The non-φ B + → J/ψ K + K − K + background in the remaining sample is small (2.1%) and neglected in the amplitude model. The related systematic uncertainty is estimated by tightening the φ mass selection window to ±7 MeV.
The mass distribution of the remaining J/ψ φK + combinations is shown in Fig. 2 together with a fit of the B + signal represented by a symmetric double-sided Crystal Ball function [39] on top of a quadratic function for the background. The fit yields 4289 ± 151 B + → J/ψ φK + events. Integration of the fit components in the 5270-5290 MeV region (twice the B + mass resolution on each side of its peak) used in the amplitude fits, gives a background fraction (β) of (23 ± 6)%. A Gaussian signal shape and a higher-order polynomial background function are used to assign systematic uncertainties which are included in, and dominate, the uncertainty given above. The B + invariant mass sidebands, 5225-5256 and 5304-5335 MeV, are used to parameterize the background in the amplitude fit.
The B + candidates for the amplitude analysis are kinematically constrained to the known B + mass [37] and to point to the closest pp interaction vertex. The measured value of m K + K − is used for the φ candidate mass, since the natural width of the φ resonance is larger than the detector resolution.

Matrix element model
We consider the three interfering processes corresponding to the following decay sequences: B + → K * + J/ψ , K * + → φK + (referred to as the K * decay chain), B + → XK + , X → J/ψ φ (X decay chain) and B + → Z + φ, Z + → J/ψ K + (Z decay chain), all followed by J/ψ → µ + µ − and φ → K + K − decays. Here, K * + , X and Z + should be understood as any φK + , J/ψ φ and J/ψ K + contribution, respectively. We construct a model of the matrix element (M) using the helicity formalism [40][41][42] in which the six independent variables fully describing the K * + decay chain are m φK , θ K * , θ J/ψ , θ φ , ∆φ K * ,J/ψ and ∆φ K * ,φ , where the helicity angle θ P is defined as the angle in the rest frame of P between the momentum of its decay product and the boost direction from the rest frame of the particle which decays to P , and ∆φ is the angle between the decay planes of the two particles (see Fig. 3). The set of angles is denoted by Ω. The explicit formulae for calculation of the angles in Ω are given in Appendix A. Figure 3: Definition of the θ K * , θ J/ψ , θ φ , ∆φ K * ,J/ψ and ∆φ K * ,φ angles describing angular correlations in B + → J/ψ K * + , J/ψ → µ + µ − , K * + → φK + , φ → K + K − decays (J/ψ is denoted as ψ in the figure).
The full six-dimensional (6D) matrix element for the K * decay chain is given by where the index j enumerates the different K * + resonances. The symbol J K * denotes the spin of the K * resonance, λ is the helicity (projection of the particle spin onto its momentum in the rest frame of its parent) and ∆λ µ ≡ λ µ + − λ µ − . The terms d J λ 1 ,λ 2 (θ) are the Wigner d-functions, R j (m φK ) is the mass dependence of the contribution and will be discussed in more detail later (usually a complex Breit-Wigner amplitude depending on resonance pole mass M 0 K * j and width Γ 0 K * j ). The coefficients A B→J/ψ K * λ J/ψ and A K * →φK λ φ are complex helicity couplings describing the (weak) B + and (strong) K * + decay dynamics, respectively. There are three independent complex A B→J/ψ K * λ J/ψ couplings to be fitted (λ J/ψ = −1, 0, 1) per K * resonance, unless J K * = 0 in which case there is only one since λ J/ψ = λ K * due to J B = 0. Parity conservation in the K * decay limits the number of independent helicity couplings A K * →φK λ φ . More generally parity conservation requires which, for the decay K * + → φK + , leads to This reduces the number of independent couplings in the K * decay to one or two. Since the overall magnitude and phase of these couplings can be absorbed in A B→J/ψ K * λ J/ψ , in practice the K * decay contributes zero or one complex parameter to be fitted per K * resonance.
The matrix element for the X decay chain can be parameterized using m J/ψ φ and the θ X , θ X J/ψ , θ X φ , ∆φ X,J/ψ , ∆φ X,φ angles. The angles θ X J/ψ and θ X φ are not the same as θ J/ψ and θ φ in the K * decay chain, since J/ψ and φ are produced in decays of different particles. For the same reason, the muon helicity states are different between the two decay chains, and an azimuthal rotation by angle α X is needed to align them as discussed below. The parameters needed to characterize the X decay chain, including α X , do not constitute new degrees of freedom since they can all be derived from m φK and Ω. The matrix element for the X decay chain also has unique helicity couplings and is given by where the index j enumerates all X resonances. To add M K * ∆λµ and M X ∆λµ coherently it is necessary to introduce the e iα X ∆λµ term, which corresponds to a rotation about the µ + momentum axis by the angle α X in the rest frame of J/ψ after arriving to it by a boost from the X rest frame. This realigns the coordinate axes for the muon helicity frame in the X and K * decay chains. This issue is discussed in Ref. [43] and at more length in Ref. [44].
The structure of helicity couplings in the X decay chain is different from the K * decay chain. The decay B + → XK + does not contribute any helicity couplings to the fit 3 , since X is produced fully polarized λ X = 0. The X decay contributes a resonance-dependent matrix of helicity couplings A X→J/ψ φ λ J/ψ ,λ φ . Fortunately, parity conservation reduces the number of independent complex couplings to one for J P X = 0 − , two for 0 + , three for 1 + , four for 1 − and 2 − , and at most five independent couplings for 2 + .
The matrix element for the Z + decay chain can be parameterized using m J/ψ K and the θ Z , θ Z J/ψ , θ Z φ , ∆φ Z,J/ψ , ∆φ Z,φ angles. The Z + decay chain also requires a rotation to align the muon frames to those used in the K * decay chain and to allow for the proper description of interference between the three decay chains. The full 6D matrix element is given by 3 There is one additional coupling, but that can be absorbed by a redefinition of X decay couplings, which are free parameters.

Parity conservation in the Z + decay requires
and provides a similar reduction of the couplings as discussed for the K * decay chain.
Instead of fitting the helicity couplings A A→B C λ B ,λ C as free parameters, after imposing parity conservation for the strong decays, it is convenient to express them by an equivalent number of independent LS couplings (B LS ), where L is the orbital angular momentum in the decay and S is the total spin of B and C, Possible combinations of L and S values are constrained via J A = L + S. The relation involves the Clebsch-Gordan coefficients (7) Parity conservation in the strong decays is imposed by Since the helicity or LS couplings not only shape the angular distributions but also describe the overall strength and phase of the given contribution relative to all other contributions in the matrix element, we separate these roles by always setting the coupling for the lowest L and S, B L min S min , for a given contribution to (1, 0) and multiplying the sum in Eq. (7) by a complex fit parameter A (this is equivalent to factoring out B L min S min ). This has an advantage when interpreting the numerical values of these parameters. The value of A j describes the relative magnitude and phase of the B L min S min j to the other contributions, and the fitted B LS j values correspond to the ratios, B LS j /B L min S min j , and determine the angular distributions.
Each contribution to the matrix element comes with its own R(m A ) function, which gives its dependence on the invariant mass of the intermediate resonance A in the decay chain (A = K * + , X or Z + ). Usually it is given by the Breit-Wigner amplitude, but there are special cases which we discuss below. An alternative parameterization of R(m A ) to represent coupled-channel cusps is discussed in Appendix D.
In principle, the width of the φ resonance should also be taken into account. However, since the φ resonance is very narrow (Γ 0 = 4.3 MeV, with mass resolution of 1.2 MeV) we omit the amplitude dependence on the invariant m K + K − mass from the φ decay.
A single resonant contribution in the decay chain B + → A..., A → ... is parameterized by the relativistic Breit-Wigner amplitude together with Blatt-Weisskopf functions, where is the Breit-Wigner amplitude including the mass-dependent width, Here, p is the momentum of the resonance A (K * + , X or Z + ) in the B + rest frame, and q is the momentum of one of the decay products of A in the rest frame of the A resonance. The symbols p 0 and q 0 are used to indicate values of these quantities at the resonance peak mass (m = M 0 ). The orbital angular momentum in B + decay is denoted as L B , and that in the decay of the resonance A as L A . The orbital angular momentum barrier factors, [45,46]: which account for the centrifugal barrier in the decay and depend on the momentum of the decay products in the rest frame of the decaying particle (p) as well as the size of the decaying particle (d). In this analysis we set this parameter to a nominal value of d = 3.0 GeV −1 , and vary it in between 1.5 and 5.0 GeV −1 in the evaluation of the systematic uncertainty.
In the helicity approach, each helicity state is a mixture of many different L values. We follow the usual approach of using in Eq. (9) the minimal L B and L A values allowed by the quantum numbers of the given resonance A, while higher values are used to estimate the systematic uncertainty.
We set BW (m) = 1.0 for the nonresonant (NR) contributions, which means assuming that both magnitude and phase have negligible m dependence. As the available phase space in the B + → J/ψ φK + decays is small (the energy release is only 12% of the B + mass) this is a well-justified assumption. We consider possible mass dependence of NR amplitudes as a source of systematic uncertainties.

Maximum likelihood fit of amplitude models
The signal PDF, P sig , is proportional to the matrix element squared, which is a function of six independent variables: m φK and the independent angular variables in the K * decay chain Ω. The PDF also depends on the fit parameters, ω, which include the helicity couplings, and masses and widths of resonances. The two other invariant masses, m φK and m J/ψ K , and the angular variables describing the X and Z + decay chains depend on m φK and Ω, therefore they do not represent independent dimensions. The signal PDF is given by: where M(m φK , Ω| ω) is the matrix element given by Eq. (5). Φ(m φK ) = p q is the phase space function, where p is the momentum of the φK + (i.e. K * ) system in the B + rest frame, and q is the K + momentum in the K * + rest frame. The function (m φK , Ω) is the signal efficiency, and I( ω) is the normalization integral, where the sum is over simulated events, which are generated uniformly in B + decay phase space and passed through the detector simulation [47] and data selection. In the simulation, pp collisions producing B + mesons are generated using Pythia [48] with a specific LHCb configuration [49]. The weights w MC j introduced in Eq. (18) contain corrections to the B + production kinematics in the generation and to the detector response to bring the simulations into better agreement with the data. Setting w MC j = 1 is one of the variations considered when evaluating systematic uncertainties. The simulation sample contains 132 000 events, approximately 30 times the signal size in data. This procedure folds the detector response into the model and allows a direct determination of the parameters of interest from the uncorrected data. The resulting log-likelihood sums over the data events (here for illustration, P = P sig ), where the last term does not depend on ω and can be dropped (N is the total number of the events in the fit). In addition to the signal PDF, P sig (m φK , Ω| ω), the background PDF, P bkg (m φK , Ω) determined from the B + mass peak sidebands, is included. We minimize the negative log-likelihood defined as where β is the background fraction in the peak region determined from the fit to the m J/ψ φK distribution ( Fig. 2), P u bkg (m φK , Ω) is the unnormalized background density proportional to the density of sideband events, with its normalization determined by 4 The equation above implies that the background term is efficiency corrected, so it can be added to the efficiency-independent signal probability expressed by |M| 2 . This way the efficiency parameterization, (m φK , Ω), becomes a part of the background description which affects only a small part of the total PDF. The efficiency parameterization in the background term is assumed to factorize as The 1 (m φK , cos θ K * ) term is obtained by binning a two-dimensional (2D) histogram of the simulated signal events. Each event is given a 1/(p q) weight, since at the generator level the phase space is flat in cos θ K * but has a p q dependence on m φK . A bi-cubic function is used to interpolate between bin centers. The 1 (m φK , cos θ K * ) efficiency and its visualization across the normal Dalitz plane are shown in Fig. 4. The other terms are again built from 2D histograms, but with each bin divided by the number of simulated events in the corresponding m φK slice to remove the dependence on this mass (Fig. 5).
The background PDF, P u bkg (m φK , Ω)/Φ(m φK ), is built using the same approach, The background function P bkg 1 (m φK , cos θ K * ) is shown in Fig. 6 and the other terms are shown in Fig. 7. The fit fraction (FF) of any component R is defined as, where in M R all terms except those associated with the R amplitude are set to zero. [MeV] . Function values corresponding to the color encoding are given on the right. By construction each function integrates to unity at each m φK value.

Background-subtracted and efficiency-corrected distributions
The background-subtracted and efficiency-corrected Dalitz plots are shown in Figs. 8-10 and the mass projections are shown in Figs. 11-13. The latter indicates that the efficiency corrections are rather minor. The background is eliminated by subtracting the scaled B + sideband distributions. The efficiency corrections are achieved by weighting events according to the inverse of the parameterized 6D efficiency given by Eq. (22). The efficiency-corrected signal yield remains similar to the signal candidate count, because we normalize the efficiency to unity when averaged over the phase space. While the m φK distribution (Fig. 11) does not contain any obvious resonance peaks, it would be premature to conclude that there are none since all K * + resonances expected in this mass range belong to higher excitations, and therefore should be broad. In fact the narrowest known K * + resonance in this mass range has a width of approximately 150 MeV [37]. Scattering experiments sensitive to K * → φK decays also showed a smooth mass distribution, which revealed some resonant activity only after partial-wave analysis [50][51][52]. Therefore, studies of angular distributions in correlation with m φK are necessary. Using full 6D correlations results in the best sensitivity.
The m J/ψ φ distribution (Fig. 12) contains several peaking structures, which could be exotic or could be reflections of conventional K * + resonances. There is no narrow X(4140) peak just above the kinematic threshold, consistent with the LHCb analysis presented in Ref. [21], however we observe a broad enhancement. A peaking structure is observed at about 4300 MeV. The high mass region is inspected with good sensitivity for the first time, with the rate having a minimum near 4640 MeV with two broad peaks on each side. The m J/ψ K distribution ( Fig. 13) peaks broadly in the middle and has a high-mass peak, which is strongly correlated with the low-mass m J/ψ φ enhancement (Fig. 10).
As explained in the previous section, the amplitude fits are performed by maximizing the unbinned likelihood on the selected signal candidates including background events and without the efficiency weights. To properly represent the fit quality, the fit projections in the following sections show the fitted data sample, i.e. including the background and without the parameterized efficiency corrections applied to the signal events.

Amplitude model with only φK + contributions
We first try to describe the data with kaon excitations alone. Their mass spectrum as predicted in the relativistic potential model by Godfrey-Isgur [53] is shown in Fig. 14 together with the experimentally determined masses of both well-established and unconfirmed K * resonances [37]. Past experiments on K * states decaying to φK [50][51][52] had limited precision, especially at high masses, gave somewhat inconsistent results, and provided evidence for only a few of the states expected from the quark model in the 1513-2182 MeV range probed in our data set. However, except for the J P = 0 + states which cannot decay to φK because of angular momentum and parity conservation, all other kaon excitations above the φK threshold are predicted to decay to this final state [54]. In B + decays, production of high spin states, like the K * 3 (1780) or K * 4 (2045) resonances, is expected to be suppressed by the high orbital angular momentum required to produce them.
We have used the predictions of the Godfrey-Isgur model as a guide to the quantum numbers of the K * + states to be included in the model. The masses and widths of all states are left free; thus our fits do not depend on detailed predictions of Ref. [53], nor on previous measurements. We also allow a constant nonresonant amplitude with J P = 1 + , since such φK + contributions can be produced, and can decay, in S-wave. Allowing the magnitude of the nonresonant amplitude to vary with m φK does not improve fit qualities.
While it is possible to describe the m φK and m J/ψ K distributions well with K * contributions alone, the fit projections onto m J/ψ φ do not provide an acceptable description of the data. For illustration we show in Fig. 15 the projection of a fit with the following composition: a nonresonant term plus candidates for two 2P 1 , two 1D 2 , and one of each of 1 3 F 3 , 1 3 D 1 , 3 3 S 1 , 3 1 S 0 , 2 3 P 2 , 1 3 F 2 , 1 3 D 3 and 1 3 F 4 states, labeled here with their intrinsic quantum numbers: n 2S+1 L J (n is the radial quantum number, S the total spin of the valence quarks, L the orbital angular momentum between quarks, and J the total angular momentum of the bound state). The fit contains 104 free parameters. The χ 2 value (144.9/68 bins) between the fit projection and the observed m J/ψ φ distribution corresponds to a p value below 10 −7 . Adding more resonances does not change the conclusion that non-K * contributions are needed to describe the data.

Amplitude model with φK + and J/ψ φ contributions
We have explored adding X and Z + contributions of various quantum numbers to the fit models. Only X contributions lead to significant improvements in the description of the data. The default resonance model is described in detail below and is summarized in Table 3, where the results are also compared with the previous measurements and the theoretical predictions forsu states [53]. The model contains seven K * + states, four X states and φK + and J/ψ φ nonresonant components. There are 98 free parameters in this fit. Projections of the fit onto the mass variables are displayed in Fig. 16. The systematic uncertainties are obtained from the sum in quadrature of the changes observed in the fit results when: the K * + and X(4140) models are varied; the Breit-Wigner amplitude parameterization is modified; only the left or right B + mass peak sidebands are used for the background parameterization; the φ mass selection window is made narrower by a factor of two (to reduce the non-φ background fraction); the signal and background shapes are varied in the fit to m J/ψ φK which determines the background fraction β; the weights assigned to simulated events, in order to improve agreement with the data on B + production characteristics and detector efficiency, are removed. More detailed discussion of the systematic uncertainties can be found in Appendix B.
The significance of each (non)resonant contribution is calculated assuming that ∆(−2ln L), after the contribution is included in the fit, follows a χ 2 distribution with the number of degrees of freedom (ndf) equal to the number of free parameters in its parameterization. The value of ndf is doubled when M 0 and Γ 0 are free parameters in the fit. The validity of this assumption has been verified using simulated pseudoexperiments. The significances of the X contributions are given after accounting for systematic variations. Combined significances of exotic contributions, determined by removing more than one exotic contribution at a time, are much larger than their individual significances given in Table 3. The significance of the spin-parity determination for each X state is determined as described in Appendix C.
The longitudinal (f L ) and transverse (f ⊥ ) polarizations are calculated for K * + contributions according to where Among the K * + states, the J P = 1 + partial wave has the largest total fit fraction (given by Eq. (24)). We describe it with three heavily interfering contributions: a nonresonant term and two resonances. The significance of the nonresonant amplitude cannot be quantified, since when it is removed one of the resonances becomes very broad, taking over its role. Evidence for the first 1 + resonance is significant (7.6σ). We include a second resonance in the model, even though it is not significant (1.9σ), because two states are expected in the quark model. We remove it as a systematic variation. The 1 + states included in our model appear in the mass range where two 2P 1 states are predicted (see Table 3), and where the K − p → φK − p scattering experiment found evidence for a 1 + state with M 0 ∼ 1840 MeV, Γ 0 ∼ 250 MeV [50], also seen in the K − p → K − π + π − p scattering data [55]. Within the large uncertainties the lower mass state is also consistent with the unconfirmed K 1 (1650) state [37], based on evidence from the K − p → φK − p scattering experiment [51].
There is also a substantial 2 − contribution to the amplitude model. When modeled as a single resonance (5.0σ significant), M 0 = 1889 ± 27 MeV and Γ 0 = 376 ± 94 MeV are obtained in agreement with the evidence from the K − p → φK − p scattering data which yielded a mass of around 1840 MeV and a width of order 250 MeV [50]. The K + p → φK + p scattering data also supported such a state at 1810 ± 20 MeV, but with a narrower width, 140 ± 40 MeV [51]. Since two closely spaced 2 − states are established from other decay modes [37], and since two 1D 2 states are predicted, we allow two resonances in the default fit. The statistical significance of the second state is 3σ. The masses and widths obtained by the fit to our data are in good agreement with the parameters of the K 2 (1770) and K 2 (1820) states and in agreement with the predicted masses of the 1D 2 states ( Table 3). The individual fit fractions are poorly defined, and not quoted, because of large destructive interferences. There is no evidence for an additional 2 − state in our data (which could be the expected 2D 2 state [53]), but we consider the inclusion of such a state among the systematic variations.
The most significant K * + resonance in our data is a vector state (8.5σ). Its mass and width are in very good agreement with the well-established K * (1680) state, which is observed here in the φK decay mode for the first time, and fits the 1 3 D 1 interpretation. When allowing an extra 1 − state (candidate for 3 3 S 1 ), its significance is 2.6σ with a mass of 1853 ± 5 MeV, but with a width of only 33 ± 11 MeV, which cannot be accommodated in thesu quark model. When limiting the width to be 100 MeV or more, the significance drops to 1.4σ. We do not include it in the default model, but consider its inclusion as a systematic variation. We also include among the considered variations the effect of an insignificant (< 0.3σ) tail from the K * (1410) resonance. Table 3: Results for significances, masses, widths and fit fractions of the components included in the default amplitude model. The first (second) errors are statistical (systematic). Errors on f L and f ⊥ are statistical only. Possible interpretations in terms of kaon excitation levels are given, with notation n 2S+1 L J , together with the masses predicted in the Godfrey-Isgur model [53]. Comparisons with the previously experimentally observed kaon excitations [37] and X → J/ψ φ structures are also given.

Contri-sign.
Fit results bution or Ref.
There is a significant (5.4σ) 2 + contribution, which we describe with one very broad resonance, consistent with the claims of a K * 2 (1980) state seen in other decays and also consistent with a broad enhancement seen in K − p → φK 0 n scattering data [52]. It can be interpreted as the 2 3 P 2 state predicted in this mass range. An extra 2 + state added to the model, as suggested e.g. by the possibility that the 1 3 F 2 state is in the probed mass range, is less than 0.7σ significant and is considered among the systematic variations. There is also 3.5σ evidence for a 0 − contribution, consistent with the previously unconfirmed K(1830) state seen in K − p → φK − p scattering data [50]. It could be a 3 1 S 0 state. An extra 0 − state added to the model (e.g. 4 1 S 0 ) is less than 0.2σ significant and is considered among the systematic variations.
We consider among the systematic variations the inclusion of several further states that are found not to be significant in the fit. These are a 3 + state (e.g. 1F 3 , < 1.8σ), a 4 + state (e.g. 1 3 F 4 , < 2.0σ or < 0.6σ if fixed to the K * 4 (2045) parameters [37]) or a 3 − state (e.g. 1 3 D 3 , < 2.0σ if fixed to the K * 3 (1780) parameters). Overall, the K * + composition of our data is in good agreement with the expectations for thesu states, and also in agreement with previous experimental results on K * states in this mass range. These results add significantly to the knowledge of K * spectroscopy.
A near-threshold J/ψ φ structure in our data is the most significant (8.4σ) exotic contribution to our model. We determine its quantum numbers to be J P C = 1 ++ at 5.7σ significance (Appendix C). When fitted as a resonance, its mass (4146.5 ± 4.5 +4.6 −2.8 MeV) is in excellent agreement with previous measurements for the X(4140) state, although the width (83 ± 21 +21 −14 MeV) is substantially larger. The upper limit which we previously set for production of a narrow (Γ = 15.3 MeV) X(4140) state based on a small subset of our present data [21] does not apply to such a broad resonance, i.e. the present results are consistent with our previous analysis. The statistical power of the present data sample is not sufficient to study its phase motion [56]. A model-dependent study discussed in Appendix D suggests that the X(4140) structure may be affected by the nearby D ± s D * ∓ s coupled-channel threshold. However, larger data samples will be required to resolve this issue.
We establish the existence of the X(4274) structure with statistical significance of 6.0σ, at a mass of 4273.3 ± 8.3 +17.2 − 3.6 MeV and a width of 56.2 ± 10.9 + 8.4 −11.1 MeV. Its quantum numbers are also 1 ++ at 5.8σ significance. Due to interference effects, the data peak above the pole mass, underlining the importance of proper amplitude analysis.
The high J/ψφ mass region also shows structures that cannot be described in a model containing only K * + states. These features are best described in our model by two J P C = 0 ++ resonances at 4506 ± 11 +12 −15 MeV and 4704 ± 10 +14 −24 MeV, with widths of 92 ± 21 +21 −20 MeV and 120 ± 31 +42 −33 MeV, and significances of 6.1σ and 5.6σ, respectively. The resonances interfere with a nonresonant 0 ++ J/ψ φ contribution that is also significant (6.4σ). The significances of the quantum number determinations for the high mass states are 4.0σ and 4.5σ, respectively.
Additional X resonances of any J P value (J ≤ 2) added to our model have less than 2σ significance. A modest improvement in fit quality can be achieved by adding Z + → J/ψ K + resonances to our model, however the significance of such contributions is too small to justify introducing an exotic hadron contribution (at most 3.1σ without accounting for systematic uncertainty). The parameters obtained for the default model components stay within their systematic uncertainties when such extra X or Z + contributions are introduced.
[MeV] )      Fig. 16 for a description of the components.

Summary
In summary, we have performed the first amplitude analysis of B + → J/ψ φK + decays. We have obtained a good description of data in the 6D phase space including invariant masses and decay angles.
Even though no peaking structures are observed in the φK + mass distributions, correlations in the decay angles reveal a rich spectrum of K * + resonances. In addition to the angular information contained in the K * + and φ decays, the J/ψ decay also helps to probe these resonances, as the helicity states of the K * + and J/ψ mesons coming from the B + decay must be equal. Unlike the earlier scattering experiments investigating K * → φK decays, we have good sensitivity to states with both natural and unnatural J P combinations.
The dominant 1 + partial wave (FF = 42 ± 8 +5 −9 %) has a substantial nonresonant component, and at least one resonance that is 7.6σ significant. There is also 2σ evidence that this structure can be better described with two resonances matching the expectations for the two 2P 1 excitations of the kaon.
Also prominent is the 2 − partial wave (FF = 10.8 ± 2.8 +1.5 −4.6 %). It contains at least one resonance at 5.0σ significance. This structure is also better described with two resonances at 3.0σ significance. Their masses and widths are in good agreement with the well established K 2 (1770) and K 2 (1820) states, matching the predictions for the two 1D 2 kaon excitations.
The 1 − partial wave (FF = 6.7 ± 1.9 +3.2 −3.9 %) exhibits 8.5σ evidence for a resonance which matches the K * (1680) state, which was well established in other decay modes, and matches the expectations for the 1 3 D 1 kaon excitation. This is the first observation of its decay to the φK final state.
The 2 + partial wave has a smaller intensity (FF = 2.9 ± 0.8 +1.7 −0.7 %), but provides 5.4σ evidence for a broad structure that is consistent with the K * 2 (1980) state observed previously in other decay modes and matching the expectations for the 2 3 P 2 state.
We also confirm the K(1830) state (3 1 S 0 candidate) at 3.5σ significance (FF = 2.6 ± 1.1 +2.3 −1.8 %), earlier observed in the φK decay by the K − p scattering experiment. We determine its mass and width with properly evaluated uncertainties for the first time.
Overall, our K * + → φK + results show excellent consistency with the states observed in other experiments, often in other decay modes, and fit the mass spectrum predicted for the kaon excitations by the Godfrey-Isgur model. Most of the K * + structures we observe were previously observed or hinted at by the Kp → φK(p or n) experiments, which were, however, sometimes inconsistent with each other.
Our data cannot be well described without several J/ψ φ contributions. The significance of the near-threshold X(4140) structure is 8.4σ with FF= 13.0 ± 3.2 +4.8 −2.0 %. Its width is substantially larger than previously determined. We determine the J P C quantum numbers of this structure to be 1 ++ at 5.7σ. This has a large impact on its possible interpretations, in particular ruling out the 0 ++ or 2 ++ D * + s D * − s molecular models [3][4][5][6][7]10]. The below-J/ψ φ-threshold D ± s D * ∓ s cusp [11,20] may have an impact on the X(4140) structure, but more data will be required to address this issue. The existence of the X(4274) structure is established (6σ) with FF= 7.1 ± 2.5 +3.5 −2.4 % and its quantum numbers are determined to be 1 ++ (5.8σ). Together, these two J P C = 1 ++ contributions have a fit fraction of 16 ± 3 +6 −2 %. Molecular bound-states or cusps cannot account for the X(4274) J P C values. A hybrid charmonium state would have 1 −+ [17,18]. Some tetraquark models expected 0 −+ , 1 −+ [13] or 0 ++ , 2 ++ [14] state(s) in this mass range. A tetraquark model implemented by Stancu [12] not only correctly assigned 1 ++ to X(4140), but also predicted a second 1 ++ state at mass not much higher than the X(4274) mass. Calculations by Anisovich et al. [15] based on the diquark tetraquark model predicted only one 1 ++ state at a somewhat higher mass. Lebed-Polosa [16] predicted the X(4140) peak to be a 1 ++ tetraquark, although they expected the X(4274) peak to be a 0 −+ state in the same model. A lattice QCD calculation with diquark operators found no evidence for a 1 ++ tetraquark below 4.2 GeV [57].
The high J/ψ φ mass region is investigated with good sensitivity for the first time and shows very significant structures, which can be described as 0 ++ contributions (FF = 28 ± 5 ± 7%) with a nonresonant term plus two new resonances: X(4500) (6.1σ significant) and X(4700) (5.6σ). The quantum numbers of these states are determined with significances of more than 4σ. The work of Wang et al. [9] predicted a virtual D * + s D * − s state at 4.48 ± 0.17 GeV. None of the observed J/ψ φ states is consistent with the state seen in two-photon collisions by the Belle collaboration [32].

A Calculation of decay angles
The decay angles are calculated in a way analogous to that documented in Appendix IX of Ref. [43]. The five angles for each decay chain are: three helicity angles of J/ψ , φ and of the resonance in question (e.g. K * ) and two angles between the decay plane of the resonance and the decay plane of either J/ψ or φ. In addition, a rotation is needed to align the muon helicity frames of the X and Z + decay chains to that of the K * in order to properly describe the interferences. The choice of K * as the reference decay chain is arbitrary. The cosine of a helicity angle of particle P , produced in two-body decay A → P B, and decaying to two particles P → C D is calculated from (Eq. (16) in Ref. [43]) where the momentum vectors are in the rest frame of the particle P . For the B + → J/ψ K * + decay, the angle between the J/ψ → µ + µ − and the K * + → φK + decay planes is calculated from 6 (Eqs. (14)- (15) in Ref. [43]) ∆φ K * ,J/ψ = atan2(sin ∆φ K * ,J/ψ , cos ∆φ K * ,J/ψ ) (29) cos with all vectors being in the B + rest frame. For the B + → Z + φ decay, the angle between the Z + → J/ψ K + and the φ → K + K − decay planes, ∆φ Z,φ , can be calculated in the same way with J/ψ → φ, µ + → K + (the K + from the φ decay) and the accompanying K + staying the same. The angle between the decay planes of two sequential decays, e.g. between the Z + → J/ψ K + and J/ψ → µ + µ − decay planes after the B + → Z + φ decay, is calculated from (Eqs. (18)- (19) in Ref. [43]) ∆φ Z,J/ψ = atan2(sin ∆φ Z,J/ψ , cos ∆φ Z,J/ψ ) (34) with all vectors being in the J/ψ rest frame. The other angles of this type are calculated in the same way, with appropriate substitutions. For example, ∆φ K * ,φ between the K * + → φK + and φ → K + K − decay planes after B + → K * + J/ψ decay, is calculated substituting φ −→ J/ψ , µ + −→ K + (K + from the φ decay), and the accompanying K + staying the same (all vectors are in the φ rest frame here). The angle aligning the muon helicity frames between the K * + and Z + decay chains is calculated from (Eqs. (20)- (21) in Ref. [43]) where the K + is the accompanying kaon and all vectors are in the J/ψ rest frame. Similarly, α X is obtained from the above equations with the K + −→ φ substitution. For the charge-conjugate B − → J/ψ φK − decays, the same formulae apply with the accompanying kaon being K − , µ + replaced by µ − and K + from the φ decay replaced by the K − from the φ decay. All azimuthal angles (∆φ and α) have their signs flipped after applying the formulae above (see the bottom of Appendix IX in Ref. [43]).

B Systematic uncertainty
Individual systematic uncertainties on masses, widths and fit fractions are presented for K * + contributions in Tables 4-5, and for X contributions in Table 6. Positive and negative deviations are summed in quadrature separately for total systematic uncertainties. The statistical uncertainties are included for comparison.
In many instances, the uncertainty in the K * + model composition is the dominant systematic uncertainty. The K * + model variations include adding the following contributions (one-by-one) to the default amplitude model: second 0 − , 1 − or 2 + states, a third 2 − state, the 3 − K * 3 (1780) state, a 3 + state, the 4 + K * 4 (2045) state, and the below threshold 1 − K * (1410) state. The variations also include omitting the second 1 + or 2 − states. The observed deviations in the fit parameters are added in quadrature and then listed in Tables 4-6.
The other sizable source of systematic uncertainty is due to the L B and L K * (or L X ) dependence of the Breit-Wigner amplitude in the numerator of Eq. (9) via Blatt-Weisskopf factors. Helicity states correspond to mixtures of allowed L values, but we assume the lowest L values in Eq. (9) in the default fit. We increase L B values by 1 for all the components (one-by-one). Values of L K * or L X can only differ by an even number because of parity conservation in strong decays. We performed such variations for states in which the higher value is allowed, except for the X states, since the fit results indicate that the higher L X amplitudes are insignificant. Again, the observed deviations in the fit parameters are added in quadrature and then listed in Tables 4-6.
The energy release in the B + → J/ψ φK + decay is small (∼ 13% on M B ), and the phase space is very limited, not offering much range for nonresonant interactions to change. In the default model the nonresonant terms are represented by constant amplitudes. When allowing them to change exponentially with mass-squared, exp(−α m 2 ), the slope parameters, α, are consistent with zero. The observed deviations in the measured parameters are included among the systematic contributions.
Replacing the Breit-Wigner amplitude for the X(4140) structure with a D ± s D * ∓ s cusp in one particular model (see Appendix D) is included among the systematic model variations.
The dependence on mass of the total resonance width (Eq. (11)) used in the default fit assumes that it is dominated by the observed decay mode. All K * + states are expected to have sizable widths to the other decay modes, Kπ, Kρ, K * (892)π etc. However, ratios of these partial widths to the φK partial width are unknown. To check the related systematic uncertainty, we perform an alternative fit (marked Γ tot in the tables) in which the mass dependence of the width is set by the lightest possible decay mode allowed: Kπ for natural spin-parity resonances and Kω for the others. This includes changing the L K * value.
The Blatt-Weisskopf factors contain the d parameter for the effective hadron size [58] (Eq. (12)), which we set to 3.0 GeV −1 in the default fit. As a systematic variation we change its value between 1.5 and 5.0 GeV −1 .
To address the systematic uncertainty in the background parameterization, we perform amplitude fits with either the left or right B + mass sideband only. The default fit uses both.
We perform fits to m J/ψ φK with alternative signal and background parameterizations to determine the systematic uncertainty on the background fraction in the signal region (β). The largest deviation in its value (∆β/β = +25%) is then used in the alternative amplitude fit to the data.
In the default fit, the simulated events used for the efficiency corrections are weighted to improve the agreement between the data and the simulation. The total Monte Carlo event weight (w M C ) is a product of weights determined as the ratios between the signal distributions in the data and in the simulated sample (generated according to the preliminary amplitude model) as functions of p T (B + ), number of charged tracks in the event, and each kaon momentum. These weights are intended to correct for any inaccuracies in simulation of pp collisions, of B + production kinematics and in kaon identification. To account for the uncertainty associated with the efficiency modelling we include among the systematic variations a fit in which the weights are not applied.
To check the uncertainty related to non-φ background, we reduce its fraction by narrowing the φ → K + K − mass selection window by a factor of two. This also accounts for any uncertainty related to averaging over this mass in the amplitude fit.
As a cross-check on both the background subtraction and the efficiency corrections the minimal value of p T for kaon candidates is changed from 0.25 GeV to 0.5 GeV, which reduces the background fraction by 54% (β = 10.4%) and the signal efficiency by 20%, as illustrated in Fig. 21. The mass projections of the fit are shown in Fig. 22. The fit results are within the assigned total uncertainties as shown at the bottom rows of Tables 4-6.
More details on the systematic error evaluations can be found in Ref. [59]. Table 4: Summary of the systematic uncertainties on the parameters of the K * + → φK + states with J P = 2 − and 1 + . The kaon p T cross-check results are given at the bottom. All numbers for masses and widths are in MeV and fit fractions in %.   )  Table 5: Summary of the systematic uncertainties on the parameters of the K * + → φK + states with J P = 0 − , 1 − and 2 + . The kaon p T cross-check results are given at the bottom. All numbers for masses and widths are in MeV and fit fractions in %.

C Spin analysis for the X → J/ψ φ states
To determine the quantum numbers of each X state, fits are done under alternative J P C hypotheses. The likelihood-ratio test is used to quantify rejection of these hypotheses. Since different spin-parity assignments are represented by different functions in the angular part of the fit PDF, they represent separate hypotheses. For two models representing separate hypotheses, assuming a χ 2 distribution with one degree of freedom for ∆(−2ln L) under the disfavored J P C hypothesis gives a lower limit on the significance of its rejection [60].
The results for the default fit approach are shown in Table 7. The J P C values of the X(4140) and X(4274) states are both determined to be 1 ++ with 7.6σ and 6.4σ significance, respectively. The quantum numbers of X(4500) and of X(4700) states are both established to be 0 ++ at 5.2σ and 4.9σ level, respectively. The separation from the alternative J P C hypothesis with likelihood closest to that for the favored quantum numbers in the default fit is studied for each state under the fit variations which have dominant effects on the resonance parameters as shown in Table 8. The lowest values are taken for the final significances of the quantum number determinations: 5.7σ for X(4140), 5.8σ for X(4274), 4.0σ for X(4500) and 4.5σ for X(4700). 7.6σ 7.2σ 5.6σ 6.8σ 2 −+ 9.6σ 6.4σ 6.5σ 6.3σ Table 8: Significance, in standard deviations, of J P C preference for the X states for dominant systematic variations of the fit model. The label "L + n" specifies which L value in Eq. (9) is increased relative to its minimal value and by how much (n). The lowest significance value for each state is highlighted.
While our 1 ++ assignment to X(4140) and its large width rule out an interpretation as a 0 ++ or 2 ++ D * + s D * − s molecule (for which 1 ++ is not allowed [3]) with large ∼ 83 MeV binding energy as suggested by many authors [3][4][5][6][7], such a structure could be formed by molecular forces in a D ± s D * ∓ s pair in S-wave [11,20]. Since the sum of D ± s and D * ∓ s masses (4080 MeV) is below the J/ψ φ mass threshold (4116 MeV), such a contribution would not be described by the Breit-Wigner function with a pole above that threshold. The investigation of all possible parameterizations for such contributions, which are model dependent, goes beyond the scope of this analysis. However, we attempt a fit with a simple threshold cusp parameterization proposed by Swanson (Ref. [61] and private communications), in which the introduction of an exponential form factor, with a momentum scale (β 0 ) characterizing the hadron size, makes the cusp peak slightly above the sum of masses of the rescattering mesons. While controversial [62], this model provided a successful description of the Z c (3900) + and Z c (4025) + exotic meson candidates with masses peaking slightly above the molecular thresholds [61].
In Swanson's model a virtual loop with two mesons A and B inside ( Fig. 1 left in Ref. [61]) contributes, in the non-relativistic near-threshold approximation, the following amplitude, where m is J/ψ φ mass, µ AB = M A M B /(M A + M B ) is the reduced mass of the pair, β 0 is a hadronic scale of order of Λ QCD , (which can be AB dependent), is a very small number ( → 0), l is the angular momentum between A and B. The lowest l values are expected to dominate. The amplitude Π(m) reflects coupled-channel kinematics. The above integral can be conveniently expressed as where −Z is the scaled mass deviation from the AB threshold. For l = 0, the integral above evaluates to For masses below the AB threshold Z > 0 and I(Z) (thus Π(Z)) has no imaginary part. For masses above the threshold Z < 0, √ Z is imaginary, which leads to both real and imaginary parts. The real and imaginary parts of −I(Z) as a function of −Z are shown in Fig. 23, while the corresponding Argand diagram is shown in Fig. 24 where it is compared to the phase motion of the Breit-Wigner function.  [61]. See the text for a more precise explanation.
The function Π(m) replaces the Breit-Wigner function BW (m|M 0 , Γ 0 ) in Eq. (10). The Blatt-Weisskopf functions in Eq. (9) still apply. Thus, the functional form of this representation, has three free parameters to determine from the data (β 0 and the complex S-wave helicity coupling). The value of β 0 obtained by the fit to the data, 297 ± 20 MeV, is close to the value of 300 MeV with which Swanson was successful in describing the other near-threshold exotic meson candidates [61]. A fit with such parameterization (see Fig. 25 for mass distributions), has a better likelihood than the Breit-Wigner fit by 1.6σ for the default model (8 free parameters in the X(4140) Breit-Wigner parameterization), and better by 3σ when only S-wave couplings are allowed (4 free parameters), providing an indication that the X(4140) structure may not be a bound state that can be described by the Breit-Wigner formula. Larger data samples will be required to obtain more insight. We have included the X(4140) cusp model among the systematic variations considered for parameters of the other fit components. The differences between the results obtained with the default amplitude model and the model in which the X(4140) structure is represented by a cusp are given in Tables 4-6.
The X(4274) mass structure can be reasonably well described by the 0 −+ cusp model for D ± s D * s0 (2317) ∓ scattering (Fig. 26). However, the multidimensional likelihood is substantially worse than for the default amplitude model (6.6σ). The likelihood remains worse for the default fit even if 1 ++ quantum numbers are assumed for such a cusp (4.4σ).   [61]. Motion with the mass is counter-clockwise. The peak amplitude is reached at threshold when the real part is maximal and the imaginary part is zero. The Breit-Wigner amplitude gives circular phase motion, also with counter-clockwise mass evolution, with maximum magnitude when zero is crossed on the real axis.
This particular cusp parameterization is not useful when trying to describe any of the higher mass J/ψ φ structures.  Figure 25: Distributions of (top left) φK + , (top right) J/ψ K + and (bottom) J/ψ φ invariant masses for the B + → J/ψ φK + data (black data points) compared with the results of the amplitude fit containing K * + → φK + and X → J/ψ φ contributions in which X(4140) is represented as a J P C = 1 ++ D + s D * − s cusp. The total fit is given by the red points with error bars. Individual fit components are also shown.
[MeV]  Figure 26: Distributions of J/ψ φ invariant mass for the B + → J/ψ φK + data (black data points) compared with the results of the amplitude fit containing K * + → φK + and X → J/ψ φ contributions in which X(4140) and X(4274) are represented as J P C = 1 ++ D ± s D * ∓ s and 0 −+ D ± s D * s0 (2317) ∓ cusps, respectively. The total fit is given by the red points with error bars. Individual fit components are also shown.