Observation of the $B^+ \rightarrow J/\psi \eta^{\prime} K^+$ decay

The $B^+ \rightarrow J/\psi \eta^{\prime} K^+$ decay is observed for the first time using proton-proton collision data collected by the LHCb experiment at centre-of-mass energies of 7, 8, and 13TeV, corresponding to a total integrated luminosity of 9fb$^{-1}$. The branching fraction of this decay is measured relative to the known branching fraction of the $B^+ \rightarrow \psi(2S) K^+$ decays and found to be $$ \frac{\mathcal{B}( B^+ \rightarrow J/\psi \eta^{\prime}K^+)}{\mathcal{B}( B^+ \rightarrow \psi(2S)K^+)} = \left(4.91\pm 0.47\pm0.29\pm0.07\right)\times10^{-2}, $$ where the first uncertainty is statistical, the second is systematic and the third is related to external branching fractions. A first look at the $J/\psi\eta^{\prime}$ mass distribution is performed and no signal of intermediate resonances is observed.

This paper reports the observation of the B + → J/ψη ′ K + decay using proton-proton (pp) collision data collected by the LHCb experiment at centre-of-mass energies of 7, 8, and 13 TeV, corresponding to a total integrated luminosity of 9 fb −1 .The measurement of its branching fraction normalised to the well-known branching fraction of the B + → ψ(2S)K + decay [37], is performed using the η ′ → ρ 0 γ decay.The observation of the signal is confirmed using the η ′ → ηπ + π − decay mode, which is also used as a cross-check.

Detector and simulation
The LHCb detector [59,60] 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 [61], a large-area silicon-strip detector located upstream of a dipole magnet with a bending power of about 4 Tm, and three stations of silicon-strip detectors and straw drift tubes [62,63] placed downstream of the magnet.The tracking system provides a measurement of the momentum, p, of charged particles with a relative uncertainty that varies from 0.5% at low momentum to 1.0% at 200 GeV/c.The minimum distance of a track to a primary pp collision vertex (PV), the impact parameter 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/c.Different types of charged hadrons are distinguished using information from two ring-imaging Cherenkov detectors [64].Photons, electrons and hadrons are identified by a calorimeter system consisting of scintillating-pad and preshower detectors, an electromagnetic and a hadronic calorimeter.Muons are identified by a system composed of alternating layers of iron and multiwire proportional chambers [65].
The online event selection is performed by a trigger [66,67], 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.At the hardware trigger stage, events are required to have a muon track with high transverse momentum or dimuon candidates in which the product of the p T of the muons has a high value.In the software trigger, two oppositely charged muons are required to form a good-quality vertex that is significantly displaced from every PV, with a dimuon mass exceeding 2.7 GeV/c 2 .
Simulated events are used to describe signal shapes and to compute the efficiencies needed to determine the branching fraction ratio.In the simulation, pp collisions are generated using Pythia [68] with a specific LHCb configuration [69].Decays of unstable particles are described by EvtGen [70], in which final-state radiation is generated using Photos [71].The interaction of the generated particles with the detector, and its response, are implemented using the Geant4 toolkit [72] as described in Ref. [73].The p T and rapidity (y) spectra of the B + mesons in simulation are corrected to match distributions in data.The correction factors are calculated by comparing the observed p T and y spectra for a high-purity data sample of reconstructed B + → J/ψK + decays with the corresponding simulated samples.In the simulation, the B + → J/ψη ′ K + decays are generated as phase-space decays and corrected using a gradient boosted decision tree reweighting algorithm [74] to reproduce the J/ψη ′ and η ′ K + mass spectra observed in data.To describe accurately the variables used for kaon identification, the corresponding quantities in simulation are resampled according to values obtained from calibration data samples of D * + → (D 0 → K − π + ) π + decays [75].The procedure accounts for correlations between the variables associated with a particular track, as well as the dependence of the kaon identification response on the track's p T and η and the multiplicity of tracks in the event.To account for imperfections in the simulation of charged-particle reconstruction, the track reconstruction efficiency is corrected using a sample of J/ψ → µ + µ − decays in data [76].Samples of the B + → J/ψK * + decays with K * + → K + (π 0 → γγ) are used to correct the photon reconstruction efficiency in simulation [77][78][79][80].

Event selection
The B + → J/ψη ′ K + candidates are reconstructed with η ′ decays to either (ρ 0 → π + π − ) γ or (η → γγ) π + π − final states.The difficulty of reconstructing photons in the η → γγ decay leads to a sample with fewer events.The B + → ψ(2S)K + normalisation decay is reconstructed using the ψ(2S) → J/ψπ + π − decay.In both signal and normalisation channels, the J/ψ meson is reconstructed in its decay to two muons.As explained in detail below, an initial loose selection is applied for both signal and normalisation channels.Subsequently, for the B + → J/ψη ′ K + candidates, where the background level is large, a multivariate estimator is used to select higher purity subset of candidates.The normali-sation channel has a high purity after the initial selection, therefore no further selection steps are applied.
To reduce systematic uncertainties, the initial selection criteria for both signal and normalisation channels are kept the same whenever possible.The selection criteria are chosen to be similar to those used in previous LHCb studies [20,21,23,29,53,77,78].The muon, pion and kaon candidates are identified by combining information from the Cherenkov detectors, calorimeters and muon detectors [81] associated with the reconstructed tracks.To reduce the combinatorial background, only tracks that are inconsistent with originating from any reconstructed PV in the event are considered.The transverse momentum of the muon candidates is required to be greater than 500 MeV/c and their momenta must exceed 6 GeV/c.Pairs of oppositely-charged muons consistent with originating from a common vertex are combined to form J/ψ → µ + µ − candidates.The reconstructed mass of the muon pair is required to be between 3.056 and 3.136 GeV/c 2 .
Tracks that are consistent with the pion or kaon hypotheses are required to have transverse momentum greater than 200 MeV/c.Photons are reconstructed from clusters in the electromagnetic calorimeter, with transverse energy above 350 MeV.The clusters must not be associated with reconstructed tracks [82,83].Photon identification is based on the combined information from electromagnetic and hadronic calorimeters, scintillation pad, preshower detectors and the tracking system.
For the reconstruction of the η ′ → ηπ + π − candidates, two photons are first combined to form an η candidate.The diphoton mass is restricted to lie within ±60 MeV/c 2 around the known mass of the η meson [37].Each η candidate is then combined with two oppositely-charged pions to form an η ′ candidate.The mass of the combination is required to lie within ±45 MeV/c 2 around the known mass of the η ′ meson [37].For an η ′ candidate reconstructed in the η ′ → ρ 0 γ decay mode, the ρ 0 candidate is formed from two oppositely-charged pions.The mass of this candidate is restricted to lie between 500 and 900 MeV/c 2 .This asymmetric region around the known mass of the ρ 0 meson [37] takes into account the shift of the ρ 0 line shape, due to the electric-dipole nature of the η ′ → ρ 0 γ transition [84][85][86][87][88][89].A photon is combined with the ρ 0 candidate in order to form the η ′ candidate, whose mass is required to lie within ±30 MeV/c 2 of the known η ′ mass [37].
Each selected J/ψ candidate is combined with a kaon track and either an η ′ candidate or two oppositely-charged pions to form a B + candidate decaying into the signal or normalisation modes, respectively.For the B + → ψ(2S)K + candidates the J/ψπ + π − mass is required to be between 3.66 and 3.71 GeV/c 2 .To improve the B + meson mass resolution a kinematic fit [90] is performed, which constrains the masses of the J/ψ, η ′ and η candidates to their known values [37], and the B + candidates to originate from its associated PV.The decay time of the B + candidates is required to be greater than 100 µm/c to suppress the large combinatorial background from tracks created in a PV.
Further selection of the B + → J/ψη ′ K + decays is based on a multivariate estimator, in the following referred to as the multi-layer perceptron (MLP) classifier.The classifier is based on an artificial neural network algorithm [91,92], configured with a cross-entropy cost estimator [93].It reduces the combinatorial background to a low level while retaining a high signal efficiency.Two MLP classifiers are trained separately for the two different η ′ meson decay modes.The list of variables used for classifiers includes the χ2 of the kinematic fit; transverse momenta of the η ′ , kaon and pion candidates; pseudorapidities of pion and kaon candidates; transverse momentum of the photon from the η ′ → ρ 0 γ decay or miminal trans-verse momentum of photons from the η ′ → (η → γγ) π + π − decay; decay time of the B + candidate; variable related to the quality of kaon identification [64,81] and, for the η ′ → ρ 0 γ decay, cosine of the angle between momenta of the π + and η ′ candidates in the rest frame of the ρ 0 candidate.The classifiers are trained using simulated samples of B + → J/ψη ′ K + decays as signal proxy, while the B + → J/ψη ′ K + candidates from data with mass above 5. 35 GeV/c 2 are used to represent the background.The B + → J/ψ (η ′ → ρ 0 γ) K + candidates with the J/ψπ + π − K + mass consistent with the known mass of the B + meson are vetoed to avoid contamination from the B + → J/ψπ + π − K + decays with a random photon added.
The requirement on each of the MLP classifiers is chosen to maximise the figure-of-merit defined as S/ √ B + S, where S represents the expected signal yield, and B is the expected background yield within a ±15 MeV/c 2 mass window centred around the known mass of the B + meson and corresponding to approximately three times the mass resolution on both sides of the peak.The background yield is calculated from fits to data, as described in Sec. 4, while the expected signal yield is estimated as S = εS 0 , where ε is the efficiency of the requirement on the response of the MLP classifier determined from simulation, and S 0 is the signal yield obtained from the fit to the data, with a loose requirement applied on the response of the MLP classifier. 2 The mass distributions for the selected B + → J/ψη ′ K + candidates are shown in Fig. 1, where clear signal peaks corresponding to B + decays are seen in data for both η ′ decay modes.

Signal yield determination
The signal yields are determined using an extended unbinned maximum-likelihood fit to the J/ψη ′ K + mass distributions with a two-component function.The signal component for both cases is modelled by a modified Gaussian function that combines a Gaussian core with power-law tails on both sides of the distribution [94,95].The background component is parameterised by a second-order positive polynomial function [96] for the η ′ → ρ 0 γ decay mode and by an exponential function in case of the η ′ → ηπ + π − decay mode.The parameters of the detector resolution function are taken from simulation, and the width of the Gaussian function is further corrected by a scale factor, s B + , that accounts for a small discrepancy between data and simulation [15,26,29,97,98].To account for the uncertainty in the tail parameters and resolution, the fit is performed simultaneously for data and simulated samples, sharing the same tail parameters, and allowing the correction factor s B + to vary.The resulting fit functions are overlaid with data distributions in Fig. 1 and the signal yields are found to be where the uncertainties are statistical only.The resolution correction factors s B + are found to be 1.08 ± 0.12 and 1.03 ± 0.16 for the η ′ → ρ 0 γ and η ′ → ηπ + π − samples, respectively.In both cases, the statistical significance of the B + → J/ψη ′ K + signal is calculated using Wilks' theorem [99] and found to exceed 17 and 12 standard deviations for the η ′ → ρ 0 γ and η ′ → ηπ + π − samples, respectively.However, as the signal yield is  much lower for the η ′ meson decays to the ηπ + π − final state, all subsequent studies are performed using only the η ′ → ρ 0 γ decay mode.
The background-subtracted J/ψη ′ , η ′ K + , and J/ψK + mass spectra from the B + → J/ψ (η ′ → ρ 0 γ) K + decays are shown in Figs.2(a-c), where the sPlot technique [100] based on the fit results is used for background subtraction.The J/ψη ′ , η ′ K + , and J/ψK + masses are calculated using a kinematic fit with J/ψ, η ′ and B + mass constraints and a PV constraint applied [90].While for the J/ψK + mass the distribution largely agrees with the shape expected from the phase-space model, for the low-mass region of the η ′ K + mass spectrum and the high-mass region of the J/ψη ′ mass spectrum a striking difference from the phase-space model is observed.These differences are potentially due to contributions from decays via intermediate heavy excited strange mesons, such as K * 0 (1430) + , K * 2 (1430) + or K * (1680) + mesons, decaying into the η ′ K + final state.The decays of the B + mesons into a J/ψ meson and heavy excited strange mesons have been studied in Refs.[13,14,16].The decays via intermediate excited kaons also contribute to the higher mass region of the J/ψη ′ mass spectrum, as shown in Fig. 2(d).The J/ψη ′ mass region below 4.7 GeV/c 2 is explicitly inspected for possible contributions from decays via excited charmonium or charmonium-like states into the J/ψη ′ final state.Fits to the background-subtracted J/ψη ′ mass distribution are performed in individual mass windows, corresponding to the well-established ψ(4160), ψ(4230), ψ(4360), ψ(4415) and ψ(4430) resonances [37].For each fit, the resonance shape is parameterised with a relativistic Breit-Wigner function convoluted with a mass resolution function.The non-resonant contribution is modelled by a first order positive polynomial function.The known masses and widths of the resonances [37] are introduced in the fits as Gaussian constraints on the corresponding parameters.The resolution function is modelled by the modified Gaussian function with parameters obtained using simulation as a function of the J/ψη ′ mass.No statistically significant signals are observed for the B + → J/ψη ′ K + decays via intermediate resonances, listed above.To probe the contribution of the resonances with higher masses, more advanced fit techniques accounnting for complicated background shape and the distortion of the signal Breit-Wigner shape are required.
For the determination of the resonant structure of the B + → J/ψη ′ K + decay, a full amplitude analysis, similar to those used in Refs.[14,16], is required.Large signal yields Figure 2: Normalised background-subtracted (a) η ′ K + , (b) J/ψη ′ , (c) J/ψK + mass spectra and (d) two-dimensional mass distribution of η ′ K + vs J/ψη ′ from the B + → J/ψη ′ K + decays.Superimposed curves are the expectations from a phase-space model.
and low background levels are important prerequisites for such analysis.The relatively large level of combinatorial background for the B + → J/ψη ′ K + signal decays with the η ′ meson reconstructed via the η ′ → ρ 0 γ decay mode, makes it difficult to carry out an amplitude analysis.With a larger data sample, expected from future data-taking periods, it will be possible to perform the full amplitude analysis using the B + → J/ψη ′ K + signal decays with the η ′ meson reconstructed via the η ′ → ηπ + π − decay mode, where the combinatorial background is smaller.
determined using an unbinned fit to this distribution with a two-component function.The component corresponding to the B + → ψ(2S)K + decays is parameterised by the modified Gaussian function.The component describing the B + → J/ψπ + π − K + decays without intermediate ψ(2S) state is modelled by a phase-space function,3 modified by a first order positive polynomial function.The tail and resolution parameters for the signal component are taken from simulation with the resolution further corrected by a scale factor, s ψ(2S) , that accounts for a small discrepancy between data and simulation [26,102].The fit is performed simultaneously to data and simulated samples, as for the signal mode described in Sec. 4. From this fit the number of B + → ψ(2S)K + decays is found to be where the uncertainty is statistical only.The correction factor s ψ(2S) is found to be 1.057 ± 0.008, which is in good agreement with results from Refs.[15,26,97,98,102,103].

Branching fraction ratio computation
The ratio of the branching fractions R, defined in Eq. ( 1), is calculated as where N B + →J/ψη ′ K + and N B + →ψ(2S)K + are the yields from Eqs. (2a) and (3), and ϵ B + →J/ψη ′ K + and ϵ B + →ψ(2S)K + are the efficiencies to reconstruct the observed final states.The efficiencies are the products of detector acceptance, reconstruction, selection and trigger efficiencies, , where x * ≡ (x − x min )/(x max − x min ), and x min , x max denote the minimal and maximal values of x, respectively [101].Here, k = 3 and n = 4 are used.

Systematic uncertainties
The signal and normalisation channels share the same set of final-state charged particles, and the same trigger and preselection requirements are applied to both.This allows many systematic uncertainties to cancel in the ratio R. The remaining nonnegligible contributions are listed in Table 1.
Several systematic uncertainties are associated with the corrections applied to the simulation.The finite size of the B + → J/ψK + signal sample used for correction of the simulated transverse momentum and rapidity spectra of B + mesons induces an uncertainty on the B + meson p T and y spectra.In turn, this uncertainty induces small changes in the ratio of efficiencies.The corresponding spread of these changes amounts to 0.1% and is taken as the systematic uncertainty related to the B + meson kinematic.
The decay model corrections for the B + → J/ψη ′ K + decay are obtained using the algorithm described in Ref. [74].The systematic uncertainty related to the correction method is estimated by varying the configuration parameters of the algorithm.The largest deviation of the efficiency value from the baseline tuning is found to be 1.1%, which is assigned as systematic uncertainty associated with the B + decay model.
There are residual differences in the reconstruction efficiency of charged-particle tracks that do not cancel completely in the ratio of total efficiencies given the slightly different kinematic distributions of the final-state particles.The track-finding efficiencies obtained from simulated samples are corrected using calibration channels [76].The uncertainties related to the efficiency correction factors are propagated to the ratios of the total efficiencies using pseudoexperiments, and are found to be 0.7%.This value is taken as the systematic uncertainty due to the tracking efficiency calibration.
Differences in the photon reconstruction efficiencies between data and simulation are studied using a large sample of B + → J/ψK * + decays, reconstructed using the K * + → K + (π 0 → γγ) decay mode [77][78][79][80].The uncertainty due to the finite size of the sample is propagated to the ratio of the total efficiencies using pseudoexperiments and is found to be less than 1.0%.The uncertainty due to the accuracy of the B + → J/ψK * + branching fraction [37] is 3.5%.These two values are added in quadrature to obtain a systematic uncertainty related to the photon reconstruction of 3.6%.
The kaon identification variable used for the MLP estimator is drawn from calibration data samples and has a dependence on the particle kinematics and track multiplicity.Systematic uncertainties in this procedure arise from the limited size of both the simulation and calibration samples, and the modelling of the particle identification variable.The limitations due to the size of the simulation and calibration samples are evaluated by using bootstrapping techniques [104,105], creating multiple samples and repeating the procedure for each of these.The impact of potential mismodelling of the kaon identification variable is evaluated by describing the corresponding distributions using density estimates with different kernel widths [75,106].For each of these cases, alternative efficiency maps are produced to determine the associated uncertainties.A systematic uncertainty of 2.8% is assigned from the observed differences with alternative efficiency maps.
A systematic uncertainty related to the knowledge of the trigger efficiencies has been previously studied using large samples of B + → (J/ψ → µ + µ − ) K + and B + → (ψ(2S) → µ + µ − ) K + decays by comparing the ratios of the trigger efficiencies in data and simulation [107].Based on this comparison, a relative uncertainty of 1.1% is assigned.
The remaining inconsistency between data and simulation, not covered by the corrections discussed in Sec. 2, is estimated by varying the requirement on the response of the MLP classifier in ranges that lead to changes in the measured signal yields as large as ±20%.The resulting difference in the data-simulation efficiency ratio is found to be 3.0%.
A different class of systematic uncertainties directly affects the fit itself, namely uncertainties associated with the fit models used to describe the J/ψ η ′ K + , ψ(2S)K + and J/ψπ + π − spectra.The systematic uncertainty is accounted for by fits with alternative models.The list of alternative models to describe signal B + and ψ(2S) components includes a modified Apollonios function [108], which has exponential instead of power-law tails, a generalised Student's t-distribution [109,110] and a modified Novosibirsk function [111].For the combinatorial background for both the signal and the normalisation channels, a third order polynomial function and an exponential function multiplied by a first order polynomial function are chosen as alternative models.For each fit only one component (either signal or background) is replaced at a time, and the same fit function is used for both signal and normalisation channel.The largest deviation of the ratio of signal yields is found to be 0.6%.As alternative background models for J/ψπ + π − candidates, a first order polynomial function and a product of an exponential function with a first order polynomial function are used.The largest deviation of the signal yield is found to be 1.5%.The two deviations are added in quadrature to obtain a 1.6% systematic uncertainty due to imperfect knowledge of the signal and background shapes.
Finally, the finite size of the simulation samples contributes an uncertainty of 0.9% on the ratio of total efficiencies.The total systematic uncertainty for the ratio of branching fractions R is calculated as the sum in quadrature of all the values listed above and is found to be 6.0%.
The statistical significance for the B + → J/ψη ′ K + decay is recalculated using Wilks' theorem for each alternative fit model, and the smallest values of 17 and 10 standard deviations for the η ′ → ρ 0 γ and η ′ → ηπ + π − cases, respectively, are taken as the significance including the systematic uncertainty.

Results and summary
The B + → J/ψη ′ K + decays are observed for the first time using proton-proton collision data collected by the LHCb experiment at centre-of-mass energies of 7, 8, and 13 TeV, corresponding to a total integrated luminosity of 9 fb −1 .In the analysis, the η ′ meson is reconstructed from the two ρ 0 γ and ηπ + π − final states.The signal significance exceeds 10 standard deviations for both modes.The branching fraction of the B + → J/ψη ′ K + decay is measured for the η ′ → ρ 0 γ sample through normalisation to the known branching fraction of the B + → ψ(2S)K + decay [37].The ratio of branching fractions is found to be B(B + → J/ψη ′ K + ) B(B + → ψ(2S)K + ) = (4.91 ± 0.47 ± 0.29 ± 0.07) where the first uncertainty is statistical, the second is systematic and the third is related to the uncertainties on the branching fractions of the intermediate resonances.
The absolute branching fraction is determined using the known branching fraction of B + → ψ(2S)K + decays, B(B + → ψ(2S)K + ) = (6.24± 0.20) × 10 −4 [37], and is found to be B(B + → J/ψη ′ K + ) = (3.06± 0.29 ± 0.18 ± 0.04) × 10 −5 , where the first uncertainty is statistical, the second is systematic and the third is due to external branching fractions uncertainties.The measured branching fraction is consistent with the upper limit previously set by the Belle collaboration [58].An inspection of the J/ψη ′ mass spectrum shows no significant contributions from the decays via intermediate charmonium or charmonium-like resonances.

Table 1 :
Summary of systematic uncertainties on the ratio of branching fractions R. The overall systematic uncertainty is calculated as a sum in quadrature of all the sources.