Observation of the B$_\mathrm{c}^+$ meson in PbPb and pp collisions at $\sqrt{s_{\mathrm{NN}}} =$ 5.02 TeV and measurement of its nuclear modification factor

The B$_\mathrm{c}^+$ meson is observed for the first time in heavy ion collisions. Data from the CMS detector are used to study the production of the B$_\mathrm{c}^+$ meson in lead-lead (PbPb) and proton-proton (pp) collisions at a center-of-mass energy per nucleon pair of $\sqrt{s_{\mathrm{NN}}} =$ 5.02 TeV, via the B$_\mathrm{c}^+$ $\to$ (J/$\psi$ $\to$ $\mu^+\mu^-$)$\mu^+\nu_\mu$ decay. The B$_\mathrm{c}^+$ nuclear modification factor, derived from the PbPb-to-pp ratio of production cross sections, is measured in two bins of the trimuon transverse momentum and of the PbPb collision centrality. The B$_\mathrm{c}^+$ meson is shown to be less suppressed than quarkonia and most of the open heavy-flavor mesons, suggesting that effects of the hot and dense nuclear matter created in heavy ion collisions contribute to its production. This measurement sets forth a promising new probe of the interplay of suppression and enhancement mechanisms in the production of heavy-flavor mesons in the quark-gluon plasma.


1
The very high temperatures and densities reached in high-energy heavy ion collisions allow quarks and gluons to form a deconfined state of matter, often referred to as the quark-gluon plasma (QGP). Among the various signatures of the QGP formation, the suppression of heavy quarkonia (e.g. J/ψ or Υ mesons) due to the screening of the heavy-quark potential has been intensely studied since its proposal by Matsui and Satz [1]. The strong J/ψ suppression observed in heavy ion collisions at the Super Proton Synchrotron [2] and Relativistic Heavy Ion Collider [3] was indeed qualitatively consistent with color screening effects. However, lead-lead (PbPb) collisions at the Large Hadron Collider (LHC) reach higher temperatures, but show less J/ψ suppression [4]. This observation is interpreted as arising from the formation of bound states of charm quarks originating from different hard scatterings, a mechanism referred to as recombination [5,6]. By contrast, the bottomonium LHC data show no evidence for recombination, consistent with the relatively small b-quark production cross section. In addition, the suppression of Υ(nS) states is sequential (more suppression for smaller binding energies) [7], as expected from color screening.
While understanding the effects of the QGP on heavy-quark bound states is progressing [8], quantifying them still remains a key challenge. The B + c meson is composed of a charm quark and a bottom antiquark, and is intermediate in size and binding energy between the J/ψ and Υ(1S) mesons [9]. Since it bridges the gap between charmonia and bottomonia, the modification of B + c yields contributes to probing heavy quark interactions at extreme temperatures. The recombination of heavy quarks can also affect the B + c yield. The production rates of J/ψ and B + c mesons from this process should scale quadratically and linearly, respectively, with the large number of charm quarks in the QGP, up to hundreds in central PbPb collisions. However, as the B + c production cross section in proton-proton (pp) scatterings is small, owing to the dominant production mechanism requiring both bb and cc pairs, the enhancement of B + c production in heavy ion collisions may greatly exceed that of the J/ψ meson [9,10]. These recombination effects would manifest more strongly at low transverse momentum (p T ) [9].
For p T m(B + c ), B + c mesons are produced predominantly via heavy-quark fragmentation [11,12], and are therefore sensitive to the energy loss of a massive color triplet charge in the QGPpossibly causing the suppression observed for other B mesons [13,14], J/ψ mesons from B decays [15], and D mesons [16,17]. Conversely, the modification of prompt J/ψ meson production for p T m(J/ψ) [15,18] probes the energy loss of a massive color octet state [19]. Therefore, comparing the B + c yield with that of other heavy flavor mesons at large p T can probe both the mass dependence of energy loss (from a possible dead-cone effect [20]) and its color charge dependence.
The B + c meson was first observed in proton-antiproton collisions at the Tevatron in the B + c → J/ψ + ν decay mode [21]. Its ground and excited states were then studied in pp collisions at the LHC [22][23][24][25][26]. In this Letter, the first observation of B + c mesons produced in heavy ion collisions is reported, and their cross sections are measured and compared in PbPb and pp collisions. The data were collected with the CMS detector in 2017 for pp and in 2018 for PbPb collisions at the same center-of-mass energy per nucleon pair, √ s NN = 5.02 TeV, corresponding to integrated luminosities of 302 pb −1 and 1.61 nb −1 , respectively. The signal is reconstructed from the three muons in the B + c → (J/ψ → µ + µ − )µ + ν µ decay mode. While this mode features a neutrino that prevents a full reconstruction of the decay, it has a much larger branching fraction than neutrinoless decay channels [24]. In this Letter, charge-conjugate states are implied, and the quoted cross sections correspond to the sum of B + c and B − c mesons. The results are presented in two kinematic regions that are defined in terms of the vector sum of the three muon momenta, and whose limits are chosen based on the single-muon acceptance of the CMS apparatus: a low-p µµµ T bin, 6 < p µµµ T < 11 GeV with rapidity 1.3 < |y µµµ | < 2.3, and a high-p µµµ T bin, 11 < p µµµ T < 35 GeV with |y µµµ | < 2.3. In simulations, the trimuon p T is, on average, about 15% smaller than the B + c p T . In PbPb collisions, the analysis is performed in the 0-90% centrality range, where centrality refers to the fraction of the inelastic nucleusnucleus cross section, with lower values denoting a larger overlap of the nuclei [27]. The results integrated over the two kinematic regions are also presented, separated in the centrality ranges 0-20 and 20-90%. To reduce potential biases, the analysis was performed in a "blind" way: the algorithms and selection procedures were finalized and formally approved using a quarter of the PbPb data, before examining the entire sample. Tabulated results are provided in a HEPData record [28].
The central feature of the CMS apparatus [29] is a superconducting solenoid providing a magnetic field of 3.8 T. Within the solenoid volume are a silicon pixel and strip tracker, a lead tungstate crystal electromagnetic calorimeter, and a brass and scintillator hadron calorimeter, each composed of a barrel and two endcap sections. Hadron forward calorimeters extend the pseudorapidity coverage to 3 < |η| < 5, and the sum of the transverse energy deposited in them is used to estimate the collision centrality. Muons are detected in gas-ionization detectors embedded in the steel flux-return yoke outside the solenoid and covering |η| < 2.4. Muons with p T > 1. Events of interest are selected using a two-tiered trigger system [32]. The first level, composed of custom hardware processors, uses information from the calorimeters and muon detectors [33]. The high-level trigger consists of a farm of processors running a fast version of the full event reconstruction software. The events used in this analysis were selected by triggers designed to collect all events containing a J/ψ meson, hence requiring two muons, without p T requirements. Loose criteria on the single-muon quality and the dimuon mass and opening angle are also applied in the PbPb high-level trigger, as in Refs. [34,35].
Monte Carlo (MC) simulations are used for various signal and background studies, and for estimating the acceptance and efficiency of the reconstruction, triggering, and selection. The B + c mesons are generated with BCVEGPY2.2 [12], while their decays are handled with EVT-GEN1.3 [36]. The underlying event is generated with PYTHIA8.212 [37], tune CP5 [38]. PYTHIA8 is also used to generate prompt and nonprompt (from B meson decays) J/ψ background samples. To simulate PbPb collisions, the generated events are embedded into simulated PbPb collisions created using HYDJET1.8 [39]. All samples are passed to GEANT4 [40] to simulate the detector response, and then reconstructed with the same software as the collision data.
The B + c → (J/ψ → µ + µ − )µ + ν µ decay features three muons originating from the same displaced vertex, an opposite-sign muon pair consistent with the J/ψ mass, and a trimuon invariant mass m µµµ between m J/ψ + m µ 3.2 GeV and m B + c 6.3 GeV. Three main background sources can mimic this topology. Fake J/ψ events arise when neither of the opposite-sign muon pairs originate from a J/ψ decay. It is estimated by summing the trimuon mass distributions obtained in the lower and higher dimuon mass sidebands of J/ψ candidates. The second category (B decays) comes from b-hadrons (excluding B + c ) decaying to a true J/ψ meson associated with a muon (usually a misidentified hadron) from the same b-hadron decay. It is estimated via simulation, where the p T spectrum is corrected using nonprompt J/ψ production measurements [15]. Its normalization is unconstrained to cover a possible mismodeling of the muon misidentification rate. The third contribution (J/ψ + random X) combines a true J/ψ meson with a muon candidate (usually an uncorrelated misidentified hadron) from another decay. It is estimated in data by rotating the momentum and decay vertex of J/ψ candidates around the collision vertex before associating them with third muon candidates. Several azimuthal rotation angles (excluding the vicinity of the original J/ψ meson) are used, with or without inverting rapidity. In PbPb collisions, the associated muons are mostly uncorrelated with the J/ψ meson, so the distributions from various rotation angles are identical (within statistical uncertainties) and averaged, with a data-derived normalization (fixed in the fit). In pp collisions, significant residual J/ψ − µ correlations lead to different distributions for different rotation angles, which is accounted for by considering various mixes of these distributions.
The offline selection includes the same event-level and single-muon identification criteria as in Refs. [34,35]. Loose kinematic acceptance criteria are applied to the muon candidates, matching the efficient region for the two triggering muons, and even looser for the third one. At least one of the two opposite-sign dimuon combinations must have an invariant mass in the J/ψ peak region, or in the sidebands used for background estimation. The sideband and peak regions are both asymmetric to account for radiative tails, and are separated by small gaps. The total sideband width equals that of the peak region, from 180 to 260 MeV depending on the muon pseudorapidity (which affects the mass resolution). For the trimuons having two oppositesign dimuons in the studied mass regions (5-6% of the overall sample), the two corresponding trimuon candidates are kept, weighted by the probability of the chosen dimuon to be a true J/ψ meson. This probability is extracted from the dimuon mass distribution from events with only one J/ψ candidate in the signal or sideband regions.
Requirements are also set on the probability of the trimuon vertex fit, the significance of its displacement from the collision vertex, the angle between the trimuon momentum and the segment joining the collision and trimuon vertices, the invariant mass corrected for the momentum of the neutrino transverse to the B + c momentum direction, and the sum of the angular separations ∆R = √ (∆η) 2 + (∆φ) 2 between the three muon pairs.
After the selection, the simulated signal and the three background samples are used to train a boosted decision tree (BDT) using the TMVA package [41]. This combines and optimizes the discriminating power of eight variables: the five discussed in the previous paragraph, the imbalance between the p T of the J/ψ and of the third muon, the ratio of the ∆R of the J/ψ muons to the sum of the ∆R values from the other two dimuon combinations, and the significance of the displacement from the collision vertex for the non-J/ψ muon.
Candidates with very low values of the resulting discriminant BDT variable (hence very high background probability) are rejected, losing only 0.1% in signal efficiency. For each analysis bin, low, medium, and high BDT intervals are set to contain about 25, 40, and 35%, respectively, of the expected signal. The first and last intervals are dominated by background and signal, respectively. A binned likelihood fit of the pp or PbPb trimuon mass distributions provides the signal yields. Using ROOFIT [42], templates from the signal and the three backgrounds are simultaneously fitted in the three BDT intervals, and in either two kinematic bins, two centrality bins, or the whole kinematic range. The BDT distribution of the sum of the fitted templates is checked against that of data, and, in pp collisions, corrected before re-running the template fit.
The results of the fits in the three BDT intervals and integrated over the two kinematic regions are shown for pp and PbPb collisions in Fig. 1. In each BDT bin, the signal purity and the measured yield, N(B c ), are given. The wrong-sign distributions, containing three same-sign muons in data, are superimposed to illustrate that the purely combinatorial background is easily rejected. The normalizations of the fake J/ψ sample, and of the J/ψ + random X sample in PbPb collisions, are provided by the data. In PbPb collisions, the J/ψ + random X and fake J/ψ backgrounds are dominant. In pp collisions, the region above the B + c mass strongly constrains the J/ψ + random X contribution, and the remaining background comes from B decays and fake J/ψ events. The signal yields extracted from the fit are corrected for the acceptance and efficiency of the reconstruction, triggering, and selection. These are calculated in each analysis bin using the simulated signal trimuons. The simulated efficiencies of single muon reconstruction, identification and triggering are corrected by a tag-and-probe method using the J/ψ resonance, similarly to Refs. [15,34]. The acceptance and efficiency are evaluated iteratively by first performing the p µµµ T -differential analysis using the original simulation. The resulting corrected yields are fitted to correct the p µµµ T spectrum of the simulation before a second run of the analysis. This p µµµ T spectrum is then corrected again based on the second-step results, notably improving upon the initial acceptance and efficiency estimation.
The corrected yields are divided by the pp integrated luminosity [43] or by its PbPb equivalent, the number of minimum bias PbPb hadronic collisions N MB times the nuclear overlap function T PbPb from Ref. [27]. The PbPb-to-pp ratio of these pp-equivalent normalized yields then provides the nuclear modification factor, R AA . In case of no modification by the medium, R AA is expected to be equal to unity.
Uncertainties arise from different sources: statistical, background (shapes and normalizations), choice of the fit method, muon efficiency, B + c kinematics (acceptance and efficiency), contamination from other B + c decays, and overall normalization. The fit uncertainties, ranging from 5 to 9% in pp and 17 to 31% in PbPb collisions, include the purely statistical and the background uncertain-ties. The latter are implemented via nuisance parameters allowing variations of the trimuon mass templates, such as controlling their statistical uncertainties with the Barlow-Beeston procedure [44], varying the fake J/ψ background between the lower and higher dimuon sideband, or varying the rotation angles in the J/ψ + random X background. Variations of the fit method are also considered, such as changing the m µµµ or BDT bin size, neglecting the low-BDT bin, using a BDT variable whose m µµµ dependence is subtracted, or regularizing the low-statistics templates instead of using the Barlow-Beeston procedure. The resulting uncertainty remains below 7% (12%) in pp (PbPb) collisions. The uncertainty from the tag-and-probe derived muon efficiency corrections is 2-5%.
Since the B + c kinematic distributions are not precisely known, acceptance and efficiency corrections are recalculated 1500 times with p µµµ T spectra fitted on variations of the measured p µµµ Tdifferential yields within the above-mentioned uncertainties. For the p µµµ T -integrated results, the root-mean-square (RMS) of the varied acceptance and efficiency corrections, of order 7 and 24% for pp and PbPb collisions, respectively, is used as the systematic uncertainty related to the B + c kinematics. For the p µµµ T dependence, these variations are correlated with the other uncertainty sources, so the combined uncertainty is assessed as the RMS of the varied corrected yields. The correlation between the variations of the spectrum and of the acceptance and efficiency is small or negative for the PbPb high-p µµµ T bin and for both pp p T bins, so that the uncertainties with or without this systematic effect are similar. This correlation is large and positive for the PbPb low-p µµµ T bin, inducing an additional 12-31% uncertainty. The uncertainty in the PbPb-to-pp ratio is the RMS of the ratios of the relevant varied quantities.
The contamination from other B + c decays, such as B + c → J/ψ(τ + → µ + X) ν τ or B + c → (cc → J/ψX) µ + ν µ , where X denotes any decay product(s), is estimated to be below 4.5%, and to have largely cancelling pp and PbPb contributions. The overall normalization uncertainty arising from the luminosity and centrality determination is 1.9-3.8%. The leading uncertainties in the p µµµ T -differential and p µµµ T -integrated measurements are from the fit and the B + c kinematics, respectively.
The significance of the B + c signal in PbPb collisions, calculated from the fit likelihood ratio and including the fit method uncertainty, is well above 5 standard deviations. The left panel of Fig. 2 shows the measured B + c meson p µµµ T -differential cross sections in pp and (pp-equivalent) PbPb collisions. The two bins of the trimuon p T correspond to different rapidity ranges. The markers of the p µµµ T bins are placed according to the Lafferty-Wyatt prescription [45]. The bin-to-bin correlation factor ρ 1-2 is also displayed. The filled and empty rectangles show the fit and total uncertainties, respectively. The ratio between the low-p The other panels of Fig. 2 show the B + c nuclear modification factor, i.e., the ratio of the (ppequivalent) PbPb to pp cross sections, as a function of p µµµ T (middle) and of centrality (right). The markers of the p µµµ T bins are placed at the average of their values for pp and PbPb collisions, while the centrality bin markers are placed at the minimum bias average number of participants N part . The filled and empty rectangles, respectively, show the bin-to-bin-uncorrelated and total uncertainties, such that the uncertainty in the difference of the two bins is the quadratic sum of uncorrelated uncertainties.
In the high-p µµµ T region, the B + c shows a moderate suppression, while the low-p µµµ T modification factor stands above unity and above the high-p µµµ T region, respectively, by 1.2 and 1.8 standard deviations, consistent with an enhancement of the integrated production and a softening of the p T spectrum in the QGP. No significant variation is observed as a function of centrality. As , other heavy mesons in these p T ranges typically show more suppression than our measurement [4,7,13,[15][16][17], which may indicate that heavy-quark recombination is a significant B + c production mechanism in the QGP. A study based on Ref. [46], ignoring the recombination of B + c excited states and possibly underestimating initial correlations, predicts an R AA (B + c ) about an order of magnitude smaller than our measurement.
In summary, the first observation of the B + c meson in heavy ion collisions is presented, using the B + c → (J/ψ → µ + µ − )µ + ν µ decay. The production cross sections in lead-lead and protonproton collisions and the nuclear modification factor derived from their ratio are measured in two bins of the trimuon transverse momentum, and in two ranges of the heavy-ion centrality. This unique bottom-charm state can help disentangle the enhancement (possibly dominant in central events at low p T ) and suppression (dominant at high p T ) mechanisms at play in the evolution of heavy quarks through the quark-gluon plasma.

Acknowledgments
We congratulate our colleagues in the CERN accelerator departments for the excellent performance of the LHC and thank the technical and administrative staffs at CERN and at other CMS institutes for their contributions to the success of the CMS effort. In addition, we gratefully acknowledge the computing centers and personnel of the Worldwide LHC Computing Grid and other centers for delivering so effectively the computing infrastructure essential to our analyses. Finally, we acknowledge the enduring support for the construction and operation of the LHC, the CMS detector, and the supporting computing infrastructure provided by the follow-  [23] CMS Collaboration, "Measurement of the ratio of the production cross sections times branching fractions of B ± c → J/ψπ ± and B ± → J/ψK ± and [28] "HEPData record for this analysis", 2021. doi:10.17182/hepdata.111309.