# Quarkonium production in high energy *pp* collisions

Jiaxing Zhao, Pol Bernard Gossiaux, Taesoo Song, Elena Bratkovskaya, Joerg Aichelin

► To cite this version:

Jiaxing Zhao, Pol Bernard Gossiaux, Taesoo Song, Elena Bratkovskaya, Joerg Aichelin. Quarkonium production in high energy *pp* collisions. 30th International Conference on Ultrarelativstic Nucleus-Nucleus Collisions, Sep 2023, Houston, United States. pp.09014, 10.1051/epjconf/202429609014. hal-04381982

HAL Id: hal-04381982

https://hal.science/hal-04381982v1

Submitted on 3 Jun 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Quarkonium production in high energy pp collisions

Jiaxing Zhao¹, Pol Bernard Gossiaux¹, Taesoo Song², Elena Bratkovskaya²,³, and Jörg Aichelin¹,*

¹SUBATECH, Nantes University, IMT Atlantique, IN2P3/CNRS 4 rue Alfred Kastler, 44307 Nantes cedex 3, France

²GSI Helmholtzzentrum für Schwerionenforschung GmbH, Planckstr. 1, 64291 Darmstadt, Germany

³Helmholtz Research Academy Hesse for FAIR (HFHF), GSI Helmholtz Center for Heavy Ion Physics, Campus Frankfurt, 60438 Frankfurt, Germany

**Abstract.** We investigate the charmonium and bottomonium production in pp collisions using the Wigner densities formalism. The Wigner density of the quarkonia is approximated by analytical 3-D isotropic harmonic oscillator Wigner densities with the same root-mean-square radius given by the solution of the Schrödinger equation. This approach reproduces quite well the available experimental transverse momentum and rapidity distributions.

Hidden heavy flavour mesons are a useful tool to study the strongly interacting quark gluon plasma, which is created during high energy heavy-ion collisions. Due to the large quark mass $m_Q$, $m_Q \gg \Lambda_{QCD}$, with $\Lambda_{QCD}$ being the QCD cutoff, quarkonium production can be factorized into the production of a heavy quark, $Q\bar{Q}$, pair, which can be described by perturbative QCD and a subsequent soft non-perturbative process, which describes the formation of a colorless quarkonium from the $Q\bar{Q}$ pair. For the latter part many approaches have been advanced. They include the Color-Evaporation Model, the Color-Singlet Model and the Color-Octet Model. The latter two are encompassed in the NRQCD approach. For a review we refer to [1]. More recently, the rapidity and $p_T$ distributions of hidden heavy flavour mesons, produced in pp collisions, have also been well reproduced in the Wigner density matrix formalism [2–4]. Here we will extend this formalism up to 3S and use EPOS4 to generate the initial heavy quarks.

The Wigner density formalism is based on the quantal density matrix projection in which the probability that a meson $i$ is produced is given by $P_i = \text{Tr}(\rho_i\rho^{(N)})$ with $\rho_i$ being the density matrix of the meson $i$ and $\rho^{(N)}$ the density matrix of the N heavy quarks and antiquarks, produced in a pp collision. A partial Fourier transformation of the density matrices yields then

$$ \frac{dP_i}{d^3\mathbf{R}d^3\mathbf{P}} = \sum \int \frac{d^3\mathbf{r}d^3\mathbf{p}}{(2\pi)^6} W_i(\mathbf{r}, \mathbf{p}) \prod_{j>2} \int \frac{d^3\mathbf{r}_j d^3\mathbf{p}_j}{(2\pi)^{3(N-2)}} W^{(N)}(\mathbf{r}_1, \mathbf{p}_1, \mathbf{r}_2, \mathbf{p}_2, \dots, \mathbf{r}_N, \mathbf{p}_N). \quad (1) $$

$W_i$ is the two-body Wigner density of the bound heavy quark pair and $W^{(N)}(\mathbf{r}_1, \mathbf{p}_1, \mathbf{r}_2, \mathbf{p}_2, \dots, \mathbf{r}_N, \mathbf{p}_N)$ is the quantal density matrix in Wigner representation of the ensemble of N heavy quarks produced in a pp collision. $\mathbf{r}(\mathbf{R})$ and $\mathbf{p}(\mathbf{P})$ are the relative (center of mass) coordinate and momentum of the heavy quark and antiquark, which are

<footnote>*Speaker</footnote>
<table><caption>**Table 1.** The masses, root-mean-square radius, and the Gaussian widths $\sigma$ of different charmonium and bottomonium states in vacuum. The experimental data are from Ref. [7].</caption><thead><tr><th></th><th>J/ψ</th><th>χ<sub>c</sub>(1P)</th><th>ψ(2S)</th><th>Υ(1S)</th><th>χ<sub>b</sub>(1P)</th><th>χ<sub>b</sub>(1D)</th><th>Υ(2S)</th><th>χ<sub>b</sub>(2P)</th><th>Υ(3S)</th></tr></thead><tbody><tr><td>Mass Theo.(GeV)</td><td>3.071</td><td>3.483</td><td>3.652</td><td>9.390</td><td>9.870</td><td>10.109</td><td>9.959</td><td>10.208</td><td>10.288</td></tr><tr><td>Mass Exp. (GeV)</td><td>3.097</td><td>3.463</td><td>3.686</td><td>9.460</td><td>9.876</td><td>10.163</td><td>10.023</td><td>10.243</td><td>10.355</td></tr><tr><td>⟨r<sup>2</sup>⟩ (fm<sup>2</sup>)</td><td>0.182</td><td>0.453</td><td>0.714</td><td>0.042</td><td>0.153</td><td>0.284</td><td>0.236</td><td>0.410</td><td>0.520</td></tr><tr><td>σ(fm)</td><td>0.348</td><td>0.426</td><td>0.452</td><td>0.167</td><td>0.247</td><td>0.285</td><td>0.260</td><td>0.302</td><td>0.307</td></tr></tbody></table>

bound in a quarkonium. We assume that the unknown quantal N-body Wigner density can be replaced by the average of classical phase space distributions, $W^{(N)} \approx \langle W_{\text{classical}}^{(N)} \rangle$.

The classical momentum space distributions of the heavy quarks is provided by EPOS4 [5, 6]. EPOS4, however, provides only the coordinate information of the vertex where the $Q\bar{Q}$ pair is created. For a heavy quark pair, created at the same vertex, we assume that the relative distance between $Q$ and $\bar{Q}$ in their center-of-mass frame is given by a Gaussian distribution and the Wigner density for a $Q\bar{Q}$ pair can be expressed as

$$W^{(2)}(\mathbf{r}, \mathbf{p}) \sim r^2 \exp\left(-\frac{r^2}{2\sigma_{Q\bar{Q}}^2}\right) f_{Q\bar{Q}}^{\text{EPOS4}}(\mathbf{p}), \qquad (2)$$

where the distance is controlled by the effective width $\sigma_{Q\bar{Q}}$. Having now momenta and positions of the heavy quarks we can calculate the yield of charmonium and bottomonium via Eq. (1).

<figure><img src="image_2.png" /><figcaption>**Figure 1.** Wave function of different charmonium (left) and bottomonium (right) states. Solid lines are from the Schrödinger equation, dashed lines are from the 3-D isotropic harmonic oscillator (Eq. (3)).</figcaption></figure>

We come now to the construction of the quarkonium Wigner density. The Wigner density is obtained by the Wigner-Weyl transformation of the density matrix of the quarkonia. The quarkonium wave function is the solution of the two-body Schrödinger equation, which we solve for charmonium and bottomonium with the Cornell potential, $V(r) = -\alpha/r + \kappa r + c$ with $\alpha = 0.513$, $\kappa = 0.17\text{GeV}^2$, $c = -0.161$, and with the quark masses $m_c = 1.5\text{GeV}$ and $m_b = 5.2\text{GeV}$. The wave functions are shown in Fig. 1, and the masses and root-mean-square radii $\langle r^2 \rangle$ are shown in Table 1. We can see that the masses are very close to the experimental values. The wave function is, however, not analytical, and therefore the Wigner
density can only be calculated numerically, what makes the solution of Eq. (1) complicated. It is therefore convenient to approximate the wave function by a 3-D isotropic harmonic oscillator wave function for the potential $V(r) = 1/(2m_Q\sigma^4)r^2$. These wave functions are analytical and can be expressed as $\psi_{nlm}(r, \theta, \phi) = R_{nl}(r)Y_{lm}(\theta, \phi)$, where $Y_{lm}$ are the spherical harmonics. The radial part can be expressed as,

$$R_{nl}(r) = \left[ \frac{2(n!)}{\sigma^3 \Gamma(n+l+3/2)} \right]^{1/2} \left(\frac{r}{\sigma}\right)^l e^{-\frac{r^2}{2\sigma^2}} L_n^{l+1/2}\left(\frac{r^2}{\sigma^2}\right), \quad (3)$$

where $L_n^{l+1/2}$ are Laguerre polynomials. The parameters of the 3-D isotropic harmonic oscillator wave functions are chosen to match the root-mean-square radius $\langle r^2 \rangle$ of the real quarkonium wave function: $\langle r^2 \rangle = 3\sigma^2/2$ for 1S, $\langle r^2 \rangle = 5\sigma^2/2$ for 1P, $\langle r^2 \rangle = 7\sigma^2/2$ for 1D and 2S, $\langle r^2 \rangle = 9\sigma^2/2$ for 2P, and $\langle r^2 \rangle = 11\sigma^2/2$ for 3S states. The corresponding widths are shown in Table 1 and the wave functions are shown in Fig. 1 with dashed lines. We can see that the ground states and low lying excited states can be well reproduced by the 3-D isotropic harmonic oscillator, while the difference increases for higher excited states, e.g. 2S, 2P, and 3S.

The Wigner densities for different states up to 3S are

$$W_{1S}(\mathbf{r}, \mathbf{p}) = 8e^{-\xi}, \quad (4)$$

$$W_{1P}(\mathbf{r}, \mathbf{p}) = \frac{8}{3}e^{-\xi}(2\xi - 3),$$

$$W_{1D}(\mathbf{r}, \mathbf{p}) = \frac{8}{15}e^{-\xi}(15 + 4\xi^2 - 20\xi + 8[p^2r^2 - (\mathbf{p} \cdot \mathbf{r})^2]),$$

$$W_{2S}(\mathbf{r}, \mathbf{p}) = \frac{8}{3}e^{-\xi}(3 + 2\xi^2 - 4\xi - 8[p^2r^2 - (\mathbf{p} \cdot \mathbf{r})^2]),$$

$$W_{2P}(\mathbf{r}, \mathbf{p}) = \frac{8}{15}e^{-\xi}(-15 + 4\xi^3 - 22\xi^2 + 30\xi - 8(2\xi - 7)[p^2r^2 - (\mathbf{p} \cdot \mathbf{r})^2]),$$

$$W_{3S}(\mathbf{r}, \mathbf{p}) = \frac{8}{315}e^{-\xi}(315 + 42\xi^4 - 336\xi^3 + 924\xi^2 - 840\xi - [2009 + 32p^2r^2 \\ + 336r^4/\sigma^4 - 1400r^2/\sigma^2 - 896p^2\sigma^2 + 224p^4\sigma^4][p^2r^2 - (\mathbf{p} \cdot \mathbf{r})^2] - [686 \\ + 608p^2r^2 + 112r^2/\sigma^2 - 896p^2\sigma^2 + 224p^4\sigma^4 - 672(\mathbf{p} \cdot \mathbf{r})^2](\mathbf{p} \cdot \mathbf{r})^2),$$

where $\xi = r^2/\sigma^2 + p^2\sigma^2$. The transverse momentum and rapidity distribution of charmonium

<figure><img src="image_3.png" /><figcaption>**Figure 2.** $p_T$ spectra of different charmonium states (left), where red is $J/\psi$, orange is $\chi_c(1P)$, and green is $\psi(2S)$. Prompt $J/\psi$ (middle). Prompt $\psi(2S)$ (right). The experimental data are from ALICE [8, 9], ATLAS [10], and CMS [11].</figcaption></figure>

and bottomonium (see Eq. (1)) are shown in Figs. 2 and 3. We see, as far as data are available,
<figure><img src="image_4.png" /><figcaption>**Figure 3.** $p_T$ and rapidity $y$ dependence of different bottomonium states. The experimental data are from CMS [12].</figcaption></figure>

a quite good agreement with the experimental results for $c\bar{c}$ as well $b\bar{b}$ mesons when choosing $\sigma_{c\bar{c}} = 0.4$ fm and $\sigma_{b\bar{b}} = 0.2$ fm in Eq. (2).

We can conclude that the experimentally available rapidity and transverse momentum distribution of $c\bar{c}$ and $b\bar{b}$ quarkonia can be well described in the Wigner density formalism. The only parameter which enters the calculation is the width of the distribution of the relative distance of the $Q\bar{Q}$ pair at production. The relative contribution of the different states is then exclusively given by their wave function.

*Acknowledgements:* This work is funded by the European Union's Horizon 2020 research and innovation program under grant agreement No. 824093 (STRONG-2020). T.S. and E.B. acknowledge support by the Deutsche Forschungsgemeinschaft through the grant CRC-TR 211, Project N 315477589 - TRR 211.

## References

[1] A. Andronic et al., Eur. Phys. J. C **76**, 107 (2016), 1506.03981

[2] T. Song, J. Aichelin, E. Bratkovskaya, Phys. Rev. C **96**, 014907 (2017), 1705.00046

[3] D.Y.A. Villar, J. Zhao, J. Aichelin, P.B. Gossiaux, Phys. Rev. C **107**, 054913 (2023), 2206.01308

[4] T. Song, J. Aichelin, J. Zhao, P.B. Gossiaux, E. Bratkovskaya, Phys. Rev. C **108**, 054908 (2023), 2305.10750

[5] K. Werner, Phys. Rev. C **108**, 064903 (2023), 2301.12517

[6] K. Werner, B. Guiot, Phys. Rev. C **108**, 034904 (2023), 2306.02396

[7] R.L. Workman, Others (Particle Data Group), PTEP **2022**, 083C01 (2022)

[8] S. Acharya et al. (ALICE), JHEP **03**, 190 (2022), 2108.02523

[9] S. Acharya et al. (ALICE), JHEP **10**, 084 (2019), 1905.07211

[10] M. Aaboud et al. (ATLAS), Eur. Phys. J. C **78**, 171 (2018), 1709.03089

[11] A.M. Sirunyan et al. (CMS), Phys. Lett. B **790**, 509 (2019), 1805.02248

[12] A.M. Sirunyan et al. (CMS), Phys. Lett. B **790**, 270 (2019), 1805.09215