Transverse momentum distributions of charmonium states with the statistical hadronization model

TL;DR

This paper uses the statistical hadronization model to predict transverse momentum spectra of charmonium states in heavy-ion collisions, showing good agreement with data and providing new predictions for the exotic state X(3872).

Contribution

It introduces detailed predictions for charmonium transverse momentum distributions, including the exotic X(3872), within the statistical hadronization framework at LHC energies.

Findings

Good agreement with existing data for J/ψ and ψ(2S) spectra.

Predicted shape and yield for X(3872) in Pb-Pb collisions.

Enhanced low transverse momentum production for X(3872).

Abstract

Calculations and predictions are presented within the framework of the statistical hadronization model for transverse momentum spectra of the charmonium states J/$\psi$, $\psi(2S)$ and $X(3872)$ produced in nucleus-nucleus collisions at LHC energies. The results are confronted with available data and exhibit very good agreement by using particle flow profiles from state-of-the-art hydrodynamic calculations. For $X(3872)$ production in Pb-Pb collisions we predict a transverse momentum distribution similar in shape to that for J/$\psi$ with a strong enhancement at low transverse momenta and a production yield of about 1\% relative to that for J/$\psi$.

Figures

Equations

$$N_{\mathrm{c}\bar{\mathrm{c}}}=\frac{1}{2}$$

$$+$$

$$\frac{{\mathrm{d}}^{2}N}{p_{\mathrm{T}}\mathrm{d}{p}_{\mathrm{T}}\mathrm{d}y}$$

$$m_{\mathrm{T}}$$

$$-p_{\mathrm{T}}$$

$$f(p_{\mathrm{T}})=C\times \frac{p_{\mathrm{T}}}{(1+(p_{\mathrm{T}}/p_0{)}^2{)}^n},$$

Full text

1 Introduction

Ultra-relativistic heavy-ion collisions are widely used to investigate the evolution from the hot and dense phase of QCD with quarks and gluons as degrees of freedom, the quark-gluon plasma (QGP), towards hadronization into color-neutral objects as degrees of freedom. Analysis of the abundances of the resulting hadrons provides a reliable tool to characterize the phase boundary in the phase diagram of the strong interaction [1].

The statistical hadronization model (SHM) is successfully used to predict and describe hadron abundances produced in relativistic nuclear collisions [2]. The underlying assumption is that at the latest at hadronization the fireball formed in such collisions is close to thermal equilibrium such that hadron yields can be characterized by a grand canonical partition function where baryon number, the $3$-component of the isospin, and strangeness are conserved on average and a rapid hadrochemical freeze-out takes place at the phase boundary. These assumptions connect the microscopic content of a heavy-ion collision with macroscopic attributes like the chemical freeze-out temperature $T_{\rm CF}$ and the baryochemical potential $\mu_{\rm B}$.

In the light-flavor sector it is observed that, besides the frequently produced mesons and baryons, also more complex and loosely bound objects, such as light nuclei and even hypernuclei, along with their anti-particles, can be described within the SHM [3, 2]. In Figure 1, the primordial and total (anti-)particle yields computed within the SHM are shown as a function of the particle mass $m$ for a collision energy of $\sqrt{s_{\rm NN}}=2.76$ TeV111It should be noted, that throughout this letter, we use $c=k_{\rm B}=\hbar=1$.. For $m\gg T_{\rm CF}$, the primordial or thermal yields follow the relation $m^{3/2}\exp(-m/T_{\rm CF})$. The total yields include contributions from resonance decays and reproduce particle yields from data over nine orders of magnitude from the lightest mesons up to nuclei with mass number $A=4$.

A striking feature in the figure is the difference between the light-flavored loosely bound state hyper-triton ${}^{3}_{\Lambda}$H and the charmonium state ${\rm J}/\psi$. Although their masses are very similar, the yield of ${\rm J}/\psi$ is about three orders of magnitude larger. The origin of this enhancement is due to the statistical hadronization of charm quarks, formed in hard initial processes, where the number of charm quarks is conserved throughout the evolution of the collision. Technically, the conservation of charm quarks leads to a fugacity in the SHM for charmed hadrons [4] which is, however, not a free parameter but determined by the measured charm cross section. Charm quarks are not confined inside the QGP, thermalize within the QGP and hadronize at the QCD phase boundary into open and hidden charm hadrons. This SHM for charmed hadrons (SHMC) provides an excellent description of charmonium production [5, 6, 7] without any new parameters and represents compelling evidence, as demonstrated in Figure 1, for this new production mechanism. A more detailed account of the SHMC is given below.

Furthermore, a large degree of thermalization is observed in the spectra and the elliptic flow of $D$-mesons and their decay electrons [8, 9]. A number of recent measurements have established the SHMC process (sometimes dubbed ‘(re)generation’) as the dominant production mechanism of ${\rm J}/\psi$ in heavy-ion collisions at LHC energies [10, 11, 12, 13]. It is therefore appealing and important to extend the intriguing results of ${\rm J}/\psi$ production beyond yields to particle spectra and to more complex charmonium as well as open charm states to further investigate the SHMC mechanism. In the present publication we focus on charmonium states. Predictions for the open charm sector will be the subject of a future publication.

Loosely bound states such as $\psi(2\rm S)$ and, more dramatically, the potential tetra-quark charmonium state ${\rm X}(3872)$ are of particular interest. Its observation by the Belle collaboration [14] and the subsequent confirmation by the CDF [15], D0 [16], and BaBar [17] collaborations showed that it is a narrow charmonium-like resonance and the close vicinity of the particle mass to the $D^{0}\bar{D}^{*0}$ production threshold suggests that the particle could be a charm meson molecule with a very small binding energy [18]. At the LHC, the state ${\rm X}(3872)$ was first observed by the LHCb collaboration [19], which later also determined [20, 21] its quantum numbers, $J^{PC}=1^{++}$.

In addition to the ${\rm X}(3872)$, the $\psi(2\rm S)$ is a natural choice when expanding the studies on the SHMC mechanism beyond the ${\rm J}/\psi$.

In this letter, we present calculations for the yields and transverse momentum spectra of the charmonium states ${\rm J}/\psi$, $\psi(2\rm S)$, and ${\rm X}(3872)$ for heavy-ion collisions at $\sqrt{s_{\rm NN}}=5$ TeV. Results will be presented for the current collision system Pb–Pb and, in the case of ${\rm X}(3872)$, also for Kr–Kr collisions where much larger luminosities are possible at the LHC.

2 Heavy quarks in the statistical hadronization model

In the SHMC it is assumed that charm quarks222In this letter we focus on charmonium, although the studies can be extended to any hadron species including charm and/or beauty content. are produced in initial hard scatterings and that during the QGP phase the number of (anti-)charm quarks is conserved, i.e. the thermal production or annihilation are negligible at LHC energies [5]. The color charge of charm quarks is screened by the color-dense medium for $T\gtrsim T_{\rm CF}=156.5$ MeV and they do not form color-less bound states in the fireball volume $V$, where $\mu_{\rm B}$ is consistent with zero for LHC energies, as determined by thermal fits, see [2]. The quarks thermalize in the QGP before the hadronization and rapid freeze-out at the phase boundary.

The (anti-)charm hadron densities computed in canonical statistical mechanics, $n_{X}^{{\rm th}}$, are anchored to the number of produced $\rm c\bar{c}$ pairs, $N_{\rm c\bar{c}}$, by a balance equation

[TABLE]

where the quantity $N_{\rm c\bar{c}}$ is interpolated via FONLL [22, 23] from charm cross-section measurements [24, 25, 26, 27] in the corresponding rapidity region. Shadowing is taken into account when calculating $N_{\rm c\bar{c}}$ for nucleus-nucleus collisions using rapidity-dependent measurements of the nuclear modification factor in proton-nucleus collisions of $D$-mesons [28], where interpolations, if necessary, are done via model calculations [29, 30, 31]. A canonical suppression factor, $I_{1}(g_{c}n_{{\rm oc}}^{{\rm th}}V)/I_{0}(g_{c}n_{{\rm oc}}^{{\rm th}}V)$, is applied to the open charm densities computed in grand canonical statistical mechanics to obtain $n_{X}^{\rm th}$ in eq. 1. Here, $I_{i}$ with $i=\{0,1\}$ are modified Bessel functions and $n^{{\rm th}}_{{\rm oc}}$ are the thermal open charm meson densities. This correction factor gains importance for decreasing charm quark densities, i.e. $N_{\rm c\bar{c}}\lesssim 1$. In order to fulfill the balance equation (eq. 1), the charm quark fugacity $g_{c}$ is needed. As explained before, it is not a free parameter but entirely defined by $N_{\rm c\bar{c}}$ and $n_{X}^{\rm th}$. The parameter $V$ used in the equations above corresponds to the volume of one unit of rapidity and is determined in the statistical model fit of non-charm hadron yields [2].

The density profile of the colliding nuclear species is taken into account by a core-corona approach. The ‘core’ part represents the central region of the colliding nuclei where nucleons undergo many scatterings and are assumed to produce a QGP. The core fraction is normalized to the thermal yields. The ‘corona’ area includes the nucleons in the outer region of the colliding nuclei. In the overlap zone with on average one or less collisions, no QGP formation is assumed. Rather the yields from the corona fraction are modeled by proton-proton ($\rm pp$) differential distributions scaled by the number of binary nucleon-nucleon collisions in the corona. The radius at which the nucleon density is not sufficient anymore to produce more than one inelastic scattering is taken from the nuclear charge density distributions. This density is found to be approximately at $10$ % of the central nuclear density. To get a feeling for the sensitivity to this estimate, we also give the result for a value of $20$ %. The fraction of the core and the corona part is estimated by a Glauber simulation [32].

The calculation of the transverse momentum spectra is based on the following consideration: the charm quarks are assumed in local thermal equilibrium in the fireball formed in the collision. At hadronization, i.e. at $T_{\rm CF}$, the charmonia states then ‘inherit’ the random thermal motion of the charm quarks superimposed with the collective velocity of the expanding QGP fluid. This velocity is, hence, extracted from state-of-the-art viscous hydrodynamic modelling of light flavor observables.

Specifically, this modelling is based on the (3+1)D viscous hydrodynamic simulation framework MUSIC [33] with IP-Glasma as initial condition333Here it is necessary to discriminate between the purpose of IP-Glasma and the previously mentioned shadowing. IP-Glasma affects the shape of the thermal contribution in the $p_{\rm T}$-differential spectrum whereas shadowing affects the $p_{\rm T}$-integrated $\rm c\bar{c}$ cross section used to calculate the total number of produced $\rm c\bar{c}$ pairs $N_{\rm c\bar{c}}$ in a heavy-ion collision used in eq. 1. [34]. This framework does not assume boost-invariance. Rather the relevant parameters are adjusted to reproduce the charged particle multiplicity distributions measured by the ALICE collaboration as function of the pseudo-rapidity $\eta$ [35, 36]. The equation of state and the shear viscosity are parameterized from lattice QCD calculations [37]. The freeze-out hyper-surface is then taken at the chemical freeze-out temperature $T=T_{\rm CF}$.

The results from these hydrodynamic simulations for the velocity profile $n$ and the transverse flow velocity $\beta^{\rm s}_{\rm T}$ are fed into a blast-wave function formulated for a Hubble-type expansion of the hyper-surface in a boost-invariant scenario as relevant for mid-rapidity data at LHC energy. The result, inspired by earlier calculations reported in [38, 39], is:

[TABLE]

where $\rho={\rm atanh}(\beta^{\rm s}_{\rm T}(r/R)^{n})$ and $I_{i}$ and $K_{i}$ are modified Bessel functions with $i=\{0,1\}$. The transverse mass $m_{\rm T}=\sqrt{m^{2}+p_{\rm T}^{2}}$ is obtained using the corresponding particle mass and the temperature is that of chemical freeze-out $T=T_{\rm CF}$. This functional form (eq. 2) is then used to model the thermal part of the spectrum and is normalized to the core fraction. To implement the boost-noninvariance discussed above in an approximate way for the description of data at forward rapidity, we adjust the flow velocity according to the (3+1)D hydrodynamic calculations, keeping otherwise the same functional form of the blast-wave description.

We have tested numerically that using eq. 2 for large masses such as that of ${\rm J}/\psi$ yields results very close to those obtained with the standard blast wave description of [38]. The overall normalization factor is obtained from the SHMC as used in Figure 1.

The shape of the corona fraction is modeled by a fit to the transverse momentum spectra measured in $\rm pp$ collision, given by

[TABLE]

where $C$, $p_{0}$ and $n$ are fit parameters.

The recent measurement of ALICE is being used [40] for the $p_{\rm T}$-differential cross section of ${\rm J}/\psi$ in $\rm pp$ collisions for mid-rapidity at $\sqrt{s_{\rm NN}}=5$ TeV.

In case of the $\psi(2\rm S)$, the $\rm pp$ cross section is taken from measurements at the corresponding collision energy and rapidity [41]. The transverse momentum spectrum of $\psi(2\rm S)$ in $\rm pp$ collisions at low $p_{\rm T}$ is not measured at mid- or forward rapidity. The shape is approximated by the ${\rm J}/\psi$ transverse momentum spectrum from [41], which is consistent with the finding that the production ratio of $\psi(2\rm S)$ and ${\rm J}/\psi$ are constant over a wide range in transverse momentum [42]. Small differences can be expected for low transverse momentum [43, 41], which should be negligible due to the dominant thermal contribution.

The ${\rm X}(3872)$ $p_{\rm T}$-differential cross section is measured at mid-rapidity for $10<p_{\rm T}/{\rm GeV}<30$ in $\rm pp$ collisions at $\sqrt{s}=7$ TeV [44]. A universal collision energy scaling is assumed for the transverse momentum spectrum to extrapolate the spectrum from $\sqrt{s}=7$ TeV to $5$ TeV for an estimate of the corona shape of the ${\rm X}(3872)$ spectrum. ‘Hybrid’ ${\rm J}/\psi$ spectra in $\rm pp$ collisions are created from the ALICE and CMS results at $5$ and $7$ TeV to cover the low- and high-$p_{\rm T}$ region simultaneously [40, 45, 46, 47]. These spectra are fitted by the functional form eq. 3 and the ratio of the transverse momentum spectra parameterizations is used to weight the ${\rm X}(3872)$ $p_{\rm T}$ spectrum at $\sqrt{s}=7$ TeV to estimate the shape of the spectrum at $\sqrt{s}=5$ TeV. Finally, the ${\rm X}(3872)$ cross section in $\rm pp$ collisions used for the determination of the corona fraction is estimated by an extrapolation from the $p_{\rm T}$ range measured by CMS down to zero $p_{\rm T}$ using eq. 3 and the corresponding $p_{\rm T}$ spectrum at $\sqrt{s}=5$ TeV.

3 Results

In Figure 2, the calculations are shown for the transverse momentum spectrum of the ${\rm J}/\psi$ at mid-rapidity for Pb–Pb collisions at $\sqrt{s_{\rm NN}}=5$ TeV. In addition to the full spectrum calculations the two components, thermal core fraction and corona, are plotted to visualize their contributions. One can see how dramatically different in shape the predicted $p_{\rm T}$ spectrum in the hadronizing QGP is due to the large fugacity, see eq. 1 and the strong collective transverse expansion. In minimum bias $\rm pp$ collisions such effects are expected to be very small.

In the upper panel of Figure 3, calculations for the transverse momentum spectrum of ${\rm J}/\psi$ in Pb–Pb collisions at a collision energy of $\sqrt{s_{\rm NN}}=5$ TeV are compared to the measurement by the ALICE collaboration at forward rapidity [11]. The model describes the data for low $p_{\rm T}$ values but falls below the data for $p_{\rm T}\gtrsim 4.5$ GeV. This suggests additional production mechanisms such as ${\rm J}/\psi$ production from gluon fragmentation in jets which could contribute towards increasing transverse momentum. In addition, predictions for the $\psi(2\rm S)$ are shown for the same collision system and rapidity window for central collisions. The lower panel shows the production ratio of $\psi(2\rm S)/{\rm J}/\psi$.

In Figure 4, predictions for the ${\rm X}(3872)$ transverse momentum spectrum at $\sqrt{s_{\rm NN}}=5$ TeV are shown for Pb–Pb and Kr–Kr collisions. The branching ratio ${\rm BR}({\rm X}(3872)\to{\rm J}/\psi\pi^{+}\pi^{-})$ is only known approximately, but can be constrained to a range between $3.2$% and about $20$% [48]. All results are therefore shown as ${\rm BR}({\rm X}(3872)\to{\rm J}/\psi\pi^{+}\pi^{-})\times{\rm d}^{2}\sigma/{\rm d}p_{\rm T}{\rm d}y$, where the branching ratio is chosen to $10$%. Heavy-ion data taking in LHC Run4 or beyond the long shutdown 4 of the LHC, i.e. after the year $2030$, might enable these measurements. We also studied ${\rm X}(3872)$ production in collisions between lighter beams such as ${}^{84}{\rm Kr}$ nuclei to assess, in future experiments, the impact of potential increased luminosity vs reduced cross section. In this case, the yields in the Kr–Kr collisions are found to be a factor $4-5$ lower than for Pb–Pb collisions for $p_{\rm T}\lesssim 5$ GeV.

4 Summary and conclusions

Calculations and predictions using the SHMC for the transverse momentum spectra of different charmonium states are presented and, whenever measurements are available, confronted with data. The spectra calculations are based on the yields obtained from the SHMC as well as a modified blast wave function with input from state-of-the-art hydrodynamical calculations for the flow profiles. Without any further modifications or parameters a very good description of the measured ${\rm J}/\psi$ transverse momentum spectra in the region $p_{\rm T}\lesssim 5$ GeV is obtained. These results further strengthen the evidence that charmonia are produced from deconfined charm quarks at the QCD phase boundary as their yields and spectra are established at the LQCD chiral crossover temperature.

The comparison of the calculations with data for ${\rm J}/\psi$ in Pb–Pb collisions at $\sqrt{s_{\rm NN}}=5$ TeV at forward rapidity shows a good agreement for low $p_{\rm T}$ but falls below the data for increasing momenta ($p_{\rm T}\gtrsim 5$ GeV). This suggests additional production mechanisms beyond (re)generation at high charmonium momenta.

Predictions are presented for $\psi(2\rm S)$ production in Pb–Pb collisions at $\sqrt{s_{\rm NN}}=5$ TeV at forward rapidity. Statistics collected in LHC heavy-ion data taking in LHC Run 3 might be sufficient for a comparison with data. These results will be important to understand open issues such as the possible presence, in the QGP, of color-less bound states. We further presented calculations for the transverse momentum spectrum of ${\rm J}/\psi$ mesons produced in Pb-Pb collisions at mid-rapidity which can be tested soon.

Predictions for the production of ${\rm X}(3872)$ for different centrality regimes of nuclear collisions at $\sqrt{s_{\rm NN}}=5$ TeV are presented for Pb–Pb and Kr–Kr systems. The measurements might become available after the long shutdown 3 of the LHC and will provide a deeper understanding of the mechanism of charmonium and exotica production at the QGP phase boundary and consequently about the QCD dynamics of deconfinement and the hadronization process.

Acknowledgment

We thank Andrea Dubla, Klaus Reygers and Chun Shen for fruitful discussions.

This work is part of and supported by the DFG Collaborative Research Centre “SFB 1225 (ISOQUANT)”. K.R. acknowledges support by the Polish National Science Centre under Maestro grant DEC-2013/10/A/ST2/00106 and OPUS 2018/31/B/ST2/01663. K.R. also acknowledges partial support of the Polish Ministry of Science and Higher Education.