# Extracting the jet transport coefficient from hadron suppressions by
confronting current NLO parton fragmentation functions

Qing-Fei Han  Key Laboratory of Quark and Lepton Physics (MOE) and Institute
of Particle Physics, Central China Normal University, Wuhan 430079, China  Man
Xie  Key Laboratory of Quark and Lepton Physics (MOE) and Institute of
Particle Physics, Central China Normal University, Wuhan 430079, China  Han-
Zhong Zhang  Key Laboratory of Quark and Lepton Physics (MOE) and Institute of
Particle Physics, Central China Normal University, Wuhan 430079, China
Guangdong Provincial Key Laboratory of Nuclear Science, Institute of Quantum
Matter, South China Normal University, Guangzhou 510006, China  Guangdong-Hong
Kong Joint Laboratory of Quantum Matter, Southern Nuclear Science Computing
Center, South China Normal University, Guangzhou 510006, China

###### Abstract

Nuclear modification factors of single hadrons and dihadrons at large
transverse momentum ($p_{\rm T}$) in high-energy heavy-ion collisions are
studied in a next-to-leading-order (NLO) perturbative QCD parton model. Parton
fragmentation functions (FFs) in $A+A$ collisions are modified due to jet
energy loss which is proportional to the jet transport coefficient $\hat{q}$
characterizing the interaction between the parton jet and the produced medium.
By confronting 6 current sets of NLO parton FFs for large $p_{\rm T}$ hadron
productions, we extract $\hat{q}$ quantitatively via a global fit to data for
both single hadron and dihadron suppressions, and obtain
$\hat{q}/T^{3}=4.74-6.72$ at $T=370$ MeV in central $Au+Au$ collisions at
$\sqrt{s_{\rm NN}}=200$ GeV, and $\hat{q}/T^{3}=3.07-3.98$ at $T=480$ MeV in
central $Pb+Pb$ collisions at $\sqrt{s_{\rm NN}}=2.76$ TeV. The numerical
results show that the uncertainties for $\hat{q}$ extraction are brought by
the different contributions of gluon-to-hadron in the 6 sets of FFs due to
gluon energy loss being $9/4$ times of quark energy loss.

##  I Introduction

The hot and dense quark-gluon plasma (QGP) could be produced in high-energy
heavy-ion collisions performed at the Relativistic Heavy Ion Collider (RHIC)
Adare  _et al._ (2008, 2013, 2010) and the Large Hadron Collider (LHC) Abelev
_et al._ (2013); Chatrchyan  _et al._ (2012); Acharya  _et al._ (2018);
Khachatryan  _et al._ (2017); Acharya  _et al._ (2019). When hard partons
generated in the initial stage of heavy-ion collisions traverse such a nuclear
medium before fragmenting into hadrons, they will lose most of their energy
due to multiple scatterings with the strongly-interacting medium. As a
consequence, the final state hadrons at large transverse momentum $p_{\rm T}$
are suppressed. This phenomenon is known as Jet Quenching Gyulassy and Plumer
(1990); Wang and Gyulassy (1992); Qin and Wang (2015), which as a hard probe
plays an essential role in the studies of QGP properties. The nuclear
modification factor, $R_{AA}(p_{\rm T})$ for single hadrons or $I_{AA}(p_{\rm
T})$ for dihadrons, is a widely-used and appropriate observable to reflect the
suppression effect quantitatively, and is defined as the ratio of the hadron
spectrum in nucleus-nucleus ($A+A$) collisions to that in proton-proton
($p+p$) collisions normalized by the number of binary nucleon-nucleon
collisions Zhang  _et al._ (2007).

In recent years, quantitative researches have become the mainstream of a large
number of theoretical and experimental studies on jet quenching. One of the
important parts is to extract the jet transport coefficient $\hat{q}$
associated with energy loss by comparing the theoretical results with the
experimental data. The parameter $\hat{q}$ characterizes the average
transverse momentum broadening squared per unit length for a jet propagating
inside the medium Baier  _et al._ (1997),

$\hat{q}=\rho\int dq_{\rm T}^{2}\frac{d\sigma}{dq_{\rm T}^{2}}q_{\rm T}^{2},$
(1)

where $\rho$ is the medium gluon number density. A great amount of abundant
experimental data Adare  _et al._ (2008, 2013, 2010); Abelev  _et al._ (2013);
Chatrchyan  _et al._ (2012); Acharya  _et al._ (2018); Khachatryan  _et al._
(2017); Acharya  _et al._ (2019) make it possible to accurately extract jet
quenching parameter. One notable work on the extraction of jet transport
coefficient $\hat{q}$ was performed by the JET Collaboration Burke  _et al._
(2014), which compared several different energy loss models and extracted the
$\hat{q}$ values from single hadron suppressions at RHIC and the LHC energies.
Besides, phenomenological investigations have been carried out to extract the
initial jet transport coefficient and the mean free path at the initial time
simultaneously Liu  _et al._ (2016), and consider bulk matter evolution Chen
_et al._ (2010, 2011) for large $p_{\rm T}$ single hadron suppression and
other jet quenching observables, such as dihadron and $\gamma$-hadron
suppressions Zhang  _et al._ (2007, 2009); Xie  _et al._ (2019, 2021). More
recently, the jet transport coefficient extractions are continuing via the
improving theoretical approaches with opacity-resummed medium induced
radiation Feal  _et al._ (2021), the quasi-particle collection for QGP based
on the linear Boltzmann transport model Liu  _et al._ (2021), as well as the
advanced analytical technique within JETSCAPE framework Cao  _et al._ (2021).
It is important and necessary for a comprehensive and thorough understanding
of the QGP properties to evaluate the systematic uncertainty of jet transport
coefficient, which is given by different theoretical methods Burke  _et al._
(2014); Feal  _et al._ (2021); Liu  _et al._ (2021), observables Chen  _et
al._ (2011); Zhang  _et al._ (2007, 2009); Xie  _et al._ (2019), even hydro
evolution information Chen  _et al._ (2010); Liu  _et al._ (2016); Xie  _et
al._ (2021), the initial parton distribution functions and the final
fragmentation functions, etc.

A recent investigation for the inclusive charged-particle spectra with a NLO
pQCD parton model using seven recent sets of parton-to-hadron fragmentation
functions (FFs) indicates that the gluon fragmentation is obviously different
in current parton FFs and the hadron production is predominantly sensitive to
the gluon-to-hadron FFs in $p+p$ collisions d’Enterria  _et al._ (2014). Since
gluon energy loss is $9/4$ times quark energy loss Deng and Wang (2010), the
use of different fragmentation functions should bring a uncertainty for the
jet transport coefficient extracted in high-energy heavy-ion collisions.
Besides, considering that single and dihadorn productions have the same jet
energy loss mechanism but slightly different production points for the initial
jets Zhang  _et al._ (2007), the study for the two observables can further
help us understand the medium properties. In this paper, we will check the
characters of 6 sets of the current NLO parton FFs Kretzer (2000); Kniehl  _et
al._ (2000); Bourhis  _et al._ (2001); Hirai  _et al._ (2007); Albino  _et
al._ (2008); de Florian  _et al._ (2007a, b), and use a NLO pQCD parton model
with jet quenching modified FFs to study high $p_{\rm T}$ hadron productions
in heavy-ion collisions Wang  _et al._ (1996); Wang and Huang (1997); Wang
(2004a); Deng and Wang (2010); Wang and Wang (2001, 2002). Confronting 6
current NLO parton FFs, we will extract the jet transport coefficient
$\hat{q}$ quantitatively via a global fit to data for both single hadron and
dihadron suppressions and check the extraction uncertainty.

The paper is organized as follows. In Sec. II, we introduce briefly the NLO
pQCD parton model and give the spectrum expressions of single hadron and
dihadon productions. In Sec.III, we give comparisons of the 6 current sets of
NLO parton FFs variously. In Sec. IV, we extract the jet transport coefficient
with 6 current sets of FFs by fitting to the experimental data for single
hadron and dihadron suppressions at RHIC and the LHC energies. A brief summary
and discussions are given in Sec. V. The above analyses use the same scales in
the parton model for different FFs. In the Appendix we add the results and
analyses of the $\hat{q}$ extracted from the hadron suppressions with the
different scales in the model for different FFs.

##  II The NLO pQCD parton model with modified fragmentation functions

###  II.1 Cross sections of single hadron and dihadron productions

According to a NLO pQCD parton model, the differential cross section of single
hadron productions in proton-proton ($p+p$) collisions can be obtained as
follows Owens (1987); Brock  _et al._ (1995),

$\displaystyle\frac{d\sigma_{pp\rightarrow h+X}}{dyd^{2}p_{\rm T}}$
$\displaystyle=\sum_{abcd}\int
dx_{a}dx_{b}f_{a/p}(x_{a},\mu^{2})f_{b/p}(x_{b},\mu^{2})$ (2)
$\displaystyle\times\frac{1}{\pi}\frac{d\sigma_{ab\rightarrow
cd}}{d\hat{t}}\frac{D_{c}^{h}(z_{c},\mu^{2})}{z_{c}}+\mathcal{O}(\alpha_{s}^{3}),$

where $f_{a/p}(x_{a},\mu^{2})$ is the parton distribution function (PDFs) for
a parton $a$ with momentum fraction $x_{a}$ from a free nucleon, and we will
take CT18 parametrizations Hou  _et al._ (2021) in following numerical
calculations. $D_{c}^{h}(z_{c},\mu^{2})$ is the parton FFs in a vacuum, for
which we will consider 6 sets of fragmentation function parametrizations in
this work. $z_{c}=p_{\rm T}/{p_{\rm T}}_{c}$ is the transverse momentum
fraction carried by the final hadrons from the parent parton $c$. Here we will
focus on hadron productions in the middle-rapidity region.
$d\sigma_{ab\rightarrow cd}/d\hat{t}$ is the differential cross section for
parton-parton hard scattering process at leading order $\alpha_{s}^{2}$. For
NLO corrections in $\mathcal{O}(\alpha_{s}^{3})$, we consider both
$2\rightarrow 3$ real tree diagram contributions and $2\rightarrow 2$ one-loop
virtual diagram contributions. In our numerical calculations, we use two
cutoffs to handle the collinear singularities and soft singularities. The
ultraviolet divergences can be solved by renormalization. For more detailed
discussions on the NLO calculations, one can find in the references Kidonakis
and Owens (2001); Harris and Owens (2002).

Similarly, the differential cross section of dihadron productions in $p+p$
collisions can be obtained as Owens (1987),

$\displaystyle\frac{d\sigma_{pp\rightarrow h_{1}+h_{2}+X}}{dPS}$
$\displaystyle=\sum_{abcd}\int\frac{dz_{c}}{z_{c}^{2}}\frac{dz_{d}}{z_{d}^{2}}f_{a/p}(x_{a},\mu^{2})f_{a/p}(x_{b},\mu^{2})$
(3) $\displaystyle\times\frac{x_{a}x_{b}}{\pi}\frac{d\sigma_{ab\rightarrow
cd}}{d\hat{t}}D_{c}^{h_{1}}(z_{c},\mu^{2})D_{d}^{h_{2}}(z_{d},\mu^{2})$
$\displaystyle\times\delta^{2}(\frac{\overrightarrow{p}_{\rm
T}^{h1}}{z_{c}}+\frac{\overrightarrow{p}_{\rm
T}^{h2}}{z_{d}})+\mathcal{O}(\alpha_{s}^{3}),$

where $dPS=dy^{h_{1}}d^{2}p_{\rm T}^{h_{1}}dy^{h_{2}}d^{2}p_{\rm T}^{h_{2}}$
in the phase space.

In high-energy nucleus-nucleus ($A+A$) collisions, the single hadron spectra
can be written as Zhang  _et al._ (2007, 2009); Chen  _et al._ (2010); Liu
_et al._ (2016),

$\displaystyle\frac{dN_{AB\rightarrow h+X}}{dyd^{2}p_{\rm T}}$
$\displaystyle=\sum_{abcd}\int
dx_{a}dx_{b}d^{2}\vec{r}t_{A}(\vec{r})t_{B}(\vec{r}+\vec{b})$ (4)
$\displaystyle\times\
f_{a/A}(x_{a},\mu^{2},\vec{r})f_{b/B}(x_{b},\mu^{2},\vec{r}+\vec{b})$
$\displaystyle\times\frac{1}{\pi}\frac{d\sigma_{ab\rightarrow
cd}}{d\hat{t}}\frac{\tilde{D}_{c}^{h}(z_{c},\mu^{2},\Delta
E_{c})}{z_{c}}+\mathcal{O}(\alpha_{s}^{3}).$

Similarly, the dihadron spectra in $A+A$ collisions can be expressed as Zhang
_et al._ (2009),

$\displaystyle\frac{dN_{AB\rightarrow h_{1}+h_{2}+X}}{dPS}$
$\displaystyle=\sum_{abcd}\int\frac{dz_{c}}{z_{c}^{2}}\frac{dz_{d}}{z_{d}^{2}}d^{2}\vec{r}t_{A}(\vec{r})t_{B}(\vec{r}+\vec{b})$
(5) $\displaystyle\times
f_{a/A}(x_{a},\mu^{2},\vec{r})f_{b/B}(x_{b},\mu^{2},\vec{r}+\vec{b})\frac{x_{a}x_{b}}{\pi}$
$\displaystyle\times\frac{d\sigma_{ab\rightarrow
cd}}{d\hat{t}}\tilde{D}^{h_{1}}_{c}(z_{c},\mu^{2},\Delta
E_{c})\tilde{D}^{h_{2}}_{d}(z_{d},\mu^{2},\Delta E_{d})$
$\displaystyle\times\delta^{2}(\frac{\overrightarrow{p}_{\rm
T}^{h1}}{z_{c}}+\frac{\overrightarrow{p}_{\rm
T}^{h2}}{z_{d}})+\mathcal{O}(\alpha_{s}^{3}).$

In Eqs. (4) and (5), $t_{A}(\vec{r})$ is the nuclear thickness function given
by the Woods-Saxon distribution Jacobs and Cooper (2000) for nucleons in a
nucleus, and is normalized as $\int d^{2}\vec{r}t_{A}(\vec{r})=A$. $\vec{b}$
is the impact parameter in $A+A$ collisions. $f_{a/A}(x_{a},\mu^{2},\vec{r})$
is the nucleus-modified parton distribution functions, which can be factorized
into the parton distribution functions inside a free nucleon
$f_{a/N}(x_{a},\mu^{2})$ and the nuclear shadowing factor
$S_{a/A}(x_{a},\mu^{2},\vec{r})$ Wang (1997); Li and Wang (2002),

$\displaystyle f_{a/A}(x_{a},\mu^{2},\vec{r})$
$\displaystyle=S_{a/A}(x_{a},\mu^{2},\vec{r})\left[\frac{Z}{A}f_{a/p}(x_{a},\mu^{2})\right.$
(6)
$\displaystyle+\left.\left(1-\frac{Z}{A}\right)f_{a/n}(x_{a},\mu^{2})\right],$

where $Z$ is the proton number of the nucleus and $A$ is the nucleus mass
number. The nuclear shadowing factor $S_{a/A}(x_{a},\mu^{2},\vec{r})$ can be
obtained using following form Emel’yanov  _et al._ (2000); Hirano and Nara
(2004),

$\displaystyle S_{a/A}(x_{a},\mu^{2},\vec{r})$
$\displaystyle=1+[S_{a/A}(x_{a},\mu^{2})-1]\frac{At_{A}(\vec{r})}{\int{d^{2}}\vec{r}[t_{A}(\vec{r})]^{2}},\
\ \ $ (7)

where $S_{a/A}(x_{a},\mu^{2})$ is given by the EPPS16 parametrizations Eskola
_et al._ (2017). Since the parton-parton scattering cross sections are
computed up to NLO, EPPS16, CT18 and FFs parametrizations are all used at NLO.

$\tilde{D}_{c}^{h}(z_{c},\mu^{2},\Delta E_{c})$ is the medium-modified
fragmentation functions and can be calculated as follows Zhang  _et al._
(2007, 2009); Wang  _et al._ (1996); Wang and Huang (1997); Wang (2004a):

$\displaystyle\tilde{D}_{c}^{h}(z_{c},\mu^{2},\Delta{E_{c}})=(1-e^{-\langle{N_{g}}\rangle})\left[\frac{z^{\prime}_{c}}{z_{c}}D_{c}^{h}(z^{\prime}_{c},\mu^{2})\right.$
$\displaystyle\phantom{X}+\left.{\langle{N_{g}}\rangle}\frac{{z_{g}}^{\prime}}{z_{c}}D_{g}^{h}({z_{g}}^{\prime},\mu^{2})\right]+e^{-\langle{N_{g}}\rangle}D_{c}^{h}({z_{c}},\mu^{2}),$
(8)

where ${z_{c}}^{\prime}=p_{\rm T}/(p_{\rm{T}c}-\Delta{E_{c}})$ is the rescaled
transverse momentum fraction of the hadron from the quenched parton. The
parton has the initial transverse momentum $p_{\rm{T}c}$ to traverse the
medium. After losing energy $\Delta{E}_{c}$, the quenched parton is fragmented
into a hadron with momentum $p_{\rm T}$. ${z_{g}}^{\prime}=p_{\rm
T}/(\Delta{E_{c}}/\langle{N_{g}}\rangle)$ is the rescaled transverse momentum
fraction of the hadron from the fragmentation of the radiated gluon with the
initial energy $\Delta{E_{c}}/\langle{N_{g}}\rangle$. $\langle{N_{g}}\rangle$
is the averaged radiation gluon number and obeys the Poisson distribution. The
factor $e^{-\langle{N_{g}}\rangle}$ is the probability for partons escaping
the medium without suffering any inelastic scattering, while the factor
$(1-e^{-\langle{N_{g}}\rangle})$ is the probability for partons encountering
at least one inelastic scattering.

$\Delta E_{c}$ is the total parton energy loss and can be calculated by the
high-twist approach Deng and Wang (2010); Wang and Wang (2001, 2002). For a
light quark $c$ with the initial energy $E$, the radiative energy loss $\Delta
E_{c}$ can be calculated as,

$\displaystyle\frac{\Delta{E}_{c}}{E}$
$\displaystyle=\frac{2C_{A}\alpha_{s}}{\pi}\int d\tau\int\frac{dl_{\rm
T}^{2}}{l_{\rm T}^{4}}\int dz$ (9)
$\displaystyle\times\left[1+(1-z)^{2}\right]\hat{q}\sin^{2}(\frac{l_{\rm
T}^{2}(\tau-\tau_{0})}{4z(1-z)E}),$

where $C_{A}=3$ and $l_{\rm T}^{2}$ is the squared transverse momentum of the
radiated gluon. Since the colour factor of gluon-gluon vertex is $9/4$ times
that of quark-gluon vertex, the energy loss of a gluon jet is $9/4$ times that
of a quark jet Deng and Wang (2010). The average number of the gluons emitted
off a hard parton is calculated as Chang  _et al._ (2014),

$\displaystyle\langle N_{g}\rangle$
$\displaystyle=\frac{2C_{A}\alpha_{s}}{\pi}\int d\tau\int\frac{dl_{\rm
T}^{2}}{l_{\rm T}^{4}}\int\frac{dz}{z}$ (10)
$\displaystyle\times\left[1+(1-z)^{2}\right]\hat{q}\sin^{2}(\frac{l_{\rm
T}^{2}(\tau-\tau_{0})}{4z(1-z)E}).$

###  II.2 The jet transport coefficient

The total parton energy loss and the number of radiated gluons are both
controlled by jet transport coefficient $\hat{q}$ (Baier  _et al._ , 1997).
According to Eq. (1) for the definition of jet transport coefficient, one can
assume the $\hat{q}$ is proportional to the local gluon density in a QGP phase
or the hadron density in a hadronic gas. Therefore the transport coefficient
in an evolving dynamical medium can be expressed as Chen  _et al._ (2010,
2011),

$\displaystyle\hat{q}(\tau,\vec{r})$
$\displaystyle=\left[\hat{q}_{0}\frac{\rho_{\scriptscriptstyle
QGP}(\tau,\vec{r}+(\tau-\tau_{0})\vec{n})}{\rho_{\scriptscriptstyle
QGP}(\tau_{0},0)}(1-f)\right.$ (11) $\displaystyle+\left.{\hat{q}_{\rm
had}(\tau,\vec{r}+(\tau-\tau_{0})\vec{n})f}\frac{}{}\right]\cdot\frac{p^{\mu}u_{\mu}}{p_{0}},$

where $\vec{n}$ is the unit length vector in the parton jet moving direction.
The first part in the above equation represents the contribution in the QGP
phase, while the second part denotes in the hadronic phase. We also consider
the effective flow dependence of jet transport coefficient. $p^{\mu}$ is the
four-momentum of the jet, and $u_{\mu}$ is the four flow velocity of the
medium. $f(\tau,\vec{r})$ is the hadronic phase fraction at a given time-space
point, which can distinguish the contribution of jet energy loss in the QGP
and the hadronic phase,

$f(\tau,\vec{r})=\begin{cases}0&T>T_{c}\\\ 0\sim 1&T=T_{c}\\\
1&T<T_{c}.\end{cases}$ (12)

Here $T$ is the local temperature of the medium. In our studies, the
hydrodynamical time-space evolution information of medium temperature $T$ and
flow velocity $u$ are obtained by the (3+1)-dimensional ideal hydrodynamic
model Hirano (2002); Hirano and Tsuda (2002), in which there is a mixed phase
platform of the first order phase transition between the QGP and hadron phases
at the critical temperature $T_{c}=170$ MeV.

For the QGP phase in Eq. 11, $\hat{q}_{0}$ denotes the jet transport
coefficient at the center of the bulk medium in the initial time $\tau_{0}$.
$\rho_{\scriptscriptstyle QGP}$ is the parton density at a given time and
space, which is propotional to the temperature cubed Hirano (2002); Hirano and
Tsuda (2002); Chen  _et al._ (2010, 2011). Thereby in numerical caculations,
we assume that $\hat{q}$ for the QGP phase has the following form as Liu  _et
al._ (2016); Xie  _et al._ (2019),

$\displaystyle\hat{q}=\hat{q}_{0}\frac{T^{3}}{T_{0}^{3}}\frac{p^{\mu}u_{\mu}}{p_{0}}(1-f),$
(13)

where $T_{0}$ is a reference temperature taken as the highest temperature in
the center of the medium at the initial time $\tau_{0}$.

Table 1:  Characteristics of the six sets of current fragmentation function parameterizations with hadron species, fitted data, minimum value of $z$, and the scale range of $\mu^{2}$. FFs set | Species | Fitted data | $z_{min}$ |  $\mu^{2}$(GeV2)
---|---|---|---|---
KRE Kretzer (2000) | $\pi^{+},\pi^{-},K^{+},K^{-},h^{+},h^{-}$ | $e^{+}e^{-}$ | 0.01 | $0.8-10^{6}$
KKP Kniehl  _et al._ (2000) |  $\pi^{+}+\pi^{-}$, $K^{+}+K^{-}$, $p+\overline{p}$, $h^{+}+h^{-}$ | $e^{+}e^{-}$ | 0.1 | $2.0-4.0\cdot 10^{4}$
BFGW Bourhis  _et al._ (2001) | $h^{+}+h^{-}$ | $e^{+}e^{-}$ | 0.001 | $2.0-1.2\cdot 10^{4}$
HKNS Hirai  _et al._ (2007) | $\pi^{+},\pi^{-},K^{+},K^{-},p,\overline{p},h^{+},h^{-}$ | $e^{+}e^{-}$ | 0.01 | $1.0-10^{8}$
AKK08 Albino  _et al._ (2008) | $\pi^{+},\pi^{-},K^{+},K^{-},p,\overline{p},h^{+},h^{-}$ |  $e^{+}e^{-}$, $pp$ | 0.05 | $1.0-10^{6}$
DSS de Florian  _et al._ (2007a, b) | $\pi^{+},\pi^{-},K^{+},K^{-},p,\overline{p},h^{+},h^{-}$ |  $e^{+}e^{-}$, $pp$, $ep$ | 0.05 | $1.0-10^{5}$

For the hadronic phase in Eq. 11, $\hat{q}_{\rm had}$ is the jet transport
coefficient in the hadronic phase and can be written as Chen  _et al._ (2010),

$\displaystyle\hat{q}_{\rm
had}=\frac{\hat{q}_{N}}{\rho_{N}}[\frac{2}{3}\sum\limits_{M}\rho_{M}(T)+\sum\limits_{B}\rho_{B}(T)],$
(14)

where $\hat{q}_{N}$ is the extracted jet transport coefficient at the center
of a large nucleus and given by $\hat{q}_{N}\approx 0.02$ GeV${}^{2}/$fm, and
$\rho_{N}\approx 0.17$ fm-3 is the nucleon density at the center of the large
nucleus Chen  _et al._ (2010). $\rho_{M}$ and $\rho_{B}$ are the meson and
baryon density in the hadronic resonance gas at a given temperature,
respectively. The factor 2/3 denotes the ratio of constituent quark numbers in
mesons and baryons. The hadron density at a temperature $T$ and zero chemical
potential is expressed as Chen  _et al._ (2010),

$\displaystyle\sum\limits_{h}\rho_{h}(T)=\frac{T^{3}}{2\pi^{2}}\sum\limits_{h}(\frac{m_{h}}{T})^{2}\sum\limits_{n=1}^{\infty}\frac{\eta_{h}^{n+1}}{n}K_{2}(n\frac{m_{h}}{T}),$
(15)

where $\eta_{h}=\pm 1$ for meson (M)/baryon (B). In this paper, we will
consider hadron resonances, and the mass is below 1 GeV, including all 17
kinds of mesons: $\pi^{+}$, $\pi^{-}$, $\pi^{0}$, $K^{+}$, $K^{-}$, $K^{0}$,
$\overline{K^{0}}$, $\eta$, $\eta^{{}^{\prime}}$, $\rho^{+}$, $\rho^{-}$,
$\rho^{0}$, $K^{*+}$, $K^{*-}$, $K^{*0}$, $\overline{K^{*0}}$, $\omega$; and 2
kinds of baryons: $p$, $n$.

###  II.3 Nuclear modification factors

Confronting the 6 current sets of parton fragmentation functions, we will
extract the parameter $\hat{q}_{0}$ using the $\chi^{2}$ fitting method of
comparing the NLO pQCD numerical results of single hadron and dihadron
suppressions with the experimental data. To demonstrate the suppression of the
single hadron spectrum in $A+A$ collisions relative to that in $p+p$
collisions, one can define the nuclear modification factor $R_{AA}(p_{\rm T})$
as Wang (2004a),

$\displaystyle R_{AA}(p_{\rm
T})=\frac{1}{T_{AA}(\vec{b})}\frac{dN_{AA\rightarrow h+X}/dyd^{2}p_{\rm
T}}{d{\sigma}_{pp\rightarrow h+X}/dyd^{2}p_{\rm T}},$ (16)

where $T_{AA}(\vec{b})=\int d^{2}\vec{r}t_{A}(\vec{r})t_{A}(\vec{r}+\vec{b})$
is the overlap function of two colliding nuclei for a given impact parameter.

In the following numerical studies for dihadrons, for a given triggered
hadron, we will focus on the away-side associated hadrons. The nuclear
modification factor $I_{AA}$ for dihadron productions can be defined as a
function of $z_{\rm T}=p_{\rm T}^{\rm assoc}/p_{\rm T}^{\rm trig}$ or a
function of $p_{\rm T}^{\rm assoc}$,

$I_{AA}(z_{\rm T})=\frac{D_{AA}(z_{\rm T})}{D_{pp}(z_{\rm T})},$ (17)

or

$I_{AA}(p_{\rm T}^{\rm assoc})=\frac{D_{AA}(p_{\rm T}^{\rm
assoc})}{D_{pp}(p_{\rm T}^{\rm assoc})},$ (18)

where $D_{AA}(z_{\rm T})=p_{\rm T}^{\rm trig}D_{AA}(p_{\rm T}^{\rm assoc})$ is
called hadron-triggered fragmentation function, and can be calculated by Wang
(2004b),

$D_{AA}(z_{\rm T})=p_{\rm T}^{h_{1}}\frac{dN_{AA\rightarrow
h_{1}+h_{2}+X}/dy^{h_{1}}dp_{\rm T}^{h_{1}}dy^{h_{2}}dp_{\rm
T}^{h_{2}}}{dN_{AA\rightarrow h_{1}+X}/dy^{h_{1}}dp_{\rm T}^{h_{1}}}.$ (19)

We calculate the $R_{AA}$ for single hadrons and $I_{AA}$ for dihadrons with
the six sets of FFs at RHIC and the LHC energies, and compare the numerical
results to the experimental data by utilizing $\chi^{2}/d.o.f$ fitting method.
The $\chi^{2}/d.o.f$ is defined as follows,

$\chi^{2}/d.o.f=\sum_{i=1}^{N}\left[\frac{(V_{th}-V_{exp})^{2}}{\sigma_{sys}^{2}+\sigma_{stat}^{2}}\right]_{i}/N$
(20)

where $V_{th}$ represents the theoretical value, $V_{exp}$ denotes the
experimental data, $\sigma_{sys}$ and $\sigma_{stat}$ are the systematic and
statistical errors for the experimental data, and $N$ is the number of data
points which are used.

##  III Comparisons of parton fragmentation functions

In this work, we employ the six sets of commonly-used parameterizations of NLO
parton-to-hadron fragmentation functions: Kretzer (KRE) Kretzer (2000), KKP
Kniehl  _et al._ (2000), BFGW Bourhis  _et al._ (2001), HKNS Hirai  _et al._
(2007), AKK08 Albino  _et al._ (2008), and DSS de Florian  _et al._ (2007a,
b), as listed in Table 1. From the species column, we can see that the charged
hadrons ($h^{\pm}$) are constructed as a sum of the individual FFs for pions
($\pi^{\pm}$), kaons ($K^{\pm}$), and (anti)protons($p,\overline{p}$) in most
of the parton FFs sets except the BFGW FFs. $\pi^{0}$ hadron is given by the
average of pions ($\pi^{+}$ and $\pi^{-}$). In this section we will briefly
show the differences of parton-to-$\pi^{0}$ or parton-to-$h^{\pm}$ between
these six sets of fragmentation function parameterizations.

Fig. 1: The invariant cross sections of single $\pi^{0}$ for the same scale
$\mu=1.2p_{\rm T}$ with different FFs in $p+p$ collisions at $\sqrt{s_{\rm
NN}}=200$ GeV (upper panel), and ratios of data over the theoretical results
(lower panel). The data are from Ref. Adare  _et al._ (2008).  Fig. 2:
Comparisons of parton-to-$\pi^{0}$ fragmentation functions between 5 sets of
FFs, KREKretzer (2000), KKP Kniehl  _et al._ (2000), HKNS Hirai  _et al._
(2007), AKK08 Albino  _et al._ (2008), and DSS de Florian  _et al._ (2007a,
b). Left panels are for quark FFs, and right panels for gluon FFs. Upper
panels are for hadrons with $p_{\rm T}=10$ GeV, and lower panels with $p_{\rm
T}=50$ GeV. The scale is set as $\mu=1.2p_{\rm T}$. Fig. 3: The contribution
fractions of quark (solid lines) and gluon (dashed lines) fragmentations to
the inclusive $\pi^{0}$ cross sections for different FFs in $p+p$ collisions
at $\sqrt{s_{\rm NN}}=200$ GeV.

###  III.1 Characteristics of parton FFs for $\pi^{0}$ hadrons

In the NLO pQCD parton model for hadron productions, there are three
independent scales: the factorization scale $\mu_{\rm fac}$, the
renormalization scale $\mu_{\rm ren}$ and the fragmentation scale $\mu_{\rm
fra}$. In following numerical calculations we choose the equal scales
$\mu_{\rm fac}=\mu_{\rm ren}=\mu_{\rm fra}$, and let them proportional to the
physics observables, such as hadron transverse momentum $p_{\rm T}$ or
dihadron invariant mass $M$.

As shown in the upper panel of Fig. 1, we firstly give the NLO pQCD results of
differential cross sections for $\pi^{0}$ productions at large transverse
momentum $p_{\rm T}$ with the 5 sets of FFs in $p+p$ collisions at
$\sqrt{s_{\rm NN}}=200$ GeV, and compare them with the experimental data Adare
_et al._ (2008). Here we choose the scale as $\mu=1.2p_{\rm T}$. The numerical
results can describe the experimental data. There are certain differences in
the specific results from the different fragmentation function
parameterizations, as shown in the lower panel of Fig. 1 for the ratios of the
experimental data to theoretical calculations.

Note that such differences mainly result from the discrepancies of gluon-
to-$\pi^{0}$ FFs d’Enterria  _et al._ (2014). We give in Fig. 2 the
comparisons of the quark (including u, d, and s) and gluon FFs for $\pi^{0}$
hadrons at $p_{\rm T}=10$ GeV and 50 GeV, respectively. The FFs of quark-
to-$\pi^{0}$ give similar contributions, while the gluon-to-$\pi^{0}$ FFs
provide significantly different contributions during the whole range of $z$.
For examples, since KRE and HKNS parameterizations have less gluon
contributions relative to other sets, the model with KRE or HKNS FFs
underestimates $\pi^{0}$ productions in $p+p$ collisions at 200 GeV as shown
in Fig. 1. Of course, in the NLO pQCD parton model using KRE and HKNS
fragmentation function parameterizations, one can add a $K$ factor ($K>1.0$)
to account for higher-order contributions or choose a smaller scale $\mu$ to
fit data better for the $p_{\rm T}$ spectra in $p+p$ collisions.

To make clear the relative contributions of quarks and gluons to final state
hadrons in these sets of FFs, we show in Fig. 3 the contribution fractions of
quarks (solid lines) or gluons (dashed lines) to the inclusive charged-hadron
cross sections in $p+p$ collisions. One can see that with the same scale
$\mu=1.2p_{\rm T}$, the relative contributions of quarks and gluons to
$\pi^{0}$ productions are both distinctly different for all the sets of FFs
due to the apparent differences in gluon fragmentation functions.

Fig. 4: The invariant cross sections of single charged hadrons for the same
scale $\mu=1.5p_{\rm T}$ with different FFs in $p+p$ collisions at
$\sqrt{s_{\rm NN}}=2.76$ TeV (upper panel), and ratios of data over the
theoretical results (lower panel). The data are from Refs. Chatrchyan  _et
al._ (2012); Abelev  _et al._ (2013). Fig. 5: Similar to Fig. 2 but 6 sets of
FFs for charged hadrons. The scale is set as $\mu=1.5p_{\rm T}$. Fig. 6: The
contribution fractions of quark (solid lines) and gluon (dashed lines)
fragmentations to the inclusive charged-hadron cross sections for different
FFs in $p+p$ collisions at $\sqrt{s_{\rm NN}}=2.76$ TeV.

###  III.2 Characteristics of parton FFs for charged hadrons

Similar to the study at RHIC energy, we also calculate the cross sections of
charged hadrons in $p+p$ collisions at $\sqrt{s_{\rm NN}}=2.76$ TeV, and give
ratios of the experimental data to theoretical results, as shown in Fig. 4. We
choose $\mu=1.5p_{\rm T}$ as the default scale for charged hadron productions.
The numerical results for charged hadrons with different FFs can also describe
the experimental data, but the model with the most of the FFs overestimates
slightly the charged hadron productions with $\mu=1.5p_{\rm T}$.

Fig. 5 shows the comparisons for the quark and gluon FFs of charged hadrons at
$p_{\rm T}=10$ GeV (upper panels) and 50 GeV (lower panels) from the 6 sets of
FFs, respectively. Likewise to pion productions, the FFs of quark-to-$h^{\pm}$
give similar contributions, but the gluon-to-$h^{\pm}$ FFs show significantly
different contributions during the whole range of $z$. Compared to KRE gluon
fragmentation, the too-hard gluon contributions of the other sets of FFs lead
to the overestimation in Fig. 4 for the charged hadron productions in $p+p$
collisions at $\sqrt{s_{\rm NN}}=2.76$ TeV.

For a more intuitive view, the relative contributions from the quark and gluon
fragmentations are given in Fig. 6 for the single charged hadrons in $p+p$
collisions at $\sqrt{s_{\rm NN}}=2.76$ TeV. Similarly to Fig. 3, there exists
different fractions of gluon (quark) contributions to hadrons among the 6 sets
of FFs.

##  IV Jet transport coefficient extractions

Because gluon energy loss is 9/4 times of quark energy loss, the different
fractions of gluon (quark) contributions to hadrons in the model with
different sets of FFs will give different strength of jet quenching for high
$p_{\rm T}$ hadron productions in high-energy nucleus-nucleus collisions. In
this section, confronting the difference in the six sets of FFs, let’s check
the jet quenching parameter $\hat{q}_{0}$ extracted from both the single
hadron and dihadron suppressions by fitting data at RHIC and the LHC.

###  IV.1 Extracting $\hat{q}_{0}$ in central $Au+Au$ collisions at
$\sqrt{s_{\rm NN}}=200$ GeV

We choose $b=2.0$ (3.2) fm in Eq. (16) for $R_{AA}(p_{\rm T})$ and Eq. (17)
for $I_{AA}(z_{\rm T})$ in 0-5% (0-10%) Au+Au collisions. Shown in Fig. 7 are
the nuclear modification factors $R_{AA}(p_{\rm T})$ for single hadrons (left
panels) in $0-5\%$ $Au+Au$ collisions, and $I_{AA}(z_{\rm T})$ for dihadrons
(right panels) in $0-10\%$ $Au+Au$ collisions at $\sqrt{s_{\rm NN}}=200$ GeV.
KREKretzer (2000), KKP Kniehl  _et al._ (2000), HKNS Hirai  _et al._ (2007),
AKK08 Albino  _et al._ (2008), and DSS de Florian  _et al._ (2007a, b) FFs are
used with the same scale $\mu=1.2p_{\rm T}$ for single hadrons and $\mu=1.2M$
for dihadrons, respectively. Theoretical results fit data well with several
appropriate values of jet transport coefficient $\hat{q}_{0}$. The solid
curves are given by the best fitting from the next $\chi^{2}$-fitting
calculations.

Shown in Fig. 8 are the $\chi^{2}/d.o.f$ fits to nuclear modification factors
in central $Au+Au$ collisions at $\sqrt{s_{\rm NN}}=200$ GeV, (a) the fits to
only single hadron $R_{AA}(p_{\rm T})$, (b) the fits to only dihadron
$I_{AA}(z_{\rm T})$, and (c) the global fits to $R_{AA}(p_{\rm T})$ \+
$I_{AA}(z_{\rm T})$. The 5 sets of fragmentation function parameterizations
are used respectively. $p_{\rm T}>5$ GeV for single hadrons $R_{AA}(p_{\rm
T})$ and $z_{\rm T}>0.25$ for $\pi^{0}$-triggered away-side charged hadrons
$I_{AA}(z_{\rm T})$ are chosen to fit data.

The same mechanism of jet quenching leads to both suppressions of large
$p_{\rm T}$ single hadrons and $\pi^{0}$-triggered away-side charged hadrons
in heavy-ion collisions, so the best-fitting values of $\hat{q}_{0}$ for
single hadrons in panel (a) of Fig. 8 are similar to those for dihadrons in
panel (b) for most sets of FFs used in the theoretical model. From a global
$\chi^{2}$ fits to both the single hadron and dihadron suppressions in panel
(c) of Fig. 8, the best-fitting values of $\hat{q}_{0}$ for each set of FFs
are obtained as: $\hat{q}_{0}=1.7$ GeV2/fm with KRE FFs, $\hat{q}_{0}=1.4$
GeV2/fm with KKP FFs, $\hat{q}_{0}=1.7$ GeV2/fm with HKNS FFs,
$\hat{q}_{0}=1.3$ GeV2/fm with AKK08 FFs, and $\hat{q}_{0}=1.2$ GeV2/fm with
DSS FFs.

Fig. 7: The nuclear modification factor $R_{AA}$ for single hadron productions
as a function of $p_{\rm T}$ (left panels) in $0-5\%$ $Au+Au$ collisions, and
$I_{AA}$ for dihadron productions as a function of $z_{\rm T}$ (right panels)
in $0-10\%$ $Au+Au$ collisions at $\sqrt{s_{\rm NN}}=200$ GeV. The
experimental data are from Refs. Adare  _et al._ (2008, 2013); Adamczyk  _et
al._ (2016).

Fig. 8: The $\chi^{2}/d.o.f$ results from fitting to nuclear modification
factors in central $Au+Au$ collisions at $\sqrt{s_{\rm NN}}=200$ GeV, (a) the
fits to only single hadron $R_{AA}(p_{\rm T})$, (b) the fits to only dihadron
$I_{AA}(z_{\rm T})$, and (c) the global fits to $R_{AA}(p_{\rm T})$ \+
$I_{AA}(z_{\rm T})$. 5 sets of fragmentation function parameterizations are
used, respectively.

The large differences in the values $\hat{q}_{0}=1.2-1.7$ GeV2/fm are obtained
by using different sets of FFs. Considering parton-to-hadron contributions in
Fig. 2 and 3, one can see that such differences for $\hat{q}_{0}$ extraction
are mainly from and sensitive to the differences of gluon-to-hadron in FFs.
For example, as shown in Fig. 3 with the different sets of FFs, one can get
the least fraction of gluon-to-hadron contribution with KRE parametrization
while the largest fraction with DSS parametrization. Since gluon energy loss
is 9/4 times of quark energy loss, model calculations with KRE FFs give
$\hat{q}_{0}=1.7$ GeV2/fm larger than $\hat{q}_{0}=1.2$ GeV2/fm with DSS FFs
for the almost same total energy loss for $R_{AA}(p_{\rm T})$ to fit data well
in central $A+A$ collisions.

Fig. 9: The nuclear modification factor $R_{AA}$ for single hadron as a
function of $p_{\rm T}$ compared to the experimental data from ALICE Abelev
_et al._ (2013) and CMS Chatrchyan  _et al._ (2012) in $0-5\%$ $Pb+Pb$
collisions at $\sqrt{s_{\rm NN}}=2.76$ TeV.

###  IV.2 Extracting $\hat{q}_{0}$ in central $Pb+Pb$ collisions at
$\sqrt{s_{\rm NN}}=2.76$ TeV

Similar to at RHIC energy, we perfrom the same analyses for single hadron
$R_{AA}(p_{\rm T})$ and dihadron $I_{AA}(p_{\rm T}^{\rm assoc})$
(corresponding to Eq. (18)) in central $Pb+Pb$ collisions at $\sqrt{s_{\rm
NN}}=2.76$ TeV. BFGW fragmentation function parameterizations Bourhis  _et
al._ (2001) is also applied for charged hadrons besides the other 5 sets of
FFs. Shown in Fig. 9 are the nuclear modification factor $R_{AA}(p_{\rm T})$
for single hadrons in 0-5% $Pb+Pb$ collisions with the same scale
$\mu=1.5p_{\rm T}$ for the 6 sets of FFs. In the meanwhile, the nuclear
modification factor $I_{AA}(p_{\rm T}^{\rm assoc})$ for dihadrons with the
same scale $\mu=1.5M$ are presented in Fig. 10 and Fig. 11, which contain the
results for different $p_{\rm T}^{\rm trig}$ ranges of the triggered hadrons,
respectively. As shown in Fig. 10, the trigger transverse momentum is chosen
with a relatively small range ($8<p_{\rm T}^{\rm trig}<16$ GeV), and the
theoretical results for $I_{AA}(p_{\rm T}^{\rm assoc})$ are compared to the
ALICE experimental data Aamodt  _et al._ (2012). The solid curve denotes the
numerical result with the best fitting $\hat{q}_{0}$ obtained by the
$\chi^{2}$ fitting to the experimental data. In Fig. 11, each panel with each
set of FFs contains four sub-panels in which theoretical results are compared
with the CMS data Conway (2013) for $I_{AA}(p_{\rm T}^{\rm assoc})$ with four
different large $p_{\rm T}^{\rm trig}$ ranges. The theoretical results fit
data well with several appropriate values of $\hat{q}_{0}$ at the LHC energy.

Fig. 10: The nuclear modification factor $I_{AA}$ as a function of $p_{\rm
T}^{\rm assoc}$ for dihadrons compared to the experimental data from ALICE
Aamodt  _et al._ (2012) in $0-10\%$ $Pb+Pb$ collisions at $\sqrt{s_{\rm
NN}}=2.76$ TeV.

Shown in Fig. 12 are the $\chi^{2}/d.o.f$ fits to nuclear modification factors
in central $Pb+Pb$ collisions at $\sqrt{s_{\rm NN}}=2.76$ TeV, (a) the fits to
only single hadron $R_{AA}(p_{\rm T})$, (b) the fits to only dihadron
$I_{AA}(p_{\rm T}^{\rm assoc})$, and (c) the global fits to $R_{AA}(p_{\rm
T})$ \+ $I_{AA}(p_{\rm T}^{\rm assoc})$. The fitted data include ALICE Abelev
_et al._ (2013); Aamodt  _et al._ (2012) and CMS Chatrchyan  _et al._ (2012);
Conway (2013) data. The 6 sets of fragmentation function parameterizations are
used, respectively. $p_{\rm T}>5$ GeV for single hadrons $R_{AA}(p_{\rm T})$
and $p_{\rm T}^{\rm assoc}>3.5$ GeV for $h^{\pm}$-triggered away-side charged
hadrons $I_{AA}(p_{\rm T}^{\rm assoc})$ are chosen to fit data. From the
global-fit results of Fig. 12 (c), we can read the best-fitting values of
$\hat{q}_{0}$ for each set of FFs as: $\hat{q}_{0}=2.2$ GeV${}^{2}/$fm with
KRE FFs, $\hat{q}_{0}=2.1$ GeV${}^{2}/$fm with KKP FFs, $\hat{q}_{0}=1.8$
GeV${}^{2}/$fm with BFGW FFs, $\hat{q}_{0}=1.7$ GeV${}^{2}/$fm with HKNS FFs,
$\hat{q}_{0}=1.9$ GeV${}^{2}/$fm with AKK08 FFs, $\hat{q}_{0}=1.8$
GeV${}^{2}/$fm with DSS FFs. The difference in the extracted values of
$\hat{q}_{0}=1.7-2.2$ GeV${}^{2}/$fm from different sets of FFs is caused by
the different contributions of gluon-to-hadron FFs, as shown in Fig. 5 and 6,
similarly to the case of RHIC energy.

Fig. 11: The nuclear modification factor $I_{AA}$ as a funtion of $p_{\rm
T}^{\rm assoc}$ for dihadrons with different $p_{\rm T}^{\rm trig}$ ranges
compared to the experimental data from CMS Conway (2013) in $0-10\%$ $Pb+Pb$
collisions at $\sqrt{s_{\rm NN}}=2.76$TeV.

Fig. 12: The $\chi^{2}/d.o.f$ results from fitting to nuclear modification
factors in central $Pb+Pb$ collisions at $\sqrt{s_{\rm NN}}=2.76$ TeV, (a) the
fits to only single hadron $R_{AA}(p_{\rm T})$ Abelev  _et al._ (2013);
Chatrchyan  _et al._ (2012), (b) the fits to only dihadron $I_{AA}(p_{\rm
T}^{\rm assoc})$ Aamodt  _et al._ (2012); Conway (2013) data, and (c) the
global fits to $R_{AA}(p_{\rm T})$ \+ $I_{AA}(p_{\rm T}^{\rm assoc})$. The 6
sets of fragmentation function parameterizations are used, respectively.

Note that in panels (a) and (b) of Fig. 12, with the same set of FFs, the
$\hat{q}_{0}$ value from dihadron suppressions is larger than that from single
hadron suppressions, which is mainly due to the different production positions
of the initial parton jets for single hadrons and dihadrons. On one hand, the
single hadrons are mainly contributed from the surface emissions of single
jets perpendicular to the surface of the hot medium, which means that a large
number of single jets produced at the center of the medium are melted by the
medium Zhang  _et al._ (2007). In contrast, the dihadrons are contributed from
a combination of tangential dijets generated at the medium surface and
punching-through dijets created at the matter center with limited energy loss
Zhang  _et al._ (2007). Thus in an $A+A$ event the average jet energy loss for
dihadrons is smaller than single hadrons. If one wants to obtain similar
hadron suppressions shown as the experimental measurements, larger average
energy loss is needed for dihadrons than single hadrons. On the other hand,
the fraction of the punching-through jets becomes larger with the increasing
$p_{\rm T}^{\rm trig}$ Zhang  _et al._ (2008). That means as the $p_{\rm
T}^{\rm trig}$ increases, the suppression of dihadrons should become less,
which is not manifested obviously by the experimental data with relatively
considerable uncertainty Conway (2013). With the larger $p_{\rm T}^{\rm
trig}$, one need to adjust $\hat{q}_{0}$ to be larger for $I_{AA}$ to fit the
experimental measurements. Therefore, in the global fit for $I_{AA}$ with all
$p_{\rm T}^{\rm trig}$ ranges, the values of $\hat{q}_{0}$ extracted from
dihadrons will be enlarged. Such significant difference of $\hat{q}_{0}$
between single hadrons and dihadrons can be attenuated by the abundant and
accurate single hadron data, as shown in panel (c) of Fig. 12.

###  IV.3 The uncertainty for $\hat{q}_{0}$ extraction due to FFs

With the same scale in the NLO pQCD parton model including different sets of
FFs for large $p_{\rm T}$ hadron productions in central $A+A$ collisions at
RHIC and the LHC, we calculate both the nuclear modification factors for
single hadron and dihadron suppressions due to jet quenching, and make a
global fit to data to extract the jet energy loss parameter $\hat{q}_{0}$, as
shown in Fig. 8 (c) and Fig. 12 (c). To illustrate clearly the uncertainties
from different FFs, we summarize all the extracted values for $\hat{q}_{0}$ to
show in the left panel of Fig. 13. The blue squares are for the scaled jet
transport coefficient of $\hat{q}_{0}/T_{0}^{3}$ in central $Au+Au$ collisions
at $\sqrt{s_{\rm NN}}=200$ GeV, and the center temperature of the medium in
the initial time is chosen as $T_{0}=370$ MeV. The pink circles denote the
results for central $Pb+Pb$ collisions at $\sqrt{s_{\rm NN}}=2.76$ TeV, and
$T_{0}=480$ MeV. The dotted-curve boxes indicate the uncertainties for
$\hat{q}_{0}/T_{0}^{3}$ from the different sets of FFs, and we can read that
with the same scales in the model, $\hat{q}/T^{3}=4.74-6.72$ at $T=370$ MeV
and $\hat{q}/T^{3}=3.07-3.98$ at $T=480$ MeV.

Fig. 13: Comparisons between the best-fitting values for the scaled
$\hat{q}_{0}/T_{0}^{3}$ extracted with different sets of FFs in the
theoretical model at RHIC and the LHC, left panel for the case of the same
scales in the model with different sets of FFs, and right panel for the
different scales. These results are also compared to those from JET
collaboration Burke  _et al._ (2014) and a previous theoretical study Xie  _et
al._ (2021) including the VISH (2+1)-dimension hydrodynamics Song and Heinz
(2008); Qiu  _et al._ (2012).

In contrast to the case of the same scale $\mu$ in the theoretical model for
different sets of FFs, we also accomplish the $\hat{q}_{0}$ extraction with
different scales in the model for different sets of FFs, as shown in the right
panel of Fig. 13. The suitable scale is included in the model for each set of
FFs to make numerical hadron spectra fit data well in $p+p$ collisions. The
best fitting values of $\hat{q}_{0}/T_{0}^{3}$ with different scales are
summarized as, $\hat{q}/T^{3}=4.74-5.53$ at $T=370$ MeV in central $Au+Au$
collisions at $\sqrt{s_{\rm NN}}=200$ GeV, and $\hat{q}/T^{3}=3.98-4.52$ at
$T=480$ MeV in central $Pb+Pb$ collisions at $\sqrt{s_{\rm NN}}=2.76$ TeV,
respectively.

The $\hat{q}_{0}/T_{0}^{3}$ uncertainty from different sets of FFs is reduced
significantly when the hadron-spectrum baseline of $p+p$ collisions is well
adjusted with different scales. For a given set of FFs used in the NLO pQCD
parton model for large $p_{\rm T}$ hadron spectra in $p+p$ collisions, the
ratio of gluon over quark contributions to hadrons decreases with the
increasing of the chosen scale, as shown in Fig. 3 and 6, or Fig. 20 and 21 of
Appendix B. The extracted $\hat{q}_{0}$ in central $A+A$ collisions increases
with the increasing of the chosen scale due to gluon energy loss being $9/4$
times of quark energy loss, so the $\hat{q}_{0}/T_{0}^{3}$ uncertainty from
different sets of FFs is balanced to some extent relative to the case of the
same scale for different sets of FFs. The more detailed discussions are
included in the following Appendix.

Finally, we also compare our results to those from JET collaboration Burke
_et al._ (2014) (black triangles) and a previous theoretical study Xie  _et
al._ (2021) (green stars) including the VISH (2+1)-dimension hydrodynamics
Song and Heinz (2008); Qiu  _et al._ (2012). All the results are consistent,
although different theoretical models or hydrodynamics cause some systematic
uncertainty for the accurate extraction of jet transport coefficient.

##  V SUMMARY

Based on the next-to-leading-order perturbative QCD parton model incorporating
modified fragmentation functions, we have studied single hadron and dihadron
productions in high-energy heavy-ion collisions at both RHIC and the LHC
energies. The fragmentation functions are modified due to jet quenching in
central $A+A$ collisions, the strength of which is characterized by the jet
transport coefficient $\hat{q}$. The six current sets of NLO fragmentation
function parameterizations are used in actual calculations for the nuclear
modification factors $R_{AA}$ for single hadrons and $I_{AA}$ for dihadrons.
We perform a global $\chi^{2}$ fitting to both $R_{AA}$ and $I_{AA}$ data to
extract the jet transport coefficient in the initial time at the center of the
created QGP medium, and obtain $\hat{q}/T^{3}=4.74-6.72$ at $T=370$ MeV in
central $Au+Au$ collisions at $\sqrt{s_{\rm NN}}=200$ GeV, and
$\hat{q}/T^{3}=3.07-3.98$ at $T=480$ MeV in central $Pb+Pb$ collisions at
$\sqrt{s_{\rm NN}}=2.76$ TeV.

The numerical results show that the significant uncertainties for
$\hat{q}/T^{3}$ extraction are mainly brought by the different contributions
of gluon-to-hadron in the different sets of fragmentation function
parameterizations due to gluon energy loss being $9/4$ times of quark energy
loss. The uncertainties are reduced, if the suitable scale $\mu$ is chosen in
the NLO pQCD parton model with each set of current fragmentation function
parameterizations to fit data well for large $p_{\rm T}$ hadron spectra in
$p+p$ collisions. However, the accurate parameterizations from a forthcoming
global fit of parton-to-hadron fragmentation functions d’Enterria  _et al._
(2014) will help to constrain the uncertainties for jet quenching parameter
extractions.

## ACKNOWLEDGMENTS

This work is supported in by National Natural Science Foundation of China
under Grants No. 11935007, Guangdong Major Project of Basic and Applied Basic
Research No. 2020B0301030008, and Science and Technology Program of Guangzhou
No. 2019050001.

## APPENDIX

###  V.1 The $\hat{q}_{0}$ extractions with different scales for different
sets of FFs

To check the nuclear or medium effects in $A+A$ collisions,
phenomenologically, one need a suitable baseline of $p+p$ collisions by
adjusting the scale $\mu$ to fit data well in $p+p$ collisions. Shown in the
upper panel of Fig. 14 are the numerical results for $\pi^{0}$ hadron spectra
at large transverse momentum $p_{\rm T}$ with the 5 sets of FFs in $p+p$
collisions at $\sqrt{s_{\rm NN}}=200$ GeV. Here the suitable scale $\mu$ is
chosen in the theoretical model with each set of FFs for fitting to the
experimental data Adare  _et al._ (2008), respectively. The lower panel of
Fig. 14 are for the ratios of the experimental data over theoretical
calculations. One can see that with the appropriate scales $\mu$ in the model
with different sets of FFs, the numerical results fit the data better relative
to those in Fig. 1.

Fig. 14: $\pi^{0}$ spectra given by the NLO pQCD parton model with the
befitted scale $\mu$ for each set of FFs in $p+p$ collisions at $\sqrt{s_{\rm
NN}}=200$ GeV (upper panel), and spectrum ratios of data over the theoretical
results (lower panel). The data are from Ref. Adare  _et al._ (2008).  Fig.
15: Charged hadron spectra given by the NLO pQCD parton model with the
befitted scale $\mu$ in $p+p$ collisions at $\sqrt{s_{\rm NN}}=2.76$ TeV
(upper panel), and spectrum ratios of data over the theoretical results (lower
panel). The data are from Refs. Chatrchyan  _et al._ (2012); Abelev  _et al._
(2013).

Similarly, in $p+p$ collisions at $\sqrt{s_{\rm NN}}=2.76$ TeV the charged
hadron spectra are also calculated with different scales in the NLO pQCD
parton model for the 6 sets of FFs, as shown in the upper panel of Fig. 15.
The ratios of the experimental data over theoretical calculations are shown in
the lower panel. With the appropriate scales, the theoretical results fit data
very well.

Using the above hadron spectra in $p+p$ collisions as baselines, we extract
the jet quenching parameter from single hadron and dihadron suppressions with
suitable scale in each set of FFs. Shown in Fig. 16 are the $\chi^{2}/d.o.f$
fits to nuclear modification factors in central $Au+Au$ collisions at
$\sqrt{s_{\rm NN}}=200$ GeV: panel (a) the fits to only single hadron
$R_{AA}(p_{\rm T})$, panel (b) the fits to only dihadron $I_{AA}(z_{\rm T})$,
and panel (c) the global fits to $R_{AA}(p_{\rm T})$ \+ $I_{AA}(z_{\rm T})$.
From the panel (c), we can read the best-fitting values of jet transport
coefficient: $\hat{q}_{0}=1.4$ GeV2/fm with KRE FFs at $\mu=0.6p_{\rm T}$,
$\hat{q}_{0}=1.4$ GeV2/fm with KKP FFs at $\mu=1.2p_{\rm T}$,
$\hat{q}_{0}=1.4$ GeV2/fm with HKNS FFs at $\mu=0.6p_{\rm T}$,
$\hat{q}_{0}=1.2$ GeV2/fm with AKK08 FFs at $\mu=1.0p_{\rm T}$, and
$\hat{q}_{0}=1.2$ GeV2/fm with DSS FFs at $\mu=1.2p_{\rm T}$. The difference
between $\hat{q}_{0}=1.2\sim 1.4$ GeV${}^{2}/$fm is narrow relative to the
same scale case.

Fig. 16: The $\chi^{2}/d.o.f$ results from fitting to nuclear modification
factors in central $Au+Au$ collisions at $\sqrt{s_{\rm NN}}=200$ GeV, (a) the
fits to only single hadron $R_{AA}(p_{\rm T})$, (b) the fits to only dihadron
$I_{AA}(z_{\rm T})$, and (c) the global fits to $R_{AA}(p_{\rm
T})+I_{AA}(z_{\rm T})$. 5 sets of fragmentation function parameterizations are
used in the model with different scales, respectively.

Fig. 17: The $\chi^{2}/d.o.f$ results from fitting to nuclear modification
factors in central $Pb+Pb$ collisions at $\sqrt{s_{\rm NN}}=2.76$ TeV, (a) the
fitts to only single hadron $R_{AA}(p_{\rm T})$ Abelev  _et al._ (2013);
Chatrchyan  _et al._ (2012), (b) the fits to only dihadron $I_{AA}(p_{\rm
T}^{\rm assoc})$ Aamodt  _et al._ (2012); Conway (2013), and (c) the global
fits to $R_{AA}(p_{\rm T})+I_{AA}(p_{\rm T}^{\rm assoc})$. 6 sets of
fragmentation function parameterizations are used in the model with different
scales, respectively.

Fig. 17 shows the $\chi^{2}/d.o.f$ fits to nuclear modification factors in
central $Pb+Pb$ collisions at $\sqrt{s_{\rm NN}}=2.76$ TeV: (a) the fits to
only single hadron $R_{AA}(p_{\rm T})$, (b) the fits to only dihadron
$I_{AA}(p_{\rm T}^{\rm assoc})$, and (c) the global fits to $R_{AA}(p_{\rm
T})$ \+ $I_{AA}(p_{\rm T}^{\rm assoc})$. We use 6 sets of fragmentation
function parameterizations with different scales in our calculations,
respectively. We get the best fitting values of the jet transport parameter
as: $\hat{q}_{0}=2.2$ GeV${}^{2}/$fm with KRE FFs at $\mu=1.5p_{\rm T}$,
$\hat{q}_{0}=2.4$ GeV${}^{2}/$fm with KKP FFs at $\mu=3.5p_{\rm T}$,
$\hat{q}_{0}=2.3$ GeV${}^{2}/$fm with BFGW FFs at $\mu=5.0p_{\rm T}$,
$\hat{q}_{0}=2.4$ GeV${}^{2}/$fm with HKNS FFs at $\mu=6.5p_{\rm T}$,
$\hat{q}_{0}=2.5$ GeV${}^{2}/$fm with AKK08 FFs at $\mu=7.5p_{\rm T}$,
$\hat{q}_{0}=2.4$ GeV${}^{2}/$fm with DSS FFs at $\mu=5.5p_{\rm T}$. The
difference between $\hat{q}_{0}$ drops to $2.2\sim 2.5$ GeV${}^{2}/$fm.

The above transport coefficients $\hat{q}/T^{3}$ extracted in the model with
different scales in different sets of FFs from hadron suppressions at RHIC and
the LHC are summarized in the right panel of Fig. 13.

###  V.2 Characteristics of fragmentation function contributions for parton to
hadron with different scales for different sets of FFs

Shown in Fig. 18 are the quark and gluon FFs of $\pi^{0}$ hadrons at $p_{\rm
T}=10$ GeV and 50 GeV for all the available FFs with different scales,
respectively. Fig. 19 are the similar plot for charged-hadron FFs of quark and
gluon with different scales in each FFs. The difference between gluon FFs are
a bit reduced with befitted scales compared with the same scale case.

Fig. 18: Comparisons of parton-to-$\pi^{0}$ fragmentation functions between 5
sets of FFs: left panels are for quark FFs, and right panels for gluon FFs.
Upper panels are for the hadrons with $p_{\rm T}=10$ GeV, and lower panels
with $p_{\rm T}=50$ GeV. The scale $\mu$ is different for each FFs. Fig. 19:
Comparisons of parton-to-$h^{\pm}$ FFs between 6 sets of FFs: left panels are
for quark FFs, and right panels for gluon FFs. Upper panels are for the
hadrons with $p_{\rm T}=10$ GeV, and lower panels with $p_{\rm T}=50$ GeV. The
scale $\mu$ is different for each FFs. Fig. 20: The contribution fractions of
quark (solid lines) and gluon (dashed lines) fragmentations to the inclusive
$\pi^{0}$ cross sections given by the model with only AKK08 FFs with
$\mu=1.0p_{\rm T}$ and $\mu=1.2p_{\rm T}$ at $\sqrt{s_{\rm NN}}=200$ GeV,
respectively. Fig. 21: The contribution fractions of quark (solid lines) and
gluon (dashed lines) fragmentations to the inclusive charged-hadron cross
sections given by the model with only AKK08 FFs with $\mu=1.5p_{\rm T}$ and
$\mu=7.5p_{\rm T}$ at $\sqrt{s_{\rm NN}}=2.76$ TeV, respectively.

Taking the AKK08 FFs as examples, we show the relative contributions of quark
and gluon to hadrons with the changed scales in $p+p$ collisions at RHIC and
the LHC energies, as shown in Fig. 20 and 21. One can see that both at RHIC
and the LHC, with $\mu$ increasing, the gluon contribution to final state
hadrons will decrease, thus we need a larger jet quenching parameter to
compensate for the total jet energy loss.

For a given set of FFs in a NLO pQCD parton model, the different scales leads
to the different fractions of gluon (quark) contributions to hadrons in $p+p$
collisions. Although the scale change can also affect the parton distribution
functions and the hard cross sections, the fraction change is mainly
contributed by parton fragmentation functions, as shown in Fig. 20 and 21. In
details, when the scale $\mu$ decreases from 1.2$p_{\rm T}$ to 1.0$p_{\rm T}$
at $\sqrt{s_{\rm NN}}=200$ GeV, the contribution of gluon-to-hadron becomes
larger, as illustrated in Fig. 20, so a relatively smaller $\hat{q}_{0}$ is
needed for the case of $\mu=1.0p_{\rm T}$ shown in Fig. 13. Meanwhile, with
the $\mu$ increasing from 1.5$p_{\rm T}$ to 7.5$p_{\rm T}$ at $\sqrt{s_{\rm
NN}}=2.76$ TeV, the contribution of gluon-to-hadron reduces, as presented in
Fig. 21, thus a relatively larger $\hat{q}_{0}$ is needed for the case of
$\mu=7.5p_{\rm T}$ shown in Fig. 13. In a word, the different fraction of
gluon (quark) contribution to hadrons will give different energy loss
parameters due to gluon energy loss being $9/4$ times of quark energy loss.

## References

  * Adare  _et al._ (2008) A. Adare _et al._ (PHENIX), Phys. Rev. Lett. 101, 232301 (2008), arXiv:0801.4020 [nucl-ex] .
  * Adare  _et al._ (2013) A. Adare _et al._ (PHENIX), Phys. Rev. C 87, 034911 (2013), arXiv:1208.2254 [nucl-ex] .
  * Adare  _et al._ (2010) A. Adare _et al._ (PHENIX), Phys. Rev. Lett. 105, 142301 (2010), arXiv:1006.3740 [nucl-ex] .
  * Abelev  _et al._ (2013) B. Abelev _et al._ (ALICE), Phys. Lett. B 720, 52 (2013), arXiv:1208.2711 [hep-ex] .
  * Chatrchyan  _et al._ (2012) S. Chatrchyan _et al._ (CMS), Eur. Phys. J. C 72, 1945 (2012), arXiv:1202.2554 [nucl-ex] .
  * Acharya  _et al._ (2018) S. Acharya _et al._ (ALICE), JHEP 11, 013 (2018), arXiv:1802.09145 [nucl-ex] .
  * Khachatryan  _et al._ (2017) V. Khachatryan _et al._ (CMS), JHEP 04, 039 (2017), arXiv:1611.01664 [nucl-ex] .
  * Acharya  _et al._ (2019) S. Acharya _et al._ (ALICE), Phys. Lett. B 788, 166 (2019), arXiv:1805.04399 [nucl-ex] .
  * Gyulassy and Plumer (1990) M. Gyulassy and M. Plumer, Phys. Lett. B 243, 432 (1990).
  * Wang and Gyulassy (1992) X.-N. Wang and M. Gyulassy, Phys. Rev. Lett. 68, 1480 (1992).
  * Qin and Wang (2015) G.-Y. Qin and X.-N. Wang, Int. J. Mod. Phys. E 24, 1530014 (2015), arXiv:1511.00790 [hep-ph] .
  * Zhang  _et al._ (2007) H. Zhang, J. F. Owens, E. Wang, and X.-N. Wang, Phys. Rev. Lett. 98, 212301 (2007), arXiv:nucl-th/0701045 .
  * Baier  _et al._ (1997) R. Baier, Y. L. Dokshitzer, A. H. Mueller, S. Peigne, and D. Schiff, Nucl. Phys. B 484, 265 (1997), arXiv:hep-ph/9608322 .
  * Burke  _et al._ (2014) K. M. Burke _et al._ (JET), Phys. Rev. C 90, 014909 (2014), arXiv:1312.5003 [nucl-th] .
  * Liu  _et al._ (2016) Z.-Q. Liu, H. Zhang, B.-W. Zhang, and E. Wang, Eur. Phys. J. C 76, 20 (2016), arXiv:1506.02840 [nucl-th] .
  * Chen  _et al._ (2010) X.-F. Chen, C. Greiner, E. Wang, X.-N. Wang, and Z. Xu, Phys. Rev. C 81, 064908 (2010), arXiv:1002.1165 [nucl-th] .
  * Chen  _et al._ (2011) X.-F. Chen, T. Hirano, E. Wang, X.-N. Wang, and H. Zhang, Phys. Rev. C 84, 034902 (2011), arXiv:1102.5614 [nucl-th] .
  * Zhang  _et al._ (2009) H. Zhang, J. F. Owens, E. Wang, and X.-N. Wang, Phys. Rev. Lett. 103, 032302 (2009), arXiv:0902.4000 [nucl-th] .
  * Xie  _et al._ (2019) M. Xie, S.-Y. Wei, G.-Y. Qin, and H.-Z. Zhang, Eur. Phys. J. C 79, 589 (2019), arXiv:1901.04155 [hep-ph] .
  * Xie  _et al._ (2021) M. Xie, X.-N. Wang, and H.-Z. Zhang, Phys. Rev. C 103, 034911 (2021), arXiv:2003.02441 [hep-ph] .
  * Feal  _et al._ (2021) X. Feal, C. A. Salgado, and R. A. Vazquez, Phys. Lett. B 816, 136251 (2021), arXiv:1911.01309 [hep-ph] .
  * Liu  _et al._ (2021) F.-L. Liu, W.-J. Xing, X.-Y. Wu, G.-Y. Qin, S. Cao, and X.-N. Wang, (2021), arXiv:2107.11713 [hep-ph] .
  * Cao  _et al._ (2021) S. Cao _et al._ (JETSCAPE), Phys. Rev. C 104, 024905 (2021), arXiv:2102.11337 [nucl-th] .
  * d’Enterria  _et al._ (2014) D. d’Enterria, K. J. Eskola, I. Helenius, and H. Paukkunen, Nucl. Phys. B 883, 615 (2014), arXiv:1311.1415 [hep-ph] .
  * Deng and Wang (2010) W.-t. Deng and X.-N. Wang, Phys. Rev. C 81, 024902 (2010), arXiv:0910.3403 [hep-ph] .
  * Kretzer (2000) S. Kretzer, Phys. Rev. D 62, 054001 (2000), arXiv:hep-ph/0003177 .
  * Kniehl  _et al._ (2000) B. A. Kniehl, G. Kramer, and B. Potter, Nucl. Phys. B 582, 514 (2000), arXiv:hep-ph/0010289 .
  * Bourhis  _et al._ (2001) L. Bourhis, M. Fontannaz, J. P. Guillet, and M. Werlen, Eur. Phys. J. C 19, 89 (2001), arXiv:hep-ph/0009101 .
  * Hirai  _et al._ (2007) M. Hirai, S. Kumano, T. H. Nagai, and K. Sudoh, Phys. Rev. D 75, 094009 (2007), arXiv:hep-ph/0702250 .
  * Albino  _et al._ (2008) S. Albino, B. A. Kniehl, and G. Kramer, Nucl. Phys. B 803, 42 (2008), arXiv:0803.2768 [hep-ph] .
  * de Florian  _et al._ (2007a) D. de Florian, R. Sassot, and M. Stratmann, Phys. Rev. D 75, 114010 (2007a), arXiv:hep-ph/0703242 .
  * de Florian  _et al._ (2007b) D. de Florian, R. Sassot, and M. Stratmann, Phys. Rev. D 76, 074033 (2007b), arXiv:0707.1506 [hep-ph] .
  * Wang  _et al._ (1996) X.-N. Wang, Z. Huang, and I. Sarcevic, Phys. Rev. Lett. 77, 231 (1996), arXiv:hep-ph/9605213 .
  * Wang and Huang (1997) X.-N. Wang and Z. Huang, Phys. Rev. C 55, 3047 (1997), arXiv:hep-ph/9701227 .
  * Wang (2004a) X.-N. Wang, Phys. Rev. C 70, 031901 (2004a), arXiv:nucl-th/0405029 .
  * Wang and Wang (2001) E. Wang and X.-N. Wang, Phys. Rev. Lett. 87, 142301 (2001), arXiv:nucl-th/0106043 .
  * Wang and Wang (2002) E. Wang and X.-N. Wang, Phys. Rev. Lett. 89, 162301 (2002), arXiv:hep-ph/0202105 .
  * Owens (1987) J. F. Owens, Rev. Mod. Phys. 59, 465 (1987).
  * Brock  _et al._ (1995) R. Brock _et al._ (CTEQ), Rev. Mod. Phys. 67, 157 (1995).
  * Hou  _et al._ (2021) T.-J. Hou _et al._ , Phys. Rev. D 103, 014013 (2021), arXiv:1912.10053 [hep-ph] .
  * Kidonakis and Owens (2001) N. Kidonakis and J. F. Owens, Phys. Rev. D 63, 054019 (2001), arXiv:hep-ph/0007268 .
  * Harris and Owens (2002) B. W. Harris and J. F. Owens, Phys. Rev. D 65, 094032 (2002), arXiv:hep-ph/0102128 .
  * Jacobs and Cooper (2000) P. Jacobs and G. Cooper, (2000), arXiv:nucl-ex/0008015 .
  * Wang (1997) X.-N. Wang, Phys. Rept. 280, 287 (1997), arXiv:hep-ph/9605214 .
  * Li and Wang (2002) S.-y. Li and X.-N. Wang, Phys. Lett. B 527, 85 (2002), arXiv:nucl-th/0110075 .
  * Emel’yanov  _et al._ (2000) V. Emel’yanov, A. Khodinov, S. R. Klein, and R. Vogt, Phys. Rev. C 61, 044904 (2000), arXiv:hep-ph/9909427 .
  * Hirano and Nara (2004) T. Hirano and Y. Nara, Phys. Rev. C 69, 034908 (2004), arXiv:nucl-th/0307015 .
  * Eskola  _et al._ (2017) K. J. Eskola, P. Paakkinen, H. Paukkunen, and C. A. Salgado, Eur. Phys. J. C 77, 163 (2017), arXiv:1612.05741 [hep-ph] .
  * Chang  _et al._ (2014) N.-B. Chang, W.-T. Deng, and X.-N. Wang, Phys. Rev. C 89, 034911 (2014), arXiv:1401.5109 [nucl-th] .
  * Hirano (2002) T. Hirano, Phys. Rev. C 65, 011901 (2002), arXiv:nucl-th/0108004 .
  * Hirano and Tsuda (2002) T. Hirano and K. Tsuda, Phys. Rev. C 66, 054905 (2002), arXiv:nucl-th/0205043 .
  * Wang (2004b) X.-N. Wang, Phys. Lett. B 579, 299 (2004b), arXiv:nucl-th/0307036 .
  * Adamczyk  _et al._ (2016) L. Adamczyk _et al._ (STAR), Phys. Lett. B 760, 689 (2016), arXiv:1604.01117 [nucl-ex] .
  * Aamodt  _et al._ (2012) K. Aamodt _et al._ (ALICE), Phys. Rev. Lett. 108, 092301 (2012), arXiv:1110.0121 [nucl-ex] .
  * Conway (2013) R. Conway (CMS), Nucl. Phys. A 904-905, 451c (2013).
  * Zhang  _et al._ (2008) H.-z. Zhang, J. F. Owens, E. Wang, and X. N. Wang, J. Phys. G 35, 104067 (2008), arXiv:0804.2381 [hep-ph] .
  * Song and Heinz (2008) H. Song and U. W. Heinz, Phys. Rev. C 77, 064901 (2008), arXiv:0712.3715 [nucl-th] .
  * Qiu  _et al._ (2012) Z. Qiu, C. Shen, and U. Heinz, Phys. Lett. B 707, 151 (2012), arXiv:1110.3033 [nucl-th] .