# Evaluating the blast-wave model as a description of 5 TeV p-Pb $\bf p_{t}$
spectra

Thomas A. Trainor  University of Washington, Seattle, WA 98195

###### Abstract

The blast-wave (BW) spectrum model is interpreted to reveal relativistic
motion (collective flow) of the hadron emission system relative to the center-
of-momentum (CM) frame in high-energy A-B collisions. In essence, any spectrum
deviation in the CM frame from a reference distribution (e.g. Boltzmann
distribution on transverse mass $m_{t}$) is interpreted to reveal a flowing
particle source. The ALICE collaboration has applied the BW model to
identified hadron (PID) spectra for four hadron species from 5 TeV $p$-Pb
collisions. From model fits BW parameters $T_{kin}$ (freeze-out temperature)
and $\langle\beta_{t}\rangle$ (transverse speed) are inferred that suggest
strong radial expansion in more-central $p$-Pb collisions. Such results from
the small $p$-Pb collision system are counterintuitive given that strong
radial expansion should be driven by large density gradients. The present
study is intended to address that problem. Several methods are employed to
evaluate the quality of the BW model data description, including logarithmic
derivatives and the Z-score statistic. The stability of the BW model
definition across several applications to data is investigated. The BW model
data description is compare to that of the two-component (soft+hard) model
(TCM) that has been previously applied to the same $p$-Pb PID spectra. The
general conclusion is that the BW model is falsified by $p$-Pb PID spectrum
data according to standard statistical measures and that the fitted parameter
values do not convey the intended meaning. Statistically acceptable data
descriptions provided by the TCM indicate that other collision mechanisms
(projectile-nucleon dissociation, dijet production), that are consistent with
conventional QCD, are more likely responsible for observed spectrum
characteristics.

††preprint: version 1.8

##  I Introduction

The blast-wave (BW) spectrum model has been widely applied to identified-
hadron (PID) spectra from A-B collisions in the context of the high-energy
heavy ion program whose central goal has been formation of a quark-gluon
plasma (QGP). The basic premise is that BW spectrum analysis might confirm the
presence of radial flow – radial expansion of a conjectured fluid assumed to
be the dominant source of final-state particle emission. Initial application
of the BW model to SPS S-S collisions (with $\sqrt{s_{NN}}\approx 19$ GeV)
proved inconclusive ssflow . However, early results from the relativistic
heavy ion collider (RHIC) ($\sqrt{s_{NN}}\approx 130$ and 200 GeV) were
interpreted to indicate the presence of significant radial flow (and hence
QGP?), at least in more-central Au-Au collisions phenradflow .

In motivating the RHIC experimental program it had been assumed that only for
more-central A-A collisions involving the largest nuclei might the requisite
energy and matter densities be generated to form a QGP. However, application
of the BW model to smaller collision systems (e.g. $p$-A, $d$-A and even
$p$-$p$) has lead to claims that significant radial flow appears in all A-B
collisions starradflow . A recent example of the latter is BW analysis of PID
spectra from 5 TeV $p$-Pb collisions aliceppbpid . Given the novelty and
significance of such claims several questions emerge regarding BW model fits
to PID spectra:

How well-defined is the theoretical basis for the BW model? That question goes
beyond the basic Cooper-Frye formulation cooper  relating particle emission in
a relativistic moving frame to observed particle momenta in the center-of-
momentum (CM) or laboratory (lab) frame. It concerns basic assumptions about
the particle emission process – what does “particle emission” mean – and what
constraints does the emission environment (presumably high-density matter of
uncertain nature) impose? Also, how consistent is the BW model formulation
from one application to another? Can different research groups apply the same
model to similar data and get the same result, or is each application unique
and therefore unverifiable?

Regarding BW model fits to data what are the standards for concluding that a
BW model (which one?) describes spectrum data with sufficient accuracy to
confirm a flowing particle-source scenario? Does a well-defined BW model
describe all available relevant data within statistical uncertainties. Is
there a competing spectrum model that describes spectrum data qualitatively
more accurately? Within the context of prevailing practice for applying the
model is it possible to falsify the BW model?

To address those questions and other related considerations the following
strategy is adopted: A reference model for comparison with BW results is
introduced in the form of a two-component (soft + hard) model (TCM) applied to
PID spectrum data from 5 TeV $p$-Pb collisions pidpart1 ; pidpart2 . A BW
model description of the same PID data, as reported in Ref. aliceppbpid , is
reviewed. Several model-independent measures are employed to compare BW and
TCM data descriptions. Limiting cases of the BW model are compared with a
Boltzmann distribution (assuming thermal emission from a stationary particle
source) and the TCM soft component (a limiting case corresponding to zero
particle density and no jet contribution to spectra). Model fit quality is
assessed with Z-scores based on published statistical and systematic data
uncertainties.

To provide fuller context for model comparisons the concept of “two cultures”
is introduced relating an approach based on perturbative QCD and parton
cascades or showers on the one hand and an approach based on fluid-dynamic
descriptions emerging from preQCD collision models on the other. Finally, the
concept of “elementary collisions” as a reference is briefly considered with
reference to particle data from $p$-$p$ and $e^{+}$-$e^{-}$ collisions.
Appendix A defines relativistic quantities appearing in BW derivations and
reviews several BW model versions.

The goals of this study are: (a) understand the assumptions and formulations
of the BW model as manifested by several versions, (b) assess whether
underlying BW model assumptions are physically reasonable given present
knowledge of QCD theory and improved data quality, (c) determine whether the
BW model describes available data within uncertainties (statistical rather
than systematic) or is falsified by data, and (d) consider an alternative
spectrum model that appears superior as a data model and more physically
interpretable.

This article is arranged as follows: Section II reviews a TCM description of
PID spectrum data from 5 TeV $p$-Pb collisions reported in Ref. aliceppbpid .
Section III briefly presents the BW model as it is applied to the same PID
spectrum data. Section IV compares the shapes of the two spectrum models via
several model-independent shape measures. Section V evaluates data fit quality
for BW and TCM models using the Z-score statistical measure. Section VI
provides a survey of BW model evolution over several decades. Section VII
reviews systematic uncertainties. Sections VIII and IX present discussion and
summary. Appendix A presents a detailed review of relativistic quantities and
several BW model derivations.

##  II TCM $\bf vs$ $\bf p$-$\bf Pb$ PID Spectrum data

The two-component model (TCM) of PID hadron spectra described here provides a
simple reference based on conventional soft and hard QCD processes that may
furnish a better understanding of the BW model and its relation to spectrum
data and collision dynamics. The PID spectrum data for seven charge-
multiplicity $n_{ch}$ classes from 5 TeV $p$-Pb collisions are as reported in
Ref. aliceppbpid . Previous TCM analysis of those data is reported in Refs.
ppbpid ; pidpart1 ; pidpart2 . A version of the PID TCM in which two
parameters are varied with $n_{ch}$ to accommodate spectrum data describes all
data within their statistical uncertainties as demonstrated in Ref. pidpart2 .

###  II.1 p-Pb spectrum TCM for identified hadrons

The $y_{t}$ spectrum TCM for unidentified hadrons is ppbpid

$\displaystyle\bar{\rho}_{0}(y_{t},n_{s})$ $\displaystyle\approx$
$\displaystyle\frac{d^{2}n_{ch}}{y_{t}dy_{t}d\eta}$ $\displaystyle=$
$\displaystyle S(y_{t},n_{s})+H(y_{t},n_{s})$ $\displaystyle\approx$
$\displaystyle\bar{\rho}_{s}(n_{s})\hat{S}_{0}(y_{t})+\bar{\rho}_{h}(n_{s})\hat{H}_{0}(y_{t}),$

where $y_{t}\equiv\ln[(m_{t}+p_{t})/m_{0}]$ is transverse rapidity with
$m_{0}\rightarrow m_{\pi}$ assumed. $\hat{S}_{0}(y_{t})$ and
$\hat{H}_{0}(y_{t})$ are unit-integral fixed model functions (independent of
A-B centrality or charge multiplicity $n_{ch}$), and $\bar{\rho}_{s}(n_{s})$
and $\bar{\rho}_{h}(n_{s})$ are soft and hard charge densities with
$\bar{\rho}_{x}\equiv n_{x}/\Delta\eta$ and
$\bar{\rho}_{0}=\bar{\rho}_{s}+\bar{\rho}_{h}$. Soft yield $n_{s}$ serves as a
centrality index. A-B charge densities are defined in terms of N-N quantities
by $\bar{\rho}_{s}=(N_{part}/2)\bar{\rho}_{sNN}$ and
$\bar{\rho}_{h}=N_{bin}\bar{\rho}_{hNN}$ with $x\equiv\rho_{hNN}/\rho_{sNN}$
and $\nu\equiv 2N_{bin}/N_{part}$ and with
$\bar{\rho}_{h}/\bar{\rho}_{s}=x(n_{s})\nu(n_{s})$. For linear superposition
of N-N collisions in A-B collisions
$\bar{\rho}_{hNN}\approx\alpha(\sqrt{s_{NN}})\bar{\rho}_{sNN}^{2}$ ppprd ;
ppbpid .

Given the A-B spectrum TCM for unidentified hadron spectra a corresponding TCM
for identified hadrons can be generated by assuming that each hadron species
$i$ comprises certain fractions of soft and hard TCM components denoted by
$z_{si}(n_{s})$ and $z_{hi}(n_{s})$ assumed independent of $y_{t}$. The PID
spectrum TCM is then expressed as pidpart1 ; pidpart2

$\displaystyle\bar{\rho}_{0i}(y_{t},n_{s})$ $\displaystyle=$ $\displaystyle
S_{i}(y_{t},n_{s})+H_{i}(y_{t},n_{s})$ $\displaystyle\approx$ $\displaystyle
z_{si}(n_{s})\bar{\rho}_{s}\hat{S}_{0i}(y_{t})+z_{hi}(n_{s})\bar{\rho}_{h}\hat{H}_{0i}(y_{t},n_{s})$
$\displaystyle\frac{\bar{\rho}_{0i}(y_{t},n_{s})}{z_{si}(n_{s})\bar{\rho}_{s}}$
$\displaystyle\equiv$ $\displaystyle X_{i}(y_{t},n_{s})$
$\displaystyle\approx$
$\displaystyle\hat{S}_{0i}(y_{t})+\tilde{z}_{i}(n_{s})x(n_{s})\nu(n_{s})\hat{H}_{0i}(y_{t}),$

where $\tilde{z}_{i}(n_{s})\equiv z_{hi}(n_{s})/z_{si}(n_{s})$ and unit-
integral model functions $\hat{S}_{0i}(y_{t})$ and $\hat{H}_{0i}(y_{t})$
depend on hadron species $i$. Soft fraction $z_{si}(n_{s})$ is expressed
within the TCM by

$\displaystyle z_{si}(n_{s})$ $\displaystyle=$
$\displaystyle\left[\frac{1+x(n_{s})\nu(n_{s})}{1+\tilde{z}_{i}(n_{s})x(n_{s})\nu(n_{s})}\right]{z_{0i}},$
(4)

where $z_{hi}(n_{s})=\tilde{z}_{i}(n_{s})z_{si}(n_{s})$ and
$\bar{\rho}_{0i}\equiv z_{0i}\bar{\rho}_{0}$ defines $z_{0i}$. Thus, if
$\tilde{z}_{i}(n_{s})$ and $z_{0i}$ are specified for relevant hadron species
then all of the PID TCM is determined. Model functions $\hat{S}_{0i}(y_{t})$
are defined on proper $m_{ti}$ for a given hadron species $i$ and then
transformed to $y_{t\pi}$. $\hat{H}_{0i}(y_{t})$ are always defined on
$y_{t\pi}$. The basis for the $\hat{S}_{0i}(y_{t})$ model definitions is the
limit of normalized spectra $X_{i}(y_{t},n_{s})$ as $n_{ch}\rightarrow 0$
(zero particle density). Soft-component models are thus always derived from
data spectra. Inferred data spectrum hard components are then by definition
the complement of model soft components so defined.

Properties of TCM soft-component models $\hat{S}_{0i}(m_{t})$ are consistent
with nuclear transparency in $p$-A collisions busza  and A-B collisions
stopping  and with the wounded-nucleon model bialas . Properties of TCM hard-
component models $\hat{H}_{0i}(y_{t})$ are consistent with measured properties
of minimum-bias jets within A-B collisions if they consist of linear
superpositions of N-N collisions hardspec ; fragevo ; eeprd ; jetspec2 ;
mbdijets .

###  II.2 5 TeV $\bf p$-Pb TCM model parameters

Table 1 presents TCM geometry parameters for 5 TeV $p$-Pb collisions inferred
from the analysis in Ref. tomglauber . Those geometry parameters, derived from
$p$-Pb $p_{t}$ spectrum and ensemble $\bar{p}_{t}$ data for unidentified
hadrons alicempt ; tommpt ; tomglauber , are assumed valid for each
identified-hadron species (confirmed in Ref. ppbpid ).
$\bar{\rho}_{0}=n_{ch}/\Delta\eta$ charge densities are measured quantities
inferred from Fig. 16 of Ref. aliceglauber . Relations $N_{bin}=N_{part}+1$
and $\nu=2N_{bin}/N_{part}$ involve number of nucleon N participants and N-N
binary collisions. $\bar{\rho}_{sNN}$ is the mean soft-component charge
density per participant pair averaged over all pairs.
$x\equiv\bar{\rho}_{hNN}/\bar{\rho}_{sNN}\approx\alpha\bar{\rho}_{sNN}$ is the
N-N hard/soft density ratio. Column $\sigma^{\prime}/\sigma_{0}$ presents
nominal centralities (bin centers) reported by Ref. aliceppbpid  in connection
with measured charge densities $\bar{\rho}_{0}$ whereas column
$\sigma/\sigma_{0}$ presents values inferred in Ref. tomglauber . The
remaining values in the table are results of the latter analysis.

Table 1: TCM fractional cross sections $\sigma/\sigma_{0}$ (bin centers) and geometry parameters, midrapidity charge density $\bar{\rho}_{0}$, TCM N-N soft component $\bar{\rho}_{sNN}$ and N-N hard/soft ratio $x(n_{s})$ used for 5 TeV $p$-Pb PID spectrum analysis tomglauber . Centrality parameters are from Ref. tomglauber . $\sigma^{\prime}/\sigma_{0}$ values are from Table 1 of Ref. aliceppbpid  as determined in Ref. aliceglauber .  $n$ | $\sigma^{\prime}/\sigma_{0}$ | $\sigma/\sigma_{0}$ | $N_{bin}$ | $\nu$ | $\bar{\rho}_{0}$ | $\bar{\rho}_{sNN}$ | $x(n_{s})$
---|---|---|---|---|---|---|---
1 | 0.025 | 0.15 | 3.20 | 1.52 | 44.6 | 16.6 | 0.188
2 | 0.075 | 0.24 | 2.59 | 1.43 | 35.9 | 15.9 | 0.180
3 | 0.15 | 0.37 | 2.16 | 1.37 | 30.0 | 15.2 | 0.172
4 | 0.30 | 0.58 | 1.70 | 1.26 | 23.0 | 14.1 | 0.159
5 | 0.50 | 0.80 | 1.31 | 1.13 | 15.8 | 12.1 | 0.137
6 | 0.70 | 0.95 | 1.07 | 1.03 | 9.7 | 8.7 | 0.098
7 | 0.90 | 0.99 | 1.00 | 1.00 | 4.4 | 4.2 | 0.047

Table 2 presents hard-component model parameters based on optimized
descriptions of spectrum hard components for the most-central ($n=1$) $p$-Pb
event class as described in Ref. pidpart1 . TCM soft-component parameters
remain unchanged from those reported in Table II of Ref. pidpart1 , with
$(T,n)=$ (145 MeV,8.5), (200 MeV,14) and (210 MeV,14) for pions, kaons and
baryons respectively. More-massive hadrons are less sensitive to exponent $n$.
Those parameters and these in Tables 1 and 2 define a fixed PID TCM
established as a reference. In Ref. pidpart2  a variable TCM is defined
wherein certain hard-component parameters (Gaussian widths $\sigma_{y_{t}}$
below or above the hard-component mode or centroids $\bar{y}_{t}$) are varied
linearly with centrality measure $x\nu$ to accommodate data (see Fig. 4 of
Ref. pidpart2 ). The variable TCM describes all spectra within statistical
uncertainties as in Fig. 1.

Table 2: Revised PID TCM hard-component model parameters $(\bar{y}_{t},\sigma_{y_{t}},q)$ for identified hadrons from 5 TeV $p$-Pb collisions derived from the differential spectrum analysis in Rev. pidpart1 . Parameter $z_{0i}$ values, inferred as limiting values of $z_{si}(n_{s})$ centrality trends (Fig. 9 of Ref. pidpart1 ), are also included. Uncertainties are determined as one half the parameter change that would produce an obvious variation in $z_{hi}(y_{t},n_{s})$ ratios (e.g. Fig. 3 of Ref. pidpart2 ). Values with no uncertainties are duplicated from a related particle type.  | $\bar{y}_{t}$ | $\sigma_{y_{t}}$ | $q$ | $z_{0i}$
---|---|---|---|---
$\pi^{\pm}$ | $2.46\pm 0.005$ | $0.575\pm 0.005$ | $4.1\pm 0.5$ | 0.82 $\pm$0.01
$K^{\pm}$ | $2.655$ | $0.568$ | $4.1$ | 0.128 $\pm$0.002
$K_{\text{S}}^{0}$ | $2.655\pm 0.005$ | $0.568\pm 0.003$ | $4.1\pm 0.1$ | 0.064 $\pm$0.002
$p$ | $2.99\pm 0.005$ | $0.47\pm 0.005$ | $5.0$ | 0.065 $\pm$0.002
$\Lambda$ | $2.99\pm 0.005$ | $0.47\pm 0.005$ | $5.0\pm 0.05$ | 0.034 $\pm$0.002

###  II.3 5 TeV $\bf p$-Pb PID spectrum data

Identified-hadron spectrum data obtained from Ref. aliceppbpid  for the
present analysis were produced by the ALICE collaboration at the LHC.
Collision events were divided into seven charge-multiplicity $n_{ch}$ or
$p$-Pb centrality classes based on yields in a VZERO-A (V0A) counter
subtending $2.8<\eta_{lab}<5.1$ in the Pb direction. Hadron species include
charged pions $\pi^{\pm}$, charged kaons $K^{\pm}$, neutral kaons
$K^{0}_{\text{S}}$, protons $p,~{}\bar{p}$ and Lambdas
$\Lambda,~{}\bar{\Lambda}$.

Figure 1 shows PID spectrum data (densities on $p_{t}$) from Ref. aliceppbpid
(points) plotted vs logarithmic variable $y_{t}$. That format (log-log with
respect to $p_{t}$) provides detailed access to low-$p_{t}$ structure (where
most jet fragments appear) and clearly shows power-law trends at higher
$p_{t}$. The curves are TCM parametrizations as reported above and
demonstrated in Ref. pidpart2  to describe spectrum data within their
statistical uncertainties. As in Ref. aliceppbpid  the spectra are scaled up
by powers of 2 according to $2^{n-1}$ where $n\in[1,7]$ is the centrality
class index and $n=1$ is least central (following the usage in Ref.
aliceppbpid ). In this paper $n=1$ denotes the most central data as in Table
1.

(a)(b)
(c)(d)

Figure 1:  $p_{t}$ spectra for identified hadrons from 5 TeV $p$-Pb collisions
aliceppbpid  plotted vs pion transverse rapidity $y_{t}$ (default) for: (a)
pions, (b) neutral kaons, (c) corrected protons, (d) Lambdas. Solid curves
represent the PID spectrum TCM from Ref. pidpart2 . Proton inefficiency
corrections are described in Ref. pidpart1 . The statistical uncertainties
(0.1 to 1% of data values, e.g. see Fig. 17) are all smaller than the point
size.

Note that in the format of Fig. 1 a Boltzmann exponential on transverse mass
$m_{t}$ appears as $A-m_{i}\cosh(y_{t})/T_{i}$, where $m_{i}$ is the mass for
hadron species $i$ and $A$ is a constant, and a power law on $m_{t}$ appears
as a straight line. The exception is for pion spectra including a substantial
resonance contribution for $y_{t}<2$ ($p_{t}<0.5$ GeV/c) pidpart1 .

##  III BW model $\bf vs$ $\bf p$-$\bf Pb$ PID spectrum data

Reference aliceppbpid  applies a blast-wave model to $p_{t}$ spectrum data for
5 TeV V0A $p$-Pb spectra, with fit parameters reported in its Fig. 6. The BW
fit parameters are mean radial speed
$\langle\beta_{T}\rangle\rightarrow\bar{\beta}_{t}$, “kinetic freezeout”
temperature $T_{kin}$ and radial-speed $\beta_{t}(r)$ profile parameter $n$ as
in Table 5 of Ref. aliceppbpid  (and see Table 3 below).

Application of certain fit models (e.g. BW models) to spectrum data for small
collision systems has been invoked in a number of studies to claim presence of
collectivity (flows) in such systems. For example, in Ref. cleymans
(reporting BW analysis of unidentified hadrons from $p$-$p$ collisions) the BW
model is said to be a “more standard description” and is “based on collective
flow in small systems.” The BW model is said to be “quite good in explaining
the bulk part of the system, however it fails at low-$p_{T}$ below 0.5 GeV/c
which could possibly be due to the decays of hadronic resonances.” “The
applicability of [the BW model] is verified by fitting the transverse momentum
spectra of the bulk part ($\sim$ 2.5 GeV/c)….”

Table 3: Blast-wave parameters for simultaneous fits of pion, charged-kaon, neutral kaon, proton and Lambda spectra from 5 TeV $p$-Pb collisions aliceppbpid . The $\sigma^{\prime}/\sigma_{0}$ values are the cross-section fractions (bin centers) reported in Ref. aliceppbpid . Table 1 shows alternative $\sigma/\sigma_{0}$ values determined in Ref. tomglauber .  $n$ | $\sigma^{\prime}/\sigma_{0}$ | $\langle\beta_{T}\rangle$ |  $T_{kin}$ (GeV) | $n$ |  $\chi^{2}$/ndf
---|---|---|---|---|---
1 | 0.025 | 0.547 | 0.143 | 1.07 | 0.27
2 | 0.075 | 0.531 | 0.147 | 1.14 | 0.33
3 | 0.15 | 0.511 | 0.151 | 1.24 | 0.36
4 | 0.30 | 0.478 | 0.157 | 1.41 | 0.35
5 | 0.50 | 0.428 | 0.164 | 1.73 | 0.43
6 | 0.70 | 0.36 | 0.169 | 2.4 | 0.54
7 | 0.90 | 0.26 | 0.166 | 3.9 | 0.84

The argument presented in Ref. aliceppbpid  for applying the BW model to PID
spectrum data and interpreting the result in a hydrodynamic context proceeds
as follows:

(a) Low-$p_{t}$ regions of hadron spectra convey important information: “The
$p_{T}$ distributions and yields of particles of different mass at low and
intermediate momenta of less than a few GeV/c (where the vast majority of
particles is produced) can provide important information about the system
created in high-energy hadron reactions.”

(b) Hydrodynamics-based spectrum models (e.g. blast-wave model) provide useful
data descriptions: “The measured $p_{T}$ distributions are compared
to…hydrodynamic models.” “Several [spectrum] parametrizations have been
tested, among which the blast-wave function…gives the best description of the
data over the full $p_{T}$ range… [model applied individually to different
hadron species].”

(c) Spectrum trends for heavy-ion collisions appear to buttress validity of a
blast-wave spectrum model: “In heavy-ion collisions, the flattening of
transverse momentum distribution and its mass ordering find their natural
explanation in the collective radial expansion of the system. This picture can
be tested in a blast-wave framework with a simultaneous fit to all particles
for each multiplicity bin. This parameterization assumes a locally thermalized
medium, expanding collectively with a common velocity field and undergoing an
instantaneous common freeze-out.” “This [collective hydrodynamic flow] results
in a characteristic dependence of the [spectrum] shape which can be described
with a common kinetic freeze-out temperature parameter $T_{kin}$ and a
collective average expansion velocity $\langle\beta_{t}\rangle$ [citing Ref.
ssflow ].”

In referring to simultaneous BW fits to PID spectra Ref. aliceppbpid  presents
no BW-data fit residuals to support its very low $\chi^{2}$ values (as
reproduced in Table 3). Claimed good agreement between individual BW fits and
data appears to extend (for neutral kaons and Lambdas) up to 7 GeV/c where one
might expect parton fragmentation to jets to produce essentially all detected
hadrons. Why should the BW model have any relation to that $p_{t}$ region?

The BW model used in Ref. aliceppbpid  to fit spectrum data is adopted from
Ref. ssflow  that introduced a hydrodynamics-based formula to describe pion
spectra from 200 GeV fixed-target ($\sqrt{s_{NN}}\approx 19$ GeV) S-S
collisions at the CERN SPS. The relevant formula is Eq. (7) of Ref. ssflow

$\displaystyle\frac{dn}{m_{t}dm_{t}}$ $\displaystyle\propto$ $\displaystyle
m_{t}\int_{0}^{R}rdrI_{0}\left[\frac{p_{t}\sinh(\rho)}{T}\right]K_{1}\left[\frac{m_{t}\cosh(\rho)}{T}\right],~{}~{}~{}~{}$
(5)

with boost $\rho=\tanh^{-1}(\beta_{t})$ and transverse speed
$\beta_{t}(r)=\beta_{s}(r/R)^{n}$, where $n=1$ corresponds to Hubble
expansion, $\beta_{s}$ is the expansion speed at the emitter surface and
$I_{0}$ and $K_{1}$ are modified Bessel functions. The mean transverse speed
is $\langle\beta_{t}\rangle\rightarrow\bar{\beta}_{t}=2\beta_{s}/(n+2)$.
Equation (5) then in effect represents a thermal (Boltzmann) energy spectrum
in the boost (comoving) frame convoluted with a source boost (speed)
distribution on source radius to describe the particle spectrum measured in
the lab frame.

(a)(b)
(c)(d)

Figure 2:  $p_{t}$ spectra for identified hadrons from 5 TeV $p$-Pb collisions
aliceppbpid  plotted vs pion transverse rapidity $y_{t}$ (default) for: (a)
pions, (b) neutral kaons, (c) corrected protons, (d) Lambdas. Solid curves
represent blast-wave model fits from Ref. aliceppbpid  with amplitudes
adjusted here to best describe data within reported BW fit intervals (arrows).

Figure 2 shows the same PID spectrum data from 5 TeV $p$-Pb collisions
appearing in Fig. 1 (points) accompanied in this case by simultaneous BW fits
to four hadron species (solid) with fit values corresponding to Table 3. The
various $n_{ch}$ classes are scaled up by powers of 2 just as in Fig. 1. The
BW model spectra have been scaled in this case to best agree with data spectra
within the indicated BW fit intervals (arrows). Cursory inspection suggests
that the BW description of pion spectra is not relevant to interpretation. The
BW descriptions of spectra for more-massive hadrons appear to accommodate data
reasonably well below $y_{t}$ = 3-3.6 ($p_{t}\approx$ 1.4-2.6 GeV/c), but see
Sec. V for detailed assessment of fit quality.

##  IV Spectrum model comparisons

The TCM provides absolute predictions for hadron yields as well as spectrum
shapes. TCM predictions span $n_{ch}$ or centrality variation within a given
A-B collision system as well as variations from system to system (e.g. from
$p$-$p$ to $p$-Pb to Pb-Pb). In contrast, the BW model is not expected to
provide such absolute predictions: “…the normalization of the spectrum…we will
always adjust for a best fit to the data, because we are only interested in
the shape of the spectra to reveal the dynamics of the collision zone at
freeze-out [emphasis added]” ssflow . Given those limitations any direct model
comparisons must be based on shape measures that are not model dependent so as
to provide unbiased results. Setting aside absolute yields how best can one
demonstrate the evolution of spectrum shape, thus how best test model
validity?

For what follows groups of spectra under comparison (e.g. $n_{ch}$ variation
for given hadron species) are rescaled so as to agree in amplitude near
$p_{t}$ = 0. The specific common value is not relevant. Successive differences
between spectra are considered in relation to TCM hard-component model shapes.
Integrals of high-$p_{t}$ and low-$p_{t}$ intervals (absolute yields) for
variable-TCM representations of data spectra are examined. Given that TCM
context, ratios of high-$p_{t}$ to low-$p_{t}$ intervals for both TCM and BW
model are considered as a function of $n_{ch}$ or centrality. Logarithmic
derivatives (curvature measures) are applied to TCM and BW models to provide
differential spectrum shape information. Spectrum shape variations for each BW
model parameter (with others held fixed) are examined in turn. Finally, the BW
model is compared with a Boltzmann exponential and with TCM
$\hat{S}_{0}(m_{t})$.

###  IV.1 Spectrum differences

Figure 3 shows data spectra (solid) for pions (left) and protons (right). The
$p_{t}$ acceptances are rather limited and thus not pursued further in this
section. The TCM (dashed) describes data within their statistical
uncertainties as demonstrated in Ref. pidpart2  Sec. III E. TCM soft
components $\hat{S}_{0}(y_{t})$ (dotted) represent the limiting case for
spectrum data with charge density $\bar{\rho}_{0}\rightarrow 0$ (i.e. zero
particle density). For pions (left) a second $\hat{S}_{0}(y_{t})$ curve below
the data at low $y_{t}$ includes no correction for conjectured resonance
contributions (see Sec. III A of Ref. pidpart1 ). Arrows in these and other
panels indicate $y_{t}$ intervals used for BW fits as reported in Ref.
aliceppbpid .

Figure 3:  Data spectra (densities on $p_{t}$) for identified pions and
protons from 5 TeV $p$-Pb collisions aliceppbpid  (solid curves) and
corresponding TCM curves (dashed) for four centrality classes plotted on
logarithmic transverse rapidity $y_{t}$. The spectra are rescaled to quantity
$X_{i}(y_{t})$ as defined in Eq. (II.1). For clarity only curves for
centralities $n=$ 1, 3, 5, 7 are plotted. The lower dotted
$\hat{S}_{0}(y_{t})$ curve for pions is without resonances.

Figure 4 shows spectrum data for neutral kaons (a) and Lambdas (c) with the
same line styles as described above. Again the TCM describes data within their
statistical uncertainties. The $p_{t}$ acceptances are large in both cases,
providing meaningful model tests below.

Figure 4 (c,d) show differences between spectra for a given $n_{ch}$ class
($n=$ 1-6) and spectra for the lowest $n_{ch}$ class ($n=7$) for each of
neutral kaons (b) and Lambdas (d), and for data (solid) and TCM (dashed). The
point of this exercise is to isolate spectrum “hard components” without
resorting to a priori physical assumptions, thereby repeating the empirical
procedure employed in Ref. ppprd  first used to establish the TCM for 200 GeV
$p$-$p$ collisions. As noted in that $p$-$p$ analysis data structures so
isolated are approximately independent of $n_{ch}$ (e.g. $p$-Pb centrality),
thus demonstrating that full spectra can be represented as superpositions of
two fixed model functions. The amplitude of one (hard) varies as the square of
the other (soft) for each N-N collision.

(a)(b)
(c)(d)

Figure 4:  Data $p_{t}$ spectra for identified neutral kaons (a,b) and Lambdas
(c,d) from 5 TeV $p$-Pb collisions aliceppbpid  (solid curves) and
corresponding TCM curves (dashed) for four centrality classes. Panels (a,c)
compare rescaled data spectra (solid) to TCM (dashed). Panels (b,d) show
differences between spectra for $n=1-6$ and peripheral spectrum $n=7$.

More specifically, modes of the “hard components” move from lower to higher
$y_{t}$ with increasing hadron mass such that spectra for various hadron
species tend to coincide at high $y_{t}$ (e.g. near $y_{t}$ = 5 or
$p_{t}\approx 10$ GeV/c). Baryon hard components also tend to fall more
rapidly below the mode than mesons. Those trends are consistent with
fragmentation functions for identified hadrons from $e^{+}$-$e^{-}$ collisions
eeprd . Systematic uncertainties below the mode for difference spectra in
panels (b,d) are substantially greater than reported in Refs. pidpart1 ;
pidpart2  because of the simpler method invoked here to reduce model
dependence. Similar results are obtained by subtracting $\hat{S}_{0}(y_{t})$
model functions (dotted) in the left panels.

Figure 5 shows pion (left) and proton (right) spectra from the BW model in Eq.
(5) using fitted parameter values in Table 3, both as reported in Ref.
aliceppbpid . Soft-component model $\hat{S}_{0}(y_{t})$ (dashed) at left does
not include the conjectured resonance contribution that appears in Fig. 3
(left, upper dotted) and describes pion data. The BW curves (solid) have been
rescaled to coincide near $y_{t}$ = 0 with the $\hat{S}_{0}(y_{t})$ soft-
component curves (dashed).

Figure 5:  Blast-wave model $m_{t}$ spectra for identified pions and protons
from 5 TeV $p$-Pb collisions aliceppbpid  for seven centrality classes. The
bold dashed curves are TCM soft-component model $\hat{S}_{0}(y_{t})$. The
curves are scaled to coincide at $y_{t}$ = 0.

Figure 6 (a,c) shows neutral kaon (a) and Lambda (c) spectra from the BW model
in Eq. (5) using fitted parameter values in Table 3, both as reported in Ref.
aliceppbpid .

(a)(b)
(c)(d)

Figure 6:  Blast-wave model $m_{t}$ spectra (solid) for identified neutral
kaons (a,b) and Lambdas (c,d) from 5 TeV $p$-Pb collisions aliceppbpid . The
bold dashed curves in (a,c) are TCM soft-component model $\hat{S}_{0}(y_{t})$.
Curves in (a,c) are scaled to coincide at $y_{t}$ = 0. Panels (b,d) show
differences (solid) between BW spectra for $n=1-6$ and peripheral spectrum
$n=7$. Dashed curves in (b,d) are TCM dashed curves from Fig. 4 (b,d) that
represent the data spectrum differences there. The TCM curves are all rescaled
by the same factor to best correspond with BW curves within fit intervals
(arrows).

Figure 6 (b,d) shows differences between BW spectra for a given $n_{ch}$ class
($n=$ 1-6) and spectra for the lowest $n_{ch}$ class ($n=7$) for each of
neutral kaons (b) and Lambdas (d). The dashed curves are TCM differences
repeated from panels (b,d) of Fig. 4 to provide a reference. In either case
the dashed curves have all been rescaled by the same factor to accommodate the
solid BW curves.

This comparison highlights two issues: (a) The BW curves below the mode remain
high while the TCM data representations fall more rapidly. That is especially
obvious for the Lambda spectra since baryon hard components have modes at
higher $y_{t}$ and are narrower compared to meson hard components pidpart2 .
(b) The BW curves above the mode fall off rapidly while the TCM data
representations follow a power-law trend that is more gradual. That is
especially obvious for the kaon spectra for the same reasons given above. The
contrast is amplified by the different BW fit intervals used for the two
cases.

###  IV.2 High/low ratio centrality trends

Certain centrality trends are more interpretable when presented in the context
of 5 TeV $p$-Pb centrality variation as reported in Ref. tomglauber  based on
accurate TCM description of ensemble $\bar{p}_{t}$ data for that collision
system. That centrality model then provided the basis for PID spectrum studies
reported in Refs. ppbpid , pidpart1  and pidpart2 ,

Figure 7 [Fig. 17 of Ref. pidpart1  in relation to its Eqs. (26) and (27)]
shows soft fraction $\bar{\rho}_{s}/\bar{\rho}_{0}$ (left) and hard/soft ratio
$\bar{\rho}_{h}/\bar{\rho}_{s}\equiv x(n_{s})\nu(n_{s})$ (right) vs charge-
density soft component $\bar{\rho}_{s}$, with $\bar{\rho}_{x}\equiv
n_{x}/\Delta\eta$ and $\bar{\rho}_{0}=\bar{\rho}_{s}+\bar{\rho}_{h}$. The
straight line at right represents the relation
$\bar{\rho}_{h}\approx\alpha\bar{\rho}_{s}^{2}$ (as first reported in Ref.
ppprd ) with $\alpha\approx 0.013$ for 5 TeV N-N collisions alicetomspec ;
pidpart1 . It is notable that (a) the lowest four of seven $p$-Pb $n_{ch}$
classes are effectively equivalent to single peripheral $p$-N collisions and
(b) for higher $n_{ch}$ classes $p$-Pb centrality does increase significantly
but jet production (measured by $\bar{\rho}_{h}/\bar{\rho}_{s}$) increases
less rapidly.

Figure 7:  Soft fraction $\bar{\rho}_{s}/\bar{\rho}_{0}$ (left) and hard/soft
ratio $\bar{\rho}_{h}/\bar{\rho}_{s}\equiv x(n_{s})\nu(n_{s})$ (right) vs
charge-density soft component $\bar{\rho}_{s}$. The straight line at right
represents the relation $\bar{\rho}_{h}\approx\alpha\bar{\rho}_{s}^{2}$ (with
$\alpha=0.013$) that is observed for single N-N collisions ppprd .

The comparison of TCM and BW spectrum variation with centrality in the
previous subsection, while model-independent, is not quantitative and is
therefore not conclusive. One can alternatively compare yields at low and high
$p_{t}$ for data (via variable TCM) and the BW model.

Figure 8 compares yields within $p_{t}$ intervals for $p_{t}<0.15$ GeV/c (low)
to those for $p_{t}>5$ GeV/c (high). The solid points are derived from data
extrapolations defined by the variable TCM (which represents PID spectrum data
accurately as established in Ref. pidpart2 ). A reference high-vs-low trend
may be derived from the fixed PID TCM from Ref. pidpart1  as follows. Since
pidpart1

$\displaystyle\bar{\rho}_{si}$ $\displaystyle\approx$ $\displaystyle
z_{si}(n_{s})\bar{\rho}_{s}~{}~{}\text{and}$ (6)
$\displaystyle\bar{\rho}_{hi}$ $\displaystyle\approx$ $\displaystyle
z_{hi}(n_{s})\bar{\rho}_{h}\approx\tilde{z}_{i}x(n_{s})\nu(n_{s})\bar{\rho}_{si}$

for hadron species $i$ then
$\text{high}\approx\tilde{z}_{i}x(n_{s})\nu(n_{s})\times\text{low}$ assuming
fixed TCM model functions (solid curves). Values for $x(n_{s})$ and
$\nu(n_{s})$ are taken from Table I and for ratios $\tilde{z}_{i}$ are derived
from Tables IV and V of Ref. pidpart1 . Note that for those tables $z_{si}$
are determined near 0.15 GeV/c ($y_{t}\approx 1$) while $z_{hi}$ are obtained
near 1 GeV/c ($y_{t}\approx 2.7$).

(a)(b)
(c)(d)

Figure 8:  High vs low spectrum-integral trends (solid dots) for PID data
spectra from 5 TeV $p$-Pb collisions and for (a) pions, (b) neutral kaons, (c)
corrected protons, (d) Lambdas. The $y_{t}$ integrals have the form
$(1/2\pi)dn_{ch}/dy_{z}$. The dashed lines represent high $\propto$ low2
reference trends (expected for N-N collisions). The solid triangles are
$x\nu\times$low trends. The solid curves are $\tilde{z}x\nu\times$low trends.
Refer to text for further details.

The solid dots represent accurate variable-TCM representations of data spectra
integrated over stated $p_{t}$ intervals. The dashed lines are high $\approx$
low2 $\rightarrow$ $x(n_{s})\times$low (reference trends) that would be
expected for single peripheral $p$-N collisions and fixed TCM model functions.
The triangles are $x(n_{s})\nu(n_{s})\times$low scaled to agree with the $p$-N
(dashed) trend for peripheral collisions. The relation of triangles to dashed
lines for pions may be compare with Fig. 7 (right). The solid curves are
$\tilde{z}_{i}x(n_{s})\nu(n_{s})\times$low scaled to agree with dashed trends
for peripheral collisions again assuming fixed TCM model functions.

Deviations of data trends (solid dots) from
$\tilde{z}_{i}x(n_{s})\nu(n_{s})\times$low trends (solid curves) indicate the
effect of hard-component model variations with $n_{ch}$ as reported in Sec.
III of Ref. pidpart2 . Meson trends fall well below fixed-TCM trends because
meson hard components above the mode are transported to lower $y_{t}$ with
increasing $p$-Pb $n_{ch}$. While baryon hard components also vary
substantially with increasing $n_{ch}$ (modes transported to higher $y_{t}$)
the variation is mainly below the mode. Baryon solid-dot and solid-curve high
vs low trends in panels (c,d) are thus different but much closer.

Those result rely on the TCM to provide accurate descriptions of absolute data
yields and spectra. As typically utilized, the BW model is not expected to
predict absolute yields: “…the normalization of the spectrum…we will always
adjust for a best fit to the data…” ssflow . However, spectrum high/low ratios
can be compared between TCM and BW models.

Figure 9 (left) shows high/low (hi/lo) ratios derived from integrals of BW
model spectra integrated over $p_{t}<0.15$ GeV/c (low) and $p_{t}>4$ GeV/c
(high) for four hadron species. The lower limit for “high” is reduced compared
to the 5 GeV/c limit in Fig. 8 to accommodate the restricted $y_{t}$ range of
the BW model spectra.

Figure 9:  Left: High/low (hi/lo) ratios derived from integrals of BW model
spectra integrated over $p_{t}<0.15$ GeV/c (low) and $p_{t}>4$ GeV/c (high)
for four hadron species. Right: High/low ratios for data represented by the
variable TCM. The lower limit for the “high” integral is 4 GeV/c to match the
left panel. The curves represent the $\tilde{z}_{i}x(n_{s})\nu(n_{s})$ trend
expected for fixed spectrum hard components. Note the systematic difference
between mesons and baryons as found in Ref. pidpart2 .

Figure 9 (right) shows high/low ratios for data (points) represented by the
variable TCM. The lower limit for “high” is here 4 GeV/c to match the left
panel. The curves are equivalent to solid curves in Fig. 8 representing the
$\tilde{z}_{i}x(n_{s})\nu(n_{s})$ trend expected for fixed spectrum hard
components. The curves are scaled vertically to pass through data points
corresponding to event class $n=4$.

While there is some similarity in the general shapes of BW and data high/low
$n_{ch}$ trends, there are large differences in mean values which are due in
part to the very different $p_{t}$ intervals imposed on BW model fits. BW
ratios for baryons are a factor 3-4 greater than those for data/TCM. The BW
fit interval for baryons is confined to higher $p_{t}$ and extends well beyond
hard-component modes. In contrast, the BW ratio trend for neutral kaons falls
significantly below that for data/TCM, and the fit interval extends on $p_{t}$
from just above the hard-component mode down to zero. The BW ratios for pions
are comparable to those for data/TCM, but the $p_{t}$ fit interval is
negligible in comparison to the other cases.

###  IV.3 Differential spectrum curvatures

While comparisons of high/low ratios as above do provide some quantitative
information on spectrum structure they are still not definitive. A more-
differential approach is called for. One possibility is logarithmic
derivatives. In Ref. pidpart2  Sec. IV E a logarithmic derivative was applied
to investigate conjectures in Ref. aliceppbpid  concerning power-law
dependences of integrated $p_{t}$ intervals on charge density
$\bar{\rho}_{0}$. It was demonstrated there that the approximate power-law
trends on $\bar{\rho}_{0}$ of baryon/meson spectrum ratios resulted from
transport of peaked spectrum hard components (for baryons) to higher $y_{t}$
with increasing $n_{ch}$, the trends corresponding quantitatively to the
Gaussian (+ exponential tail) shapes on $y_{t}$ of spectrum hard components.
In the present context direct analysis of $y_{t}$ spectra via logarithmic
derivative provides the required information.

The logarithmic derivative applied to spectra $\bar{\rho}_{0}(y_{t})$ can be
illustrated using the TCM hard-component model [i.e.
$\bar{\rho}_{0}(y_{t})\rightarrow\hat{H}_{0}(y_{t})$, see Ref. pidpart2  Eq.
(8)]:

$\displaystyle-\frac{d\ln[\bar{\rho}_{0}(y_{t})]}{dy_{t}}$
$\displaystyle\rightarrow$
$\displaystyle\frac{y_{t}-\bar{y}_{t}}{\sigma^{2}_{y_{t}}}~{}~{}\text{near the
mode}$ $\displaystyle\approx$ $\displaystyle q~{}~{}\text{well above the
mode}.$

The second derivative leads to

$\displaystyle-\frac{d^{2}\ln[\bar{\rho}_{0}(y_{t})]}{dy_{t}^{2}}$
$\displaystyle\rightarrow$ $\displaystyle 1/\sigma_{y_{t}}^{2}~{}~{}\text{near
the mode}$ $\displaystyle\approx$ $\displaystyle 0~{}~{}\text{well above the
mode}$

The second derivative functions as a measure of local curvature or rate of
change of local slope of a spectrum.

$p$-Pb spectra can be characterized in general as follows (based on TCM
structure): At low $p_{t}$ spectra vary as $\hat{S}_{0}(m_{t})$ –
approximately as a Boltzmann exponential on $m_{ti}$ with $m_{ti}\rightarrow
m_{i}\cosh(y_{ti})$ for hadron species $i$. At high $p_{t}$ spectra
approximate a power law on $p_{t}$ $\sim p_{t}^{n}$ (with $p_{t}\sim
m_{0}\sinh(y_{t})$) or exponential on $y_{t}$. Based on those characteristics
one expects for Eq. (IV.3)

$\displaystyle-\frac{d^{2}\ln[\bar{\rho}_{0}(y_{t})]}{dy_{t}^{2}}$
$\displaystyle\sim$ $\displaystyle\cosh(y_{t})~{}~{}\text{at low $p_{t}$}$
$\displaystyle\sim$ $\displaystyle 0~{}~{}~{}\text{at high $p_{t}$}.$

Note that $\hat{S}_{0}(y_{t})$ (a Lévy distribution on $m_{t}$) itself goes to
a power law on $m_{t}$ at higher $y_{t}$. Its logarithmic second derivative
should thus also fall to zero at higher $y_{t}$.

Figure 10 (left) shows Eq. (IV.3) (left side) applied to neutral kaon spectra
from 5 TeV $p$-Pb collisions as represented by the variable TCM. The kaon data
are preferred because of the large $y_{t}$ acceptance. The line styles
proceeding down from most-central are solid, dashed, dotted and dash-dotted
followed by solid for the remainder. The trends at low $y_{t}$ correspond to
$\cosh(y_{t})$ and at high $y_{t}$ drop to zero as anticipated by Eq. (IV.3).
The hatched band corresponds to $1/\sigma^{2}_{y_{t}}$ for neutral kaons (see
Table 2). The trend expected for data hard components alone is constant
(within the hatched band) terminating near $y_{t}$ = 4 (where the hard
component transitions from Gaussian to exponential) dropping toward zero above
that point. Variation between $y_{t}$ = 2 and 4 corresponds to relative
contributions of soft and hard components at a given $y_{t}$ varying with
$p$-Pb centrality ($\sim$ product $x\nu$).

Figure 10:  Left: Logarithmic second derivative Eq. (IV.3) applied to variable
TCM spectra accurately representing data for neutral kaons. Line styles
proceed down from most-central as solid, dashed, dotted, dash-dotted. The
hatched band is limiting case $1/\sigma_{y_{t}}^{2}$ for Eq. (IV.3). Right:
Logarithmic second derivative Eq. (IV.3) applied to BW spectra with the same
line-style sequence. The bold solid curve is the $n=1$ most-central curve from
the left panel.

Figure 10 (right) shows Eq. (IV.3) (left side) applied to BW spectra (curves
of several line styles) as reported in Ref. aliceppbpid  and as in Fig. 6 (a).
The result for the variable TCM representing most-central ($n=1$) $p$-Pb data
(left panel) is included as a reference. The strong deviation between BW and
data/TCM is evident above $y_{t}$ = 2.5 ($p_{t}\approx 0.8$ GeV/c). Variation
of BW model parameters with $n_{ch}$ leads to rather small spectrum shape
variations (compared to the left panel) below $y_{t}$ = 3 where $\chi^{2}$
fits should matter. What follows is a procedure to examine the effect of
varying each of three BW parameters in turn.

###  IV.4 BW model response to individual parameters

Figure 11 (a) shows the effect of varying “kinetic freezeout” temperature
$T_{kin}$ over the values reported in Ref. aliceppbpid  while holding the
other two parameters fixed at their values for middle centrality $n=4$. The
result is relatively little variation of BW spectrum shapes.

Figure 11 (b) shows the effect of varying mean radial speed $\bar{\beta}_{t}$
over its reported values while holding the other two parameters fixed. The
result is large variation of curvatures, spectra being shifted strongly to
higher $y_{t}$ with increasing centrality across the entire $y_{t}$ interval.

(a)(b)
(c)(d)

Figure 11:  BW second derivatives with only one parameter varying: (a)
$T_{kin}$ varies, (b) $\bar{\beta}_{t}$ varies, (c) profile parameter $n$
varies. In each case other BW parameters are maintained at centrality $n=4$
values and the TCM reference for $n=1$ is included. (d) BW $\beta_{t}(r)$
profiles for seven $p$-Pb centralities.

Figure 11 (c) shows the effect of varying flow-profile parameter $n$ over its
reported values (per Table 3) while holding the other two parameters fixed.
The result is moderate shift of spectra to lower $y_{t}$ with increasing
centrality, but only above $y_{t}$ = 2. Figure 11 (d) shows the effect of
varying parameter $n$ on the “flow profile” or variation of radial speed
$\beta_{t}(x)$ with fractional radius $x=r/R$.

Based on Fig. 11 (a) and Z-scores in the next section the BW model as applied
to neutral kaons is not relevant above $y_{t}\approx 2.8$ ($p_{t}\approx 1.15$
GeV/c). Below that point $T_{kin}$ does not significantly affect the model
compared to the other two BW parameters (acts mainly as an overall scale
factor that is seen as not relevant). $\bar{\beta}_{t}$ affects the model over
a broad interval but $n$ is controlling only above $y_{t}$ = 1.8
($p_{t}\approx 0.4$ GeV/c), and the two parameters produce strongly opposing
variations in spectrum shape. In effect, the BW model as applied to neutral
kaons is a two-parameter model function that attempts to describe data over a
limited $p_{t}$ interval. Based on the above results physical interpretation
of model parameters directly via comparisons with data trends seems
problematic.

Figure 12 (left) shows the effect of turning off radial flow
($\bar{\beta}_{t}=0$). The shape change corresponds to faster fall-off with
decreasing $T_{kin}$ value (and increasing $p$-Pb centrality). The bold dotted
curve is the Boltzmann (exponential on $m_{t}$) limit of $\hat{S}_{0}(y_{t})$
(i.e. if Lévy exponent $n\rightarrow\infty$) with $T=145$ MeV. BW model
spectra approach that limit for central $p$-Pb collisions and $T_{kin}\approx
145$ MeV.

Figure 12:  Left: BW model spectra with $\bar{\beta}_{t}=0$ (various line
styles). $\hat{S}_{0}(y_{t})$ (bold dashed) corresponds to $T_{i}=200$ MeV for
kaons. The reference Boltzmann distribution (bold dotted) corresponds to
$T=145$ MeV. Right: BW model spectra with fixed $T_{kin}=145$ MeV (= TCM $T$
value for unidentified hadrons) with $\bar{\beta}_{t}$ and $n$ varying as in
Table 3.

Figure 12 (right) shows the effect of holding temperature $T_{kin}$ fixed at
145 MeV while the other two parameters remain consistent with Table 3. The BW
model fitted to the most-peripheral $p$-Pb data (with $\bar{\beta}_{t}\approx
0.25$) corresponds closely to TCM model $\hat{S}_{0}(y_{t})$ (bold dashed).

Note that $\hat{S}_{0}(y_{t})$ for the TCM (bold dashed) in the left panel is
defined to match data spectra as $n_{ch}\rightarrow 0$. As such it has no jet
contribution and represents zero particle density. The BW model with no jet
component and zero flow contribution should have the best opportunity to
describe spectrum data there given BW model assumptions. But the BW model
strongly deviates from $\hat{S}_{0}(y_{t})$ above $y_{t}$ = 2 ($p_{t}\approx
0.5$ GeV/c). BW curves with $\bar{\beta}_{t}=0$ and $T_{kin}\approx 145$ MeV
(central $p$-Pb in Table 3) instead closely approximate a Boltzmann
exponential (bold dotted) on transverse mass $m_{t}$ with slope parameter
$T=145$ MeV scaled to coincide with $\hat{S}_{0}(y_{t})$ at low $y_{t}$.

The power-law tail on $\hat{S}_{0}(y_{t})$ (as a Lévy distribution levywilk )
can be interpreted to represent inhomogeneity of the emitting system, i.e. the
Lévy model parameter $n$ can be expressed in the form
$1/n\sim\sigma^{2}_{T}/\bar{T}^{2}$ wilklevy . The inferred width $\sigma_{T}$
may be interpreted in terms of incomplete equilibration or in terms of $k_{t}$
broadening within a longitudinal splitting cascade. The result in the right
panel suggests that a substantial part of the $\bar{\beta}_{t}$ contribution
simply accommodates the Lévy shape of soft component $\hat{S}_{0}(y_{t})$
common to all A-B collision systems and has nothing to do with radial flow per
se. That explains why low-$n_{ch}$ $p$-$p$ and peripheral A-B collision data
require $\bar{\beta}_{t}\approx 0.25$ alicepppid  for Boltzmann-based BW
models.

##  V BW $\bf vs$ TCM Spectrum fit quality

This section quantitatively evaluates spectrum fit quality for the BW model
compared to the TCM as applied to 5 TeV $p$-Pb collisions. Three goodness-of-
fit measures are compared: (a) Z-scores, (b) $\chi^{2}$ and (c) data/model
spectrum ratios. Fit quality based on statistical errors is compared with that
based on systematic errors, both as published with spectrum data in
association with Ref. aliceppbpid .

###  V.1 Goodness-of-fit measures

A standard measure of data description or fit quality for a model is the
Z-score (for an $i^{th}$ observation) zscore

$\displaystyle Z_{i}$ $\displaystyle=$
$\displaystyle\frac{O_{i}-M_{i}}{\sigma_{i}},$ (10)

where $O_{i}$ is an observation (datum), $M_{i}$ is a model prediction and
$\sigma_{i}$ is the uncertainty (“error”) for the corresponding observation.
In some presentations $M_{i}\rightarrow E_{i}$ (expectation) and
$\sigma_{i}\rightarrow\sqrt{E_{i}}$ (assuming Poisson statistics for $O_{i}$).
For an acceptable model one expects $Z_{i}$ to have an r.m.s. value near 1.
The $\chi^{2}$ statistic is simply related:

$\displaystyle\chi^{2}$ $\displaystyle=$
$\displaystyle\sum_{i=1}^{N}Z_{i}^{2}$ (11)

for $N$ data points. An acceptable model should yield $\chi^{2}\sim\nu=N-D$
where $D$ is the number of degrees of freedom of the model (e.g. number of
model parameters).

There is a tendency to evaluate model quality in terms of data/model ratios –
e.g. Fig. 7 of Ref. aliceppbpid . The relation between data/model ratios and
Z-scores is

$\displaystyle\frac{\text{data}}{\text{model}}-1$ $\displaystyle\approx$
$\displaystyle\text{Z-score}\times\frac{\text{error}}{\text{data}},$ (12)

with error/model (exact) $\rightarrow$ error/data (approximate). The
error/data ratio may vary by orders of magnitude between different particle
types and collision systems, and even across $y_{t}$ intervals (see Figs. 17
and 18). Whereas the correct test of model validity is the relation of
Z-scores to 1 (or $\chi^{2}$ to $\nu$) the interpretation of data/model ratios
relative to 1 is problematic and typically quite misleading. The error/data
ratios on the right are typically $\ll 1$ thereby suppressing, in data/model
ratios on the left, what may be very significant data-model deviations.

###  V.2 Fit quality based on statistical uncertainties

Figures 13 and 14 (a,c) compare BW model fits (dashed) to $p$-Pb spectrum data
(solid). The spectrum data are scaled as shown in the $y$-axis labels and as
reported in Refs. pidpart1 ; pidpart2 . The BW curves correspond to the fitted
parameters in Table 3 and are scaled to best match the data, consistent with
the fitting procedure.

(a)(b)
(c)(d)

Figure 13:  Data $p_{t}$ spectra (solid) for pions (a) and neutral kaons (c)
and corresponding BW model curves (dashed). Data spectra are rescaled to
quantity $X_{i}(y_{t})$ as defined in Eq. (II.1). BW curves are scaled to
match the data within fit intervals (arrows). Panels (b,d) show corresponding
Z-scores based on statistical errors published with the spectrum data.

Figures 13 and 14 (b,d) show Z-scores for BW model vs spectrum data using
statistical uncertainties that accompany the published data. It is evident
that even within $y_{t}$ intervals chosen to produce an acceptable fit
aliceppbpid  the Z-scores are consistently substantially greater than 1.
Results for BW fits can be compared with results from a TCM analysis of the
same data reported in Refs. pidpart1 ; pidpart2 .

(a)(b)
(c)(d)

Figure 14:  Data $p_{t}$ spectra for protons (a) and Lambdas (c) (solid) and
corresponding BW model curves (dashed). Data spectra are rescaled to quantity
$X_{i}(y_{t})$ as defined in Eq. (II.1). BW curves are scaled to match the
data within fit intervals (arrows). Panels (b,d) show corresponding Z-scores
based on statistical errors published with the spectrum data.

Figures 15 and 16 compare data/TCM ratios (left panels) to Z-scores (right
panels), for mesons (pions and kaons) in Fig. 15 and for baryons (protons and
Lambdas) in Fig. 16. Based solely on data/model ratios it would seem that the
TCM description for pion data is much better than that for Lambda data. Yet a
comparison of corresponding Z-scores shows that the data descriptions are
actually of comparable quality. While Z-scores within some $y_{t}$ intervals
are consistent with statistical deviations (i.e. r.m.s. consistent with 1)
prominent localized excursions are apparent – for instance $y_{t}$ near 2 for
pions and $y_{t}$ within 2-3 for protons. Given that the significant data-TCM
deviations are consistent across $p$-Pb $n_{ch}$ classes and highly localized
on $y_{t}$ it is unlikely that the deviations result from model issues. It is
reasonable to conjecture that such structures are inevitable within the
complex procedure of processing particle data. These results suggest that the
TCM describes $p$-Pb PID spectrum data within statistical uncertainties except
for local anomalies that may be due to detector effects.

(a)(b)
(c)(d)

Figure 15:  Left: Data/TCM spectrum ratios for pions (a) and neutral kaons
(c). Right: Z-scores for pions (b) and neutral kaons (d). The model is the
variable TCM defined in Ref. pidpart2 .

(a)(b)
(c)(d)

Figure 16:  Left: Data/TCM spectrum ratios for protons (a) and Lambdas (c).
Right: Z-scores for protons (b) and Lambdas (d). The model is the variable TCM
defined in Ref. pidpart2 .

###  V.3 Fit quality based on systematic uncertainties

The BW-fit $\chi^{2}/\text{ndf}$ values appearing in Table 3, as reported in
Ref. aliceppbpid , are notable for two reasons: they are (a) consistently low
compared to the expected r.m.s. value $\approx 1$ for acceptable fits and (b)
apparently at odds with the Z-scores shown in Figs. 13 and 14 that are
consistent with $\chi^{2}/\text{ndf}$ values values substantially greater than
10. Those results suggest that $\chi^{2}/\text{ndf}$ values reported in Ref.
aliceppbpid  were obtained assuming systematic rather than statistical
uncertainties $\sigma_{i}$ in Eq. (10).

Figures 17 and 18 (a,c) show BW-fit Z-scores recalculated using systematic
errors included with published $p$-Pb spectra from Ref. aliceppbpid . Vertical
scales are the same as those in Figs. 13 and 14 (a,c) for comparison. Those
results suggest that systematic uncertainties were in fact used to calculate
$\chi^{2}/\text{ndf}$ values appearing in Table 3.

(a)(b)
(c)(d)

Figure 17:  BW model Z-scores for pions (a) and neutral kaons (c) based on
systematic uncertainties or errors published with spectrum data. Panels (b,d)
show error/data ratios for statistical errors (solid) and systematic errors
(dashed).

Figures 17 and 18 (b,d) compare published $p$-Pb systematic uncertainties
(dashed) with statistical uncertainties (solid) in ratio to corresponding
spectrum data values. Systematic uncertainties are approximately constant on
$y_{t}$ with magnitudes 5-10% of data values. The numbers are consistent with
total systematic uncertainties presented in Tables 3 and 4 of Ref. aliceppbpid
. The systematic uncertainties tend to exceed statistical uncertainties by a
factor five or more, especially within fit intervals employed for BW fits,
implying that reported $\chi^{2}/\text{ndf}$ values are at least a factor 25
less than what would be obtained using statistical uncertainties. That
description responds to items (a) and (b) at the beginning of this subsection.

(a)(b)
(c)(d)

Figure 18:  BW model Z-scores for protons (a) and Lambdas (c) based on
systematic uncertainties or errors published with spectrum data. Panels (b,d)
show error/data ratios for statistical errors (solid) and systematic errors
(dashed).

Appropriate uncertainty estimates are of central importance for realistic
assessment of model quality. The most basic estimate (where data result from
counting discrete elements – e.g. particles) is the statistical uncertainty
corresponding to Poisson statistics. Corresponding Z-scores then show all
statistically significant data-model deviations no matter what their origin.
Various sources of systematic uncertainty, whether correlated with collision
event class ($n_{ch}$) or particle properties ($y_{t}$, hadron species) or
not, may be relevant to correct assessment of model quality but should be
considered only relative to Z-scores based on statistical uncertainties.

More generally the $\chi^{2}$ statistic is not very effective for model
testing. $\chi^{2}$ is intended to measure a significant difference between
model and data but contains much less information than the differential
Z-scores it represents as in Eq. (11). Z-scores based on statistical errors
clearly show not only what data-model differences are statistically
significant but also how such deviations are correlated, e.g. across $y_{t}$
bins and $n_{ch}$ classes. In comparing two models the one with the greater
$\chi^{2}$ may even be preferred because of the detailed structure of
Z-scores.

##  VI Blast wave model evolution

In order to interpret BW model fit results it is useful to review early
theoretical developments that led to the current general form applied to
ultrarelativistic A-B collisions and some significant variations in the
structure of the model over several decades.

###  VI.1 Early BW model development

The blast-wave model arguably has its origins in the collision theories of
Landau (1950s, e.g. Ref. landau ) and Hagedorn (1960s, e.g. Ref. haged ).
Hadron production is based in the first case on hydrodynamic evolution of a
hadronic fluid and in the second case on emission from rapidly-moving
“fireballs.” In either case “collectivity” – collective motion of a continuous
hadronic fluid or discrete collection of objects (fireballs) in configuration
space – is a central feature. Detected hadrons that emerge are thereby
correlated in momentum space, such correlations interpreted to indicate a form
of collectivity.

The correct relation between the measured particle momentum distribution in
the center-of-momentum (CM) or laboratory (lab, for symmetric colliders) frame
and the hadron emission distribution in the local, comoving or boost frame was
determined by Cooper and Frye cooper  as

$\displaystyle E\frac{d^{3}N}{dp^{3}}$ $\displaystyle=$
$\displaystyle\frac{d^{3}N}{m_{t}dm_{t}dy_{z}d\phi_{p}}=\int_{\sigma}f(x,p,u)p^{\mu}d^{3}\sigma_{\mu},$
(13)

where $u(x)$ is the velocity field of an emitting fluid in the lab frame,
$p^{\mu}$ is the particle four-momentum measured in the lab frame, $\sigma$ is
the “freezeout” hypersurface ($\sim$ transition from small to large mean free
path), and $d^{3}\sigma_{\mu}$ is a differential volume four-vector. A
limiting case for $f(x,p,u)$ in the local or boost frame is isotropic
emission, according to Bose or Fermi statistics or (as an approximation) the
Boltzmann distribution. That formalism was employed to describe spherical
expansion within nucleus-nucleus collisions in Ref. siemens  wherein the term
“blast wave” was introduced in connection with A-B collisions. More complex
expansion geometries (i.e. nucleon flows) were encountered within the Bevalac
program bevalac .

###  VI.2 BW model for ultrarelativistic A-B collisions

For the Bevalac and alternating gradient synchrotron (AGS) nucleus-nucleus
programs, where the collision energy per nucleon is of order the nucleon mass-
energy, fractional deviations from spherical expansion (e.g. squeeze-out,
side-splash) were measured via the sphericity tensor or related measures
sphericity ; daniel . For A-B collisions at the super proton synchrotron
(SPS), RHIC or LHC, collision energies are much greater than the nucleon mass-
energy. As a consequence the velocity field represented by $u(x)$ above is
highly asymmetric – strongly elongated along the collision axis. Reference
ssflow  reported a version of Eq. (13) adapted to that context and applied to
SPS S-S collisions at $\sqrt{s_{NN}}\approx 19$ GeV. It is useful to examine
details of the derivation and accompanying arguments.

It is assumed that within a local (co-moving or boost) frame hadron emission
is locally thermal and isotropic (i.e. within a transparent environment):
“…the invariant distribution function [in the CM or lab frame]…we assume to be
an isotropic thermal distribution boosted by the local fluid velocity
$u^{\mu}$, and we approximate the respective Bose and Fermi distributions by
the Boltzmann distribution.” Equation (13) is modified by

$\displaystyle f(x,p,u)$ $\displaystyle\rightarrow$
$\displaystyle[g/(2\pi)^{3}]\exp[-(u^{\nu}p_{\nu}-\mu)/T],$ (14)

and potential $\mu$ is omitted below (see comment below on “normalization”).
Configuration space is represented by cylindrical coordinates with a
hypersurface defined by

$\displaystyle\sigma_{\mu}$ $\displaystyle=$
$\displaystyle[t(\eta_{s}),r\hat{e}(\phi_{x}),z(\eta_{s})]$ $\displaystyle=$
$\displaystyle[\tau\cosh(\eta_{s}),r\hat{e}(\phi_{x}),\tau\sinh(\eta_{s})]$

with $\zeta\rightarrow\eta_{s}$ and

$\displaystyle d^{3}\sigma_{\mu}$ $\displaystyle=$
$\displaystyle\left[\frac{\partial z}{\partial\eta_{s}},0,0,\frac{\partial
t}{\partial\eta_{s}}\right]d\eta_{s}rdrd\phi_{x}$ $\displaystyle\rightarrow$
$\displaystyle\tau[\cosh(\eta_{s}),0,0,\sinh(\eta_{s})]d\eta_{s}rdrd\phi_{x}.$

Space-time rapidity is $\eta_{s}=(1/2)\ln[(t+z)/(t-z)]=\ln[(t+z)/\tau]$.
Quantities $u_{\mu}$ and $p_{\mu}$ are as defined in Eqs. (27) and (A.1.3).
Whereas Ref. ssflow  invokes “flow angles” $\rho$ and $\eta$, equivalent
quantities are presented here by $\eta_{t}$ and $\eta_{z}$ respectively. The
conventional assumption $\eta_{z}\rightarrow\eta_{s}$ is consistent with
Bjorken expansion.

Given those definitions Eq. (13) then becomes

$\displaystyle\frac{d^{3}N}{m_{t}dm_{t}dy_{z}d\phi_{p}}$ $\displaystyle=$
$\displaystyle\frac{g}{(2\pi)^{3}}m_{t}$
$\displaystyle\times\int_{-H_{s}}^{H_{s}}d\eta_{s}\left[\cosh(y_{z})\frac{\partial
z}{\partial\eta_{s}}-\sinh(y_{z})\frac{\partial t}{\partial\eta_{s}}\right]$
$\displaystyle\times\int_{0}^{R}rdr\exp\left\\{-m_{t}\cosh[\eta_{t}(r)]\cosh(y_{z}-\eta_{z})/T\right\\}$
$\displaystyle\times\int_{0}^{2\pi}d\phi_{x}\exp\\{p_{t}\sinh[\eta_{t}(r)]\cos(\phi_{x}-\phi_{p})/T\\},$

where $\eta_{s}\in[-H_{s},H_{s}]$ delimits the particle source. The quantity
in square brackets then simplifies to

$\displaystyle\cosh(y_{z})\frac{\partial
z}{\partial\eta_{s}}-\sinh(y_{z})\frac{\partial t}{\partial\eta_{s}}$
$\displaystyle\rightarrow$ $\displaystyle\tau\cosh(y_{z}-\eta_{s}).$ (18)

The integral over $\phi_{x}$ reduces to $2\pi I_{0}[p_{t}\sinh(\eta_{t})/T]$.
Using Eq. (18) and integrating both sides over $y_{z}$ leads to integrand
$K_{1}[m_{t}\cosh(\eta_{t})/T]$ in the third line. “We obtain the transverse
mass spectrum $dn/m_{T}dm_{T}$ by integrating over rapidity [$y_{z}$] using
the modified Bessel function $K_{1}$….” The integral over $\eta_{s}$ then
simply leads to factor $2Z_{t}\rightarrow 2H_{s}$ consistent with the interval
over space-time rapidity $\zeta\rightarrow\eta_{s}$ assumed in Eq. (14) of
Ref. ssflow . The result is Eq. (5) above except
$\rho\rightarrow\eta_{t}=\tanh^{-1}(\beta_{t})$, with
$\beta_{t}(r)=\beta_{s}(r/R)^{n}$. That expression is used by Ref. aliceppbpid
to obtain the BW model fits summarized in Sec. III.

As noted, the absolute scale of the BW model in relation to spectrum data is
not seen as significant (in contrast to the TCM). Reference ssflow  concludes:
“…we have shown that a thermal model is perfectly possible for S+S collisions
despite (or because of?) the similarity of S+S and pp spectra. The data force
us to include resonance decays and longitudinal flow while they make no
decisive statement about the existence of transverse flow.”

###  VI.3 Alternative BW model applications

It is informative to examine some alternative routes to BW model functions and
their applications to a variety of collision data. The examples below relate
to the experimental programs at the SPS, RHIC and LHC. Further details on BW
model derivations and applications are provided in App. A.

####  VI.3.1 Tomásik, Wiedemann and Heinz tomasik

The study in Ref. tomasik  combines BW spectrum analysis of single-particle
spectra and Bose-Einstein correlation (BEC) analysis of two-particle
correlations to infer space-time properties of the particle source (emission
region, freezeout state). The models (BW + BEC) were applied to NA49 data from
158A GeV Pb-Pb collisions. The particle source is modeled by an “emission
function” $S(x,p)d^{4}x$ that includes the factor $p^{\mu}d^{3}\sigma_{\mu}$
from Eq. (13) plus densities on $r$, $\eta_{s}$ and $\tau$ that explicitly
define the integrated space-time (source) volume.

Combining Eqs. (3.1) and (3.10) from Ref. tomasik

$\displaystyle E\frac{d^{3}N}{dp^{3}}$ $\displaystyle=$ $\displaystyle\int
S(x,p)d^{4}x$ $\displaystyle=$ $\displaystyle\frac{1}{(2\pi)^{3}}\int
d\eta_{s}rdrd\phi_{x}G(r)H(\eta_{s})T(\tau)$ $\displaystyle\times
m_{t}\,\tau\cosh(y_{z}-\eta_{s})\exp[-(u^{\nu}p_{\nu}-\mu)/T],$

where $H(\eta_{s})$ and $T(\tau)$ are Gaussians on $\eta_{s}$ and $\tau$
respectively, the latter normalized to unity.
$p_{\mu}=[m_{t}\cosh(y_{z}),p_{t},0,m_{t}\sinh(y_{z})]$. $G(r)$ is modeled as
a Gaussian or square (flat) distribution on $r$. In contrast, Ref. ssflow
implicitly assumes flat distributions on $r$ and $\eta_{s}$ and a delta
function on $\tau$. Velocity field $u(x)$ is initially defined in terms of
longitudinal $\eta_{l}\rightarrow\eta_{z}$ and transverse $\eta_{t}$
“expansion rapidities” (boosts) and $u(x)$ is then as in Eq. (27). $\eta_{l}$
is then identified with space-time rapidity $\eta\rightarrow\eta_{s}$ (i.e.
Bjorken expansion). In Ref. ssflow  transverse speed is parametrized by
$\beta_{t}(r)=\beta_{s}(r/R)^{n}$ but Ref. tomasik  expresses transverse boost
as $\eta_{t}(r)=\eta_{f}(r/r_{rms})$.

####  VI.3.2 Florkowski and Broniowski wojciechs

Reference wojciechs  briefly summarizes the history of the BW model and then
presents an alternative route to the model as it appears in Eq. (5). Symbols
used in Ref. wojciechs  are related to the standard set employed in this
study.

Space-time rapidity is $\alpha_{\parallel}\rightarrow\eta_{s}$ (what
correlates $t$ and $z$ given proper time $\tau$). $\zeta$ is a parameter that
correlates proper time $\tau$ and radius $\rho\rightarrow r$. That system is
assumed boost invariant over all $\eta_{s}$ in contrast to ssflow . Thus,
instead of $[t(\zeta),z(\zeta)]$ with $\zeta\rightarrow\eta_{s}$ space-time
rapidity there is $[\tau(\zeta),r(\zeta)]$ with $\zeta\in[0,1]$ simply a
correlation parameter. $\sigma_{\mu}(x)$ is formally as in the second line of
Eq. (VI.2), but the differential volume element $d^{3}\sigma_{\mu}$ (with
$\rho\rightarrow r$) is

$\displaystyle\left[\frac{dr}{d\zeta}\cosh(\eta_{s}),\frac{d\tau}{d\zeta}\hat{e}(\phi_{x}),\frac{dr}{d\zeta}\sinh(\eta_{s})\right]d\zeta\tau
rd\eta_{s}d\phi_{x}.$ (20)

$u(x)=\cosh(\eta_{t})[\cosh(\eta_{z}),\tanh(\eta_{t})\hat{e}(\phi_{x}),\sinh(\eta_{z})]$
is the velocity field (corrected from Ref. wojciechs ) with
$\alpha_{\perp}\rightarrow\eta_{t}$ and in this case
$\alpha_{\|}\rightarrow\eta_{z}$. Longitudinal flow is represented by
$v_{z}\rightarrow\beta_{z}=\tanh(\eta_{z})$, and transverse flow by
$\beta_{t}(\zeta)=\tanh[\eta_{t}(\zeta)]$. $\beta_{t}$ is correlated with
$\tau$ and $r$ by parameter $\zeta$ as alternative to explicit
$\beta_{t}(r)=\beta_{s}(r/R)^{n}$ in Ref. ssflow . $p_{\mu}$ and
$u^{\mu}p_{\mu}$ are as in Eqs. (A.1.3) and (A.1.4).

The transverse component $d\tau/d\zeta$ in Eq. (20) would introduce additional
complexity compared to Eq. (VI.2) and is assumed zero in Ref. wojciechs  “to
achieve the simplest possible form of the model.” $d^{3}\sigma_{\mu}$ then has
the form of Eq. (VI.2) (second line) except
$d\zeta\rho(\zeta)d\rho/d\zeta\rightarrow rdr$. Equation (28) of Ref.
wojciechs  corresponds to Eq. (VI.2) with $g\rightarrow 1$, but its Eq. (29)
assumes that $\alpha_{\perp}(\zeta)\rightarrow\eta_{t}(\zeta)$ is constant.

####  VI.3.3 Rath, Khuntia, Sahoo and Cleymans cleymans

The BW model invoked in Ref. cleymans  is nominally that of Ref. ssflow  but
differs substantially as follows. The differential volume element is given as

$\displaystyle d^{3}\sigma_{\mu}$ $\displaystyle=$
$\displaystyle\tau[\cosh(\eta_{s}),0,0,[-]\sinh(\eta_{s})]d\eta_{s}rdrd\phi_{x}$
(21)

but the bracketed minus sign is incorrect – compare Eq. (VI.2) (second line).
Quantities $u_{\mu}$ and $p_{\mu}$ are as defined in Eqs. (27) and (A.1.3).
Quantity $\eta\rightarrow\eta_{s}$ is space-time rapidity. It is notable that
“Bjorken correlation in rapidity” is interpreted as
$(y=\eta)\rightarrow(y_{z}=\eta_{s})$. In that case the expression in Eq. (18)
above goes to $\tau$ and the integration over $\eta_{s}$ in Eq. (VI.2) is
trivial. Function $K_{1}(x)$ could not arise based on that assumption. It
would be replaced by an exponential with the same arguments.

Explicitly, $\eta_{s}$ is space-time rapidity parameterizing $\sigma_{\mu}(x)$
and is usually assumed equivalent to longitudinal boost $\eta_{z}$ in velocity
field $u_{\mu}(x)$ (Bjorken expansion), both of those being flows in
configuration space. The assumption $y_{z}=\eta_{s}$ equates longitudinal
rapidity $y_{z}$ of detected particles in momentum space with space-time
rapidity or boost of a particle source in configuration space. Equation (14)
of Ref. cleymans  appears to be missing a combined factor $p_{t}m_{t}$ or
$m_{t}^{2}$ compared to Eq. (VI.2). The analysis is based on taking $n=1$ for
“flow profile” $\beta_{t}(r)\propto(r/R)^{n}$.

####  VI.3.4 Ray and Jentsch lanny

Reference lanny  uses the BW model as formulated in Ref. tomasik  to study the
effect of event-wise fluctuations within A-A collisions and their
manifestations in single-particle spectra on transverse rapidity
$y_{t}=\ln[(p_{t}+m_{t})/m_{0}]$ and on two-particle correlation space
$(y_{t1},y_{t2})$. The single-particle distribution in Ref. lanny  is derived
from Ref. tomasik  [see Eq. (VI.3.1) above] but it is assumed that
$y\rightarrow y_{z}\approx 0$. That is a reasonable assumption for data near
midrapidity and for the higher energies at the RHIC and LHC in contrast to
Ref. ssflow  relying on integration over $y_{z}$ to obtain Eq. (5) for
description of SPS S-S data at $\sqrt{s_{NN}}\approx 19$ GeV. The expression
on the right in Eq. (18) then becomes $\tau\cosh(\eta_{s})$ (with
$\tau_{0}\rightarrow\tau$). In Eq. (VI.2)
$\cosh(y_{z}-\eta_{z})\rightarrow\cosh(\eta_{z})$ in the exponential.

####  VI.3.5 Lao, Liu and Ma lao

Reference lao  invokes the BW model as reported in Eq. (18) of Ref. lanny  to
produce its Eq. (1). A “blast-wave model with fluctuations” is attributed to
Ref. tomasik  but the word “fluctuations” does not appear in that document.
Reference lanny  does describe a BW model with fluctuations but their Eq.
(18), adopted from Ref. tomasik  and the basis for Eq. (1) of Ref. lao ,
appears before the subject of fluctuations (in $\beta=1/T$ and in radial boost
$\eta_{t}$) is addressed. In Ref. lao  $y\rightarrow y_{z}$ remains variable
whereas $\cosh(y_{z}-\eta_{z})\rightarrow\cosh(y_{z})\cosh(\eta_{z})$ is
assumed. The derivation then proceeds to “single source emission” wherein
$\eta_{s}\rightarrow 0$ and source distribution $H(\eta_{s})\rightarrow 1$.
The particle source is thereby modeled effectively as a thin disk at
$\eta_{s}=0$. As a result $K_{1}[m_{t}\cosh(\eta_{t})/T]$ in Eq. (5)
$\rightarrow\exp[-m_{t}\cosh(y_{z})\cosh(\eta_{t})/T]$ in Eq. (4) of Ref. lao
. The $\beta_{t}(r)$ flow-profile parameter $n_{0}\rightarrow n$ is held fixed
at 2. It is further remarked that “$n_{0}$ is not a sensitive quantity. It
does not matter if $n_{0}=1$ or $n_{0}=2$.” But refer to Fig. 11 (c) for the
strong effect of parameter $n$.

Equation (4) of lao  is then applied to a broad array of LHC A-B collision
data. The study notes “…one can see that we have used the species-dependent
parameters to fit the spectra… [i.e. BW fits to individual hadron species].
This is not the usual way to use the blast-wave model, which fits various
spectra simultaneously. We have examined the simultaneous fit of the model for
various spectra and know that narrow and different $p_{T}$ ranges have to be
used for different particles [e.g. Ref. aliceppbpid ]. We do not think that
the simultaneous fitting of various spectra can be better in wide $p_{T}$
ranges. Instead, we may use the individual fit for different spectra and
obtain better fits.”

Even with fits to individual hadron species $\chi^{2}/\text{ndf}$ numbers are
often much greater than 1. Reference lao  comments that large
$\chi^{2}/\text{ndf}$ values “…indicates that the fitting was qualitative and
approximately acceptable, and the large dispersion between the curve and data
exists.” One observation is notable in the context of the present study:
“However, the disadvantage of the [BW] model is also obvious. In some cases,
one[-]component model cannot fit the data well. In fact, it is only applicable
to the low-$p_{T}$ region, but not to the high-$p_{T}$ region. This problem is
not only a disadvantage in the blast-wave model…but also a disadvantage in all
thermal models.”

These examples serve to illustrate several BW issues: Assumptions vary
concerning the structure of the flow field and emission volume and how they
should be represented. Concerning derivation of a specific BW model function
the basic Cooper-Frye context is accepted. However, significant details of
derivations vary. Different symbol choices by authors for basic kinematic
quantities complicate interpretation. There is also substantial variation
concerning BW model application to data: e.g. a common fit to all hadron
species as opposed to individual fits to species. There is the question what
$p_{t}$ fit intervals should be used, as determined by what criteria? It is
not clear how fits to data might be used to clarify such ambiguities. As a
matter of interpretation is there any circumstance in the comparison of model
to data wherein the BW model (which version?) might be falsified?

####  VI.3.6 Tsallis-BW models tsallisbw

A basic assumption of the BW models above is that any deviation from a
Boltzmann exponential (say below $m_{t}\approx 3$ GeV/c) must arise from a
moving particle source as in fluid flow ssflow . That assumption fails for two
reasons: (a) The Boltzmann exponential is an idealization quite unlikely in
real situations, especially in relation to parton splitting cascades resulting
in heterogeneous systems (parton showers) as hadron sources. (b) The strong
contribution from minimum-bias jet fragments (from large-angle-scattered
partons) extends down to 0.5 GeV/c with a peak near 1 GeV/c. That system is
described accurately by the TCM as demonstrated for 5 TeV $p$-Pb collisions in
Secs. II and V.2 and Refs. ppbpid ; pidpart1 ; pidpart2 .

Incorporation of the Tsallis spectrum model ($q$-exponential, Tsallis
statistics) as a replacement for the Boltzmann exponential within the BW model
formulation tsallisbw  is, in effect, an attempt to deal with (a) but relies
on variation of Tsallis model parameters with A-B centrality which is never
observed for soft component $\hat{S}_{0}(m_{t})$ of the TCM. As a monolithic
spectrum model the Tsallis model has, in principle, no capacity to accommodate
item (b) above. But model parameters $T$ and $q$ are varied as an attempt to
accommodate strong increase of the jet contribution with centrality. Model
fits in Ref. tsallisbw  are restricted to the range $p_{t}\in[0.5,3]$ GeV/c.

The main difference from the Boltzmann BW versions is that $\bar{\beta}_{t}$
for $p$-$p$ and peripheral A-B collisions is now consistent with zero as one
might expect, since the Tsallis model is a variant of soft-component model
$\hat{S}_{0}(m_{t})$ that describes all A-B soft components within statistical
uncertainties. Ironically, the decrease with increasing centrality of Tsallis
parameter $q-1$ (corresponding to increase of TCM parameter $n$ indicating
softening of the $\hat{S}_{0}(m_{t})$ power-law tail) may actually respond in
part to “jet quenching.” Although it cautions that “…the physical
interpretation of this [Tsallis] statistical model in the context of high
energy nuclear collisions remains to be fully understood” Ref. tsallisbw
interprets the $q-1$ trend as “…indicating an evolution from a highly non-
equilibrated system in p+p collisions toward an almost thermalized system in
central Au+Au.” There is no acknowledgement of a jet contribution to spectra
assumed (erroneously) as confined to “high-$p_{T}$ particles.” The possibility
that parameter $\bar{\beta}_{t}$ and $T_{kin}$ values, as inferred from a BW
model effectively falsified by data, may not be physically interpretable is
also not considered.

##  VII Systematic uncertainties

Systematic uncertainties for BW model fits to spectra relate to the primary
PID spectra, various BW model versions and fit methods and the validity of any
BW model version for part or all of the spectrum $p_{t}$ acceptance.

###  VII.1 Systematic uncertainties from Ref. aliceppbpid

Regarding PID spectra themselves Ref. aliceppbpid  states that “The main
sources of systematic uncertainties for the analysis of charged and neutral
particles are summarized in Tables 3 and 4, respectively. The study of
systematic uncertainties was repeated for the different multiplicity bins in
order to separate the sources of uncertainty which are dependent on
multiplicity and uncorrelated across different bins (depicted as shaded boxes
in the figures).” From that and other language in the text there are
apparently three types of spectrum uncertainties reported: statistical,
$n_{ch}$-dependent and total denoted by solid bars, shaded boxes and unshaded
boxes respectively. Captions for Figs. 2, 3, 4 and 5 include “The empty boxes
show the total systematic uncertainty; the shaded boxes indicate the
contribution uncorrelated across multiplicity bins….”

Spectrum data reported in Ref. aliceppbpid  apparently include statistical and
total uncertainties, the latter consistent with “Total” in Tables 3 and 4 of
Ref. aliceppbpid . There is no distinction in those tables among systematic
errors that are $y_{t}$ dependent or not and $n_{ch}$ dependent (correlated vs
uncorrelated) or not, and there is no estimation of what might contribute to
point-to-point (on $y_{t}$) uncertainties. See Figs. 15 and 16 (right) for
statistically-significant systematic errors that are correlated across
multiplicity $n_{ch}$ bins but uncorrelated (i.e. vary strongly) across
$y_{t}$ bins.

Properly-estimated systematic uncertainties relevant to model descriptions
would ideally be consistent with actual data-model discrepancies in the event
of a valid model. The results in Figs. 15 and 16 indicate that actual data-
model systematic errors are rather small, highly-localized on $y_{t}$ and
closely-correlated on $n_{ch}$. As noted above, TCM data descriptions are not
obtained by fits to individual spectra. The Z-scores in those figures suggest
that the systematic errors presented in Tables 3 and 4 of Ref. aliceppbpid
greatly overestimate the uncertainties relevant to data description by an
appropriate model.

###  VII.2 Sensitivity to BW fit ranges

A particular issue for BW model fits is the $p_{t}$ or $y_{t}$ intervals over
which fits to data are imposed. Fit ranges are determined by “the available
data at low $p_{T}$” (i.e. the effective particle acceptance lower cutoff) and
“the agreement with the data at high $p_{t}$, justified considering that the
assumptions underlying the blast-wave model are not expected to be valid at
high $p_{T}$” [emphasis added]. That policy begs the question where, if
anywhere, is the BW model valid and by what criteria? Parameter systematic
uncertainties are assigned empirically based on range variations: “The [fit]
results are reported in Tab. 5 and Fig. 6. Variations of the fit range[s] lead
to large shifts ($\sim$ 10%) of the fit results (correlated across
centralities)….” What criteria are invoked to determine such variations?

Concerning utility of BW fits there is also the statement “It has be be kept
in mind, however, that the actual values of the fit parameters depend
substantially on the fit range. In spite of this limitations [sic], the blast-
wave model still provides a handy way to compare the transverse momentum
distributions and their evolution in different collision systems.” The
statement apparently refers to BW fits to individual hadron species over the
entire $y_{t}$ acceptance as they appear in Fig. 1 of Ref. aliceppbpid . The
quality of those data descriptions is not reported, and the parameter values
are presumably meaningless based on the argument for limiting fit ranges above
(“…not expected to be valid at high $p_{T}$”). As to convenient data
representation the variable-TCM description in Fig. 1 above describes all PID
spectrum data within their statistical uncertainties without relying on fits
to individual spectra and based on physically meaningful parameters.

###  VII.3 Systematic vs statistical errors and Z-scores

As noted in Sec. V, Z-scores provide a superior differential measure of fit
quality compared to the $\chi^{2}$ statistic. Spectrum data associated with
Ref. aliceppbpid  include both statistical and systematic uncertainties. As
noted in Sec. V.3 $\chi^{2}$ values presented in Table 5 of Ref. aliceppbpid
appear based on total systematic uncertainties rather than statistical
uncertainties. That is one reason to question interpretation of such
$\chi^{2}$ values. Another is the practice to limit model-fit $p_{t}$
intervals to what will provide a “good” fit, perhaps as determined by
resulting $\chi^{2}$ values.

An example of Z-score usage as standard practice for particle physics analysis
can be found in Ref. alephff  reporting properties of hadronic Z decays from
the LEP. In its Fig. 1 appears “the difference between the distributions of
the QCD models and the [$e^{+}$-$e^{-}$] data in units of the data error,”
which is just the Z-score defined in Eq. (10). The “data error” here is
“combined statistical and systematic errors.” Although that seems to
correspond formally to $\chi^{2}$ values in Ref. aliceppbpid  the two
uncertainty types have approximately the same magnitudes for the
$e^{+}$-$e^{-}$ spectrum data in Fig. 20 below, in contrast to the factor
$\approx$ 10 difference for spectrum data reported in Ref. aliceppbpid .

Part of systematic uncertainty estimation relates to the interpretability of
analysis results, but applies in this case to model uncertainty estimations
themselves as opposed to data. If a statistic is employed to determine data-
model agreement/disagreement do its resulting values provide a reliable
evaluation of model validity? In the case of $\chi^{2}$ values in Table 5 of
Ref. aliceppbpid  it is likely that the values reported are misleading, hence
uncertain.

It can be argued that Z-scores should be calculated only with statistical
uncertainties, which are usually unambiguous in a context where discrete
objects are counted. Z-scores based solely on statistical uncertainties
directly reveal systematic errors (as opposed to uncertainties) in both data
and model. Data systematic errors are apparent for instance in Fig. 16 (b)
(large narrow excursions) since the applied model (TCM) has no capacity to
describe, and the collision process itself is unlikely to generate, such
narrow structures. On the other hand, large and smooth deviations apparent in
Fig. 13 (right) unambiguously reveal model systematic errors. A more effective
venue for estimation of systematic model and data errors is differential
Z-scores based on statistical uncertainties that are unambiguous in their
definition.

##  VIII Discussion: two cultures

The blast-wave model represents one of two “cultures” relating to the dynamics
of high-energy nuclear collisions, associated with a hydrodynamics- or
hadronic-fluid-based description (flows) on the one hand and a particle or
high-energy physics (partons) version of QCD on the other. As noted in Sec.
VI.1 the former has its origins in the nuclear theory of Landau landau
(1950s) and elements of the fireball description of Hagedorn haged  (1960s)
which substantially predate the experimental discovery of quarks and
development of QCD theory (1970s) at which point it became evident that parton
cascades or showers play a major role in high-energy collision dynamics.

###  VIII.1 Parton splitting cascades or showers

One of the most popular Monte Carlo models of elementary collisions is PYTHIA
pythia1 ; pythia2 : “The [PYTHIA] program is designed to simulate the physics
processes that can occur in collisions between high-energy particles, e.g. at
the LHC collider at CERN. … A combination of perturbative results and models
for semihard and soft physics … are combined to trace the evolution towards
complex [e.g. hadronic] final states.” Initial developments that led to the
PYTHIA Monte Carlo – “reproducing the work of Field and Feynman, and extending
it to the analytical model developed in Lund” – were motivated by an early
model of jet formation reported in Refs. field1 ; field2 ; field3 . “The
‘jets’ observed in both cases [$e^{+}$-$e^{-}$ and $p$-$p$ collisions] are
thought to arise from quarks that fragment or cascade into a collection of
hadrons moving in roughly the direction of the original quark [emphasis
added]” field3 .

Concerning early development of detailed mechanisms “Both of these [early jet
models other than PYTHIA] were based on the concept of independent
fragmentation, wherein each of the $q$, ${\bar{q}}$ and $g$ jets are assumed
to fragment symmetrically around a jet axis defined by the direction of the
respective parton in the CM frame of the event. In Lund, instead another
picture had been developed, string fragmentation. Here the connecting colour
field is approximated by a massless relativistic string, with gluons
represented by pointlike momentum-carrying ‘kinks’.” “PYTHIA lacked the parton
showers that gave the other two programs [ISAJET and FIELDAJET] realistic jet
shapes.” “It was clear that parton showers would play a key role in order to
produce multijet final states…” pythia2 .

Reference pythia2  notes that “…the soft-gluon emission pattern around a
$q{\bar{q}}g$ topology could be viewed as a sum of radiation off two
independent dipoles…mimicking the nonperturbative string picture. …this
offered a starting point to formulate a [parton] shower as a successive
branching [cascade] of dipoles, an idea that today is a standard choice for
most shower algorithms….” In Ref. gosta  a Lund Dipole Cascade model,
implemented as the DIPSY Monte Carlo, is described: The model “gives a good
description of inclusive pp and ep cross-sections (including diffraction), and
a fair description of exclusive final states (min. bias and underlying event).
The model can also be applied to reactions with nuclei, with some early
results available.” In particular the model is able to describe $\eta$ charge
densities within $\eta\in[-2.5,5]$ for 0.9 and 7 TeV $p$-$\bar{p}$ collisions
at the LHC. In summary, various Monte Carlo models based on parton splitting
cascades or showers are able to describe the hadronic final states of
elementary ($e^{+}$-$e^{-}$, $p$-$p$) collisions in considerable detail.

###  VIII.2 Hadron gas in thermodynamic equilibrium

An interesting alternative approach to $e^{+}$-$e^{-}$ collision data is
reported in Ref. statmodel . The study acknowledges that $e^{+}$-$e^{-}$
collisions (as an example) include both hard processes (parton showers)
described by perturbative QCD and soft processes (parton fragmentation to
hadron jets) that must be described phenomenologically, with the Jetset Monte
Carlo (precursor to PYTHIA) as one example. Reference statmodel  comments:
“The main unsatisfactory feature of these [Monte Carlo] models…is the large
number of free parameters required in order to correctly reproduced
experimental data. As a consequence, those models have a rather poor
predictive power.”

The analysis assumes that “each jet [of an $e^{+}$-$e^{-}$ hadronic final
state] represents an independent phase in complete thermodynamic
equilibrium….” That assumption is taken to imply that “one can describe a jet
as an object defined by thermodynamic and mechanical quantities such as
temperature, volume….” The basis for analysis is “…the canonical partition
functions of systems with internal symmetries.” The analysis is applied to
final-state abundances of an array of hadron species. The study concludes that
“…this model is able to fit impressively well the average multiplicities of
light hadrons….” The model parameters are temperature $T$, volume $V$ and
$\gamma_{s}$ related to strange-quark chemical equilibrium.

While that analysis is informative and technically competent it can be argued
that “an independent phase in complete thermodynamic equilibrium” is not
demonstrated by model agreement with hadron species abundances. In discussing
the statistical properties of hadrons emerging from “freezeout” of A-A
collisions at the SPS Ref. stock  warns that “…this apparent ‘thermal’
equilibrium is a result of the decay process, the nature of which lies well
beyond the statistical model which ‘merely’ captures the apparent statistical
order, prevailing right after decay. … The observed equilibrium is, thus, not
achieved by inelastic transmutation of the various hadronic species densities,
in final hadron gas rescattering cascades, i.e. not by hadron rescattering
approaching a dynamical equilibrium [emphasis added].” One may argue by
analogy that the partons within a jet do not rescatter so as to achieve
thermodynamic equilibrium prior to hadronization. It is the hadronization
process itself, as a quantum transition following a least-action principle,
that leads to approximate agreement with the description in Ref. statmodel .

The array of species abundances among jet fragments from $e^{+}$-$e^{-}$
collisions (what is addressed by Ref. statmodel ) is only part of the
information carried by the fragment distribution from a jet ensemble. Single-
particle momentum distributions and multiparticle correlations (on angle
porter2 ; porter3  and $p_{t}$ or $y_{t}$ ytyt ), as well as local
correlations among charge, flavor and baryon number, clearly indicate that
jets are not featureless equilibrated gases. In summary, while thermal
equilibrium of a phase might imply a certain data profile, observation of such
a data profile by experiment does not demonstrate thermal equilibrium.

###  VIII.3 Hadron emission from a flowing source

Within a parton cascade the initial momentum of a single leading particle is
distributed by splitting to a shower of lower-momentum particles. Transfer of
momentum to final-state hadrons can be described as top-down. In contrast,
within a fluid-dynamics picture of high-energy collisions the initial
projectile energy is to some extent “stopped” in the CM leading to high energy
densities and pressure gradients that drive fluid flow. The resulting flow
field, convoluted with a locally-isotropic Boltzmann distribution (Cooper-Frye
formalism cooper ), is expected to describe final-state momenta of emitted
hadrons. Transfer of momentum/energy can be described as bottom-up.

Reference wong  discusses A-A collision dynamics at the LHC in the context of
the Landau hydrodynamic model of high-energy nuclear collisions: “Consider
first the case of central AA collisions with A $\gg$ 1 such that nucleons of
one nucleus collide with a large numbers of nucleons of the other nucleus and
the whole energy content is used in particle production. This is the case of
‘full stopping’ [emphasis added]. … Common to both Landau hydrodynamics and
Bjorken hydrodynamics is the basic assumption that in the [conjectured] dense
hot matter produced in high-energy heavy-ion collisions, the density of the
quanta of the medium is so high that a state of local thermal equilibrium can
be maintained through out.” The large densities imply corresponding large
density gradients that drive flows. In such models much of the hadron momentum
originates solely from pressure-induced flows combined with thermal emission.
Any small jet contribution is assumed restricted to limited $p_{t}$ intervals.

That is the context in which the BW model has been applied to particle data at
the SPS and higher energies. Reference ssflow  asserts “In this paper we want
to develop a phenomenological model for the hadronic matter by starting out
with thermalization as the basic assumption and adding more features as they
are dictated by the analysis of the measured hadronic spectra.” Referring to
nucleon flows at the Bevalac is the statement “Though this directed flow [at
the Bevalac] is of different nature as our collective expansion flow [at the
SPS], it suggests the relevance of hydrodynamics also at higher energies.”

But there immediately appears evidence of conflict between such assumptions
and what had emerged from high-energy (particle) physics study of elementary
collisions. Reference ssflow  first admits that “Our current understanding of
QCD results basically from high energy experiments with small collision
systems [e.g. $p$-$p$ collisions] suffering hard interactions which are
relatively easy to analyse.” But given the similarity of S-S and $p$-$p$
spectra at SPS energies Ref. ssflow  then argues that “It would be too
impulsive to deduce from the apparent similarity of pp and S+S spectra that
the collective interpretation [i.e. flows in S-S] is wrong, since a) pp is by
no means an elementary collision system which we understand in sufficient
detail to serve as the antipode [opposite] of a collective system and b) only
a few of the observed features of the [hadron] spectra can be fully reproduced
by these two radically different philosophies which share only a small set of
common principles as local energy-momentum conservation, relativistic space-
time picture, etc.”

One may assume that “radically different philosophies” refers to “high energy
experiments” and resulting parton-cascade descriptions (e.g. PYTHIA) on the
one hand and a hydrodynamic approach implemented via a BW model on the other.
If “pp is by no means an elementary collision system” then what is —
$e^{+}$-$e^{-}$? Point a) above essentially excludes $p$-$p$ collisions as a
reference system for the study of A-B collisions although that has
subsequently been common practice at the RHIC and LHC.

Concerning point b) is the statement “For other observables, which are not
directly tied to the dynamics, e.g. strangeness, we already see a big
enhancement compared to [minimum-bias] pp, thus indicating fundamental
differences between both [S-S and $p$-$p$] collision systems.” But that
conclusion can be strongly questioned in the context of measured jet
production in $p$-$p$ and for instance $p$-Pb collisions as reported in Refs.
ppbpid ; pidpart1 ; pidpart2 , in which minimum-bias jets are found to make
strong contributions to strangeness and baryon production. The “fundamental
differences” between S-S and $p$-$p$ are more likely simply due to increased
jet production in S-S collisions due to multiple N-N binary collisions (and
hence increased jet production) per nucleon participant with otherwise no
significant deviation from $p$-$p$ observations.

In summary, descriptions of hadron emission from a flowing source, e.g. the BW
model, compete directly with data descriptions in terms of parton cascades –
within projectile nucleons leading to (soft) hadron distributions along the
beam axis or within (hard) jets resulting from large-angle scattering of
partons from projectile nucleons.

###  VIII.4 BW model compared with elementary collisions

It is informative to consider some examples of “elementary” collisions in
comparison to model elements from the TCM. In particular, the
$\sqrt{s_{NN}}=19.4$ S-S spectrum data invoked in Ref. ssflow  are compared
with 17 GeV (SPS) $p$-$p$ data and with the soft-component model function that
describes $\sqrt{s}=200$ GeV $p$-$p$ spectra ppprd . These comparisons address
in part the question: what is an elementary collision and how should it be
modeled?

Figure 19 (left) shows a pion spectrum on $m_{t}$ (solid dots,
$2\times\pi^{-}$) for SPS fixed-target 0-2% central S-S collisions at 200A GeV
($\sqrt{s_{NN}}=19.4$ GeV) na35ss . A hydro-model (BW) analysis of those data
was used to infer radial flow with mean transverse speed
$\langle\beta_{t}\rangle\rightarrow\bar{\beta}_{t}\sim 0.25$ ssflow  (see Fig.
12, right, and associated text). The curve labeled $S_{NN}$ is the soft
component for $\sqrt{s}=$ 200 GeV NSD $p$-$p$ collisions from the RHIC ppprd .
The line labeled M-B is the Boltzmann (exponential) limit of Lévy distribution
$S_{NN}$ with $T\approx 145$ MeV. A $\pi^{+}+\pi^{-}$ spectrum from $p$-$p$
collisions at 158A GeV (open circles, $\sqrt{s}=17.3$ GeV) is included for
comparison na49pp .

Figure 19:  Left panel: $m_{t}$ spectra from 0-2% central S-S collisions at
200A GeV (solid points, $\sqrt{s_{NN}}=19.4$ GeV) na35ss  and from $p$-$p$
collisions at $\sqrt{s}=17.3$ GeV (open circles) na49pp . The dashed curve is
Lévy soft component $\hat{S}_{0}(m_{t})$ from 200 GeV $p$-$p$ collisions. The
dash-dotted curve is the Boltzmann limiting case of $S_{NN}$. Right panel:
$m_{t}$ spectrum (points) from $e^{+}$-$e^{-}$collisions at $\sqrt{s}=91$ GeV
alephff . The curves are duplicated from the left panel. The $e^{+}$-$e^{-}$
data are rescaled as in the axis label to match the $S_{NN}$ curve at larger
$m_{t}$.

Fig. 19 (right) shows an $m_{t}$ spectrum (points) from LEP
$e^{+}$-$e^{-}$collisions at $\sqrt{s}\approx 91$ GeV alephff .
$e^{+}$-$e^{-}$ hadronic events from Z0 decays are dominated by two-jet events
(see Fig. 13 of Ref. alephff ). The spectrum, derived from a sphericity
analysis of $q$-$\bar{q}$ dijets, reveals the jet fragment momentum
distribution transverse to the thrust (dijet) axis. The LEP
$dn_{ch}/p_{t}dp_{t}$ data were rescaled to overlap the $p$-$p$ soft component
$S_{NN}$ (dashed curve), and hence the SPS S-S spectrum. The shape of the
$m_{t}$ spectrum from 91 GeV Z0 decays is consistent with the soft component
of $p$-$p$ spectra at 200 GeV and the spectrum from central 19 GeV S-S
collisions. Commonality of the soft-component shape (Lévy distribution) across
energies and collision systems suggests that the TCM soft component is a
universal feature of fragmentation for any leading particle – be it parton or
hadron.

While the comparison in Fig. 19 (right) provides a hint of the relation among
dijet formation in $e^{+}$-$e^{-}$ collisions, the hadronic final state of A-B
nuclear collisions and the BW model a more direct comparison can be made.

Figure 20 (left) shows the $e^{+}$-$e^{-}$ $p_{t}$ spectrum from Fig. 19
(right) replotted on $y_{t}$ (derived from measured $p_{t}$ assuming all
hadrons are pions as for $y_{z}$ below). The points are published
$dn_{ch}/dp_{t}$ data (Fig. 18a of Ref. alephff ) divided by $p_{t}$ with no
other scaling. According to Ref. alephff  this is the “charged particle
momentum component transverse to the sphericity [dijet] axis and projected
into the event plane.” The solid curve is a Lévy distribution on $m_{t}$ [i.e.
$\hat{S}_{0}(m_{t})$] with slope parameter $T=90$ MeV and Lévy exponent
$n=7.8$. Dashed curves show the effect of varying $T$ by 10 MeV. For
comparison, soft-component $\hat{S}_{0}(y_{t})$ for pions in Fig. 1 (a) has
$T=145$ MeV and $n=8.5$. The $e^{+}$-$e^{-}$ data are described within their
published uncertainties except for the lowest three points at 0.02, 0.07 and
0.13 GeV/c.

Figure 20:  Left: The hadron distribution on $y_{t}$ (perpendicular to the
dijet axis) from 91 GeV $e^{+}$-$e^{-}$ collisions alephff . The solid curve
is a Lévy distribution [i.e. $\hat{S}_{0}(y_{t})$] with parameters ($T,n$)
noted in the figure. The dashed curves correspond to variation of $T$ by 10
MeV. Dotted curve M-B is the exponential (on $m_{t}$) limit for
$n\rightarrow\infty$ and $T=90$ MeV. Right: Hadron distribution on $y_{z}$
(parallel to the dijet axis) from 91 GeV $e^{+}$-$e^{-}$ $\rightarrow$
$q$-$\bar{q}$ collisions. The curve is Eq. (VIII.4) explained in the text.

Figure 20 (right) shows a hadron distribution on $y_{z}$ from $Z$ decays
reported by Ref. alephff . $y_{z}$ ($y_{s}$ in Fig. 16a of Ref. alephff  with
$s$ for “sphericity”) is derived from measured $p_{L}\rightarrow p_{z}$ along
the event sphericity (dijet) axis and is based on assigning the pion mass to
all hadrons. The solid curve is a generalized $q$-Gaussian distribution

$\displaystyle G(y_{z})$ $\displaystyle=$
$\displaystyle\frac{A}{[1+(|y_{z}|/\sigma_{y_{z}})^{\kappa}/n]^{n}}$
$\displaystyle\sim$ $\displaystyle
A\exp_{q}[-(|y_{z}|/\sigma_{y_{z}})^{\kappa}],$

with $\kappa=5$, $\sigma_{y_{z}}=3.35$, $n=4.2$ and $A=6.6$. The distribution
differs from a Gaussian in two ways: (a) the parameter value $\kappa=$ 5
rather than 2 leads to a more rectangular distribution and (b) the parameter
value $n=$ 4.2 controls the distribution tails. The $q$ in “$q$-Gaussian”
refers to the relation $q-1=1/n$. For $n\rightarrow\infty$ (and $\kappa=2$)
the distribution reverts to a standard Gaussian form. The kinematic limits of
the data distribution are $y_{z}\approx\pm 6.5$ (for pions) corresponding to
$\sqrt{s}/2\approx 46$ GeV. The “1/e” points at $\pm 3.35$ correspond to 2
GeV/c.

Several issues are addressed by Figs. 12, 19 and 20:

Figure 12 demonstrates the relation between the commonly-applied BW model and
(a) the Boltzmann distribution (agreement expected for $\bar{\beta}_{t}=0$ by
definition of BW model) and (b) $\hat{S}_{0}(y_{t})$ inferred from spectrum
data as the zero-density limit but requiring $\bar{\beta}_{t}\approx 0.25$
from the BW model to approximate the Lévy distribution.

Figure 19 demonstrates that the soft-component model (Lévy distribution)
$\hat{S}_{0}(y_{t})$ ($S_{NN}$ in those figures) describes spectra from SPS
$p$-$p$ and S-S collisions and (for higher $p_{t}$) LEP $e^{+}$-$e^{-}$
collisions, all of which deviate strongly from a Boltzmann-distribution
exponential for reasons related to known QCD processes. The small difference
between SPS central S-S spectra and $p$-$p$ spectra may be due an expected jet
contribution at $\sqrt{s_{NN}}\approx 19$ GeV jetspec2 . The power-law tail of
soft component $\hat{S}_{0}(y_{t})$ inferred from 200 GeV $p$-$p$ collisions
is apparently able to approximate the jet contribution to the 19 GeV S-S
spectrum, consistent with collision-energy dependence of the $p$-$p$ soft
component (see Ref. alicetomspec , Fig. 12 right).

Figure 20 reveals a direct connection between the TCM description of hadron-
hadron collisions and $e^{+}$-$e^{-}$ collisions from LEP alephff . The left
panel shows the $p_{t}$ spectrum from $e^{+}$-$e^{-}$ collisions (solid dots)
consistent in shape with soft component $\hat{S}_{0}(y_{t})$ (solid curve)
from $p$-$p$ and $p$-Pb collisions. While the soft component for pions has
$T\approx 145$ MeV that for $e^{+}$-$e^{-}$ collisions has the substantially
lower value $T\approx 90$ MeV. Whereas the parton parents of final-state
hadrons from A-B collisions experience Fermi motion within projectile nucleons
and $k_{t}$ contributions within a parton cascade the $q$-$\bar{q}$ pair in a
$e^{+}$-$e^{-}$ collision has zero transverse momentum in the CM frame leading
to a “cooler” hadron fragment spectrum. On the other hand $n=7.8$ for
$e^{+}$-$e^{-}$ collisions represents a “harder” power-law tail for that
spectrum than the $n=8.5$ value for $p$-$p$ collisions. Both observations are
consistent with $\hat{S}_{0}(y_{t})$, as the description inferred from A-B
data, reflecting hadron production via parton splitting cascades within a
conventional QCD context. The resulting Lévy data trend is inconsistent with
the Boltzmann exponential trend assumed for the BW model.

Figure 20 (right) shows the complementary dijet fragment distribution on
$y_{z}$ (from Ref. alephff  Fig. 16a) that can be compared with similar
distributions on pseudorapidity from hadronic A-B collisions. Within a
hydrodynamic model context such distributions are attributed within hadronic
collisions to longitudinal flow or Bjorken expansion. Yet the very similar
fragment distribution from $e^{+}$-$e^{-}$ collisions may be attributed to
parton splitting cascades leading to dijets. In the first case longitudinal
hadron momentum is attributed to pressure gradients, in the second case to a
parton parent in the cascade carrying a fraction $x$ of the leading-parton
momentum. One description (flows) says that a particle has been pushed while
another (partons) says that it has been “pulled.”

For both longitudinal and transverse momentum distributions from
$e^{+}$-$e^{-}$ collisions the trend at higher momentum is a power law
controlled by $n=4.2$ for $p_{z}$ and $n=7.8$ for $p_{t}$, indicated by
straight-line trends in Figure 20 on $y_{t}$ or $y_{z}$. A similar power-law
trend is observed for $p_{t}$ spectrum hard components [e.g. Fig. 4 (b,d)]
associated quantitatively in that case with large-angle-scattered partons
fragmenting to jets fragevo .

Detailed study of $p_{t}$ spectra from 5 and 13 TeV $p$-$p$ collisions
tomnewppspec  demonstrates the two-component (projectile-proton dissociation
and scattered parton fragmentation to dijets) nature of hadron production
consistent with basic QCD in such “elementary” collisions. Statistical
analysis of two models applied to the same data with intent to demonstrate
thermalization and/or flows in $p$-$p$ collisions tommodeltests  are falsified
via Z-score tests as in Sec. V. One might expect that “elementary” collisions
($e^{+}$-$e^{-}$, $p$-$p$ and even $p$-A systems) should be well understood
after fifty years of theoretical and experimental QCD-related research whereas
larger A-B collision systems might remain ambiguous. However, there seems to
be a persistent tendency to assume that A-A collisions comprise the fixed
reference and “elementary” collisions are seen as not elementary and still
returning unanticipated results.

##  IX Summary

Hadron $p_{t}$ spectrum analysis is a principal component of data
interpretation relating to high-energy nuclear collisions. A major element of
spectrum analysis is fitting models to spectrum data both as a form of data
compression (reducing spectra to a few, possibly interpretable, parameters)
and as a test of specific physical model predictions. Several models available
for spectrum description are comprised of simple functions: (a) the Boltzmann
exponential on transverse mass $m_{t}$, (b) the Tsallis model (“Tsallis
statistics”) on $m_{t}$ or $p_{t}$, (c) the blast-wave (BW) model on $m_{t}$
and (d) the two-component (soft+hard) spectrum model (TCM). Also available are
several complex Monte Carlo models (e.g. PYTHIA).

The BW model in particular, formulated to describe the effect on hadron
momenta in the center-of-momentum frame of moving (flowing) particle sources,
was intended to facilitate discovery of flows in A-A collisions as one
manifestation of QGP formation. Early results from the relativistic heavy ion
collider (RHIC) seemed to demonstrate flow manifestations in spectra.

More recently, the BW model as applied to $p$-A and even $p$-$p$ collisions
has been interpreted to support the same conclusions for small systems,
although matter and energy densities achievable in such systems are a priori
unlikely to produce a QGP. Such results cast doubt on the legitimacy of the BW
model applied to collision data.

The present study of BW model applications to identified-hadron (PID) spectra
from 5 TeV $p$-Pb collisions is a response to that situation. This study is
intended to assess basic formulation(s) of the BW model and its (their)
underlying assumptions, to determine data fit quality in relation to possible
model falsification and to compare the BW model to alternative models that may
provide superior data description and lead to different (and perhaps more
credible) physical interpretations.

The general procedure includes the following steps: (a) visually compare BW
model and TCM to data in basic semilog plots on transverse rapidity $y_{t}$
where data trends have a particularly simple structure, (b) compare data and
model shapes (what is emphasized in BW model applications) using model-
independent measures, (c) assess data-model fit quality differentially and
quantitatively using Z-scores, and (d) review BW model evolution across
several decades via a number of published examples.

The specific results are as follows: (a) As established in a previous study
the TCM describes spectra for four hadron species within point-to-point
uncertainties over the entire $p_{t}$ acceptance and the TCM is not fitted to
individual spectra. In contrast, the BW model fails to describe spectra except
over limited $p_{t}$ intervals defined by data-model agreement. That failure
is most obvious for the case of pion spectra. (b) Data-model comparisons using
model-independent methods (e.g. logarithmic derivatives) demonstrate
qualitative differences between BW model and data. On that basis the model is
falsified. (c) Application of the Z-score statistic, a standard measure of
data-model agreement, demonstrates that the BW model is excluded as a spectrum
model for this collision system, whereas the TCM is an acceptable model as
established in a previous study. (d) While the general structure of the BW
model is accepted, specific applications may be quite different in detail
(what parameters are allowed to vary over what ranges, even basic algebraic
structure) leading to difficulty in interpreting results.

If the BW model is accepted as a valid physical data model, albeit only within
artificially-restricted $p_{t}$ intervals (with doubtful justification), it is
commonly concluded that the great majority of hadron production arises from a
moving source exhibiting “collective” motion. However, the present study
establishes that for the $p$-Pb collision system at 5 TeV the BW model is
falsified by PID spectrum data. Previous studies have arrived at the same
conclusion for 5 and 13 TeV $p$-$p$ collisions.

In contrast, the TCM does provide a statistically-acceptable data description
for all available PID spectrum data from 5 TeV $p$-Pb collisions. From that
result the following conclusions may be drawn:

The TCM soft component is consistent with nuclear transparency (1976) wherein
a projectile nucleon may transit a target nucleus intact, albeit “excited”
(wounded nucleon model), with later dissociation to hadron fragments outside
the collision space-time volume. That picture is consistent with the soft-
component model inferred from data – a Boltzmann exponential with power-law
tail indicating a heterogeneous source (e.g. parton shower).

The TCM hard component is quantitatively consistent with minimum-bias dijet
production as represented by measured jet spectra and parton fragmentation
functions. In that case final-state hadrons do indeed emerge from moving
particle sources which however happen to be scattered partons as expected from
basic QCD theory.

##  Appendix A BW model derivation

Formulation of the BW model for single-particle momentum spectra requires
evaluation of the integral cooper

$\displaystyle E\frac{d^{3}N}{dp^{3}}$ $\displaystyle=$
$\displaystyle\frac{g}{(2\pi)^{3}}\int_{\sigma}e^{-(u^{\mu}p_{\mu})/T}p^{\nu}d^{3}\sigma_{\nu},$
(23)

where a factor $\exp(\mu/T)$ is omitted. Evaluation requires determination of
quantities $u^{\mu}(x)$, $p_{\mu}$, $\sigma_{\mu}$ and $d^{3}\sigma_{\mu}$,
where $u^{\mu}(x)$ is a fluid velocity field on space-time $x$ in the CM or
lab frame, $p_{\mu}$ is particle four-momentum in the lab frame and
$\sigma(r,\phi_{x},\eta_{s})$ is the volume (hypersurface) from which detected
particles emerge (freeze out). The Boltzmann exponential assumes particle
emission according to an isotropic thermal distribution in the local frame.

###  A.1 Basic four-vectors

The basic four-vectors required for evaluating Eq. (23) are defined using a
set of self-consistent symbols that overlap as much as possible the different
conventions utilized in an assortment of relevant publications.

####  A.1.1 Space-time rapidity

The position four-vector in cylindrical coordinates preferred for large
$\sqrt{s}$ is

$\displaystyle x_{\mu}$ $\displaystyle=$ $\displaystyle[t,x,y,z]$
$\displaystyle=$
$\displaystyle[\tau\cosh(\eta_{s}),r\hat{e}(\phi_{x}),\tau\sinh(\eta_{s})],$

where $\hat{e}(\phi_{x})$ is a unit vector, with space-time rapidity

$\displaystyle\eta_{s}$ $\displaystyle\equiv$
$\displaystyle\frac{1}{2}\ln\left(\frac{t+z}{t-z}\right)=\ln\left(\frac{t+z}{\tau}\right)$
(25)

where $t=\tau\cosh(\eta_{s})$, $z=\tau\sinh(\eta_{s})$ and
$v_{s}=z/t=\tanh(\eta_{s})$. Azimuth angles in configuration and momentum
spaces are distinguished by $\varphi$ vs $\phi$ or $\phi_{x}$ vs $\phi_{p}$,
with the latter notation used below to reduce ambiguity. The space-time
parameter $\zeta$ in Ref. ssflow  here goes to $\eta_{s}$.

Configuration-space rapidity $\eta_{s}$ is distinguished from (a) momentum-
space pseudorapidity $\eta$ with no subscript and (b) longitudinal flow angle
or boost $\eta_{z}$ (in configuration space). Those distinctions are
consistently maintained below unless explicitly noted. The approximation
$\eta_{z}(z)\approx\eta_{s}(z)$ (Bjorken expansion) is often assumed. In what
follows the two quantities are maintained distinct unless specifically noted
otherwise.

####  A.1.2 Fluid (particle source) velocity field

Fluid velocity fields are distributed on configuration space $x$. The
following are Lorentz hyperbolic (boost) angles for particle sources moving in
the CM frame

$\displaystyle\eta_{t}(r)$ $\displaystyle=$
$\displaystyle(1/2)\ln\left(\frac{1+\beta_{t}}{1-\beta_{t}}\right)=\ln[\gamma_{t}(1+\beta_{t})]$
(26) $\displaystyle\eta_{z}(z)$ $\displaystyle=$
$\displaystyle(1/2)\ln\left(\frac{1+\beta_{z}}{1-\beta_{z}}\right)=\ln[\gamma_{z}(1+\beta_{z})],$

with $\beta_{x}=\tanh(\eta_{x})$, $\gamma_{x}=\cosh(\eta_{x})$ and
$\gamma_{x}^{2}=1/(1-\beta_{x}^{2})$. Given those definitions the four-
velocity field distributed on configuration space is

$\displaystyle u^{\mu}(x)$ $\displaystyle=$
$\displaystyle\cosh(\eta_{t})[\cosh(\eta_{z}),\tanh(\eta_{t})\hat{e}(\phi_{x}),\sinh(\eta_{z})]$
(27) $\displaystyle=$
$\displaystyle\gamma_{t}[\gamma_{z},\beta_{t}\hat{e}(\phi_{x}),\gamma_{z}\beta_{z}]$
$\displaystyle u^{\mu}u_{\mu}$ $\displaystyle=$ $\displaystyle 1$

In Refs. ssflow ; cleymans  $\rho\rightarrow\eta_{t}$ and
$\eta\rightarrow\eta_{z}$ (“boost angles”) are functions of $r$ and $z$
respectively based on assumptions.

####  A.1.3 Particle momenta in the CM or lab frame

Given measured momenta in cylindrical coordinates $(p_{t},\theta,\phi)$ (e.g.
in a solenoidal magnetic field) pseudorapidity
$\eta=-\ln\left[\tan(\theta/2)\right]$ ($\sim y_{z}$ rapidity). The three-
vector momentum is $\vec{p}=(p_{x},p_{y},p_{z})=(p_{t},\eta,\phi_{p})$ with
magnitude $p=p_{t}\cosh(\eta)$ and $p_{z}=p_{t}\sinh(\eta)$. Symbol $\eta$
without subscript here always represents pseudorapidity.

Transverse mass $m_{t}$ is defined by $m_{t}^{2}=p_{t}^{2}+m_{0}^{2}$.
Transverse rapidity is $y_{t}=\ln[(m_{t}+p_{t})/m_{0}]$ in a longitudinally
comoving frame with $m_{t}=m_{0}\cosh(y_{t})$, $p_{t}=m_{0}\sinh(y_{t})$.
Longitudinal rapidity is $y\rightarrow y_{z}=\ln[(E+p_{z})/m_{t}]$ with
$E=m_{t}\cosh(y_{z})$, $p_{z}=m_{t}\sinh(y_{z})$. With those definitions the
particle four-momentum is

$\displaystyle p_{\mu}$ $\displaystyle=$
$\displaystyle[E,p_{t}\hat{e}(\phi_{p}),p_{z}]$ $\displaystyle=$
$\displaystyle\ [m_{t}\cosh(y_{z}),p_{t}\hat{e}(\phi_{p}),m_{t}\sinh(y_{z})]$
$\displaystyle=$ $\displaystyle
m_{t}[\cosh(y_{z}),\tanh(y_{t})\hat{e}(\phi_{p}),\sinh(y_{z})]$ $\displaystyle
p^{\mu}p_{\mu}$ $\displaystyle=$ $\displaystyle m_{0}^{2}$

####  A.1.4 Particle energy in the boost frame

With the velocity field for particle emission and the particle four-momentum
defined the product represents particle energy in a conjectured boost frame
based on measured particle momentum in the CM or lab frame cooper .

$\displaystyle u^{\mu}p_{\mu}$ $\displaystyle=$
$\displaystyle\gamma_{t}m_{t}[\cosh(y_{z})\cosh(\eta_{z})-\sinh(y_{z})\sinh(\eta_{z})]$
$\displaystyle-\gamma_{t}\beta_{t}p_{t}\cos(\phi_{x}-\phi_{p})$
$\displaystyle=$
$\displaystyle\gamma_{t}[m_{t}\cosh(y_{z}-\eta_{z})-\beta_{t}p_{t}\cos(\phi_{x}-\phi_{p})]$
$\displaystyle\rightarrow$ $\displaystyle m_{t}\cosh(y_{z})=E~{}~{}\text{in
the lab frame for no flows}.$

###  A.2 Emission hypersurface

The freezeout hypersurface $\sigma$ appearing in Eq. (23) is the
configuration-space volume from which particles are emitted (decoupled from
the source fluid) and then fly freely to a detector. In Ref. cooper  the
differential product $p^{\mu}d\sigma_{\mu}=(E-\beta p)dx=\bar{E}d\bar{x}$ for
a 1D system, where bars denote quantities in the local, comoving or boost
frame. That follows since $p_{\mu}=[E,p]$, $\sigma_{\mu}=[t,x]$,
$d\sigma_{\mu}=[dx,dt]$ and $dt=\beta dx$ if $d\bar{t}=0$. In what follows,
several examples are presented for defining $d^{3}\sigma_{\mu}$ and
determining $p^{\mu}d^{3}\sigma_{\mu}$.

####  A.2.1 Schnedermann et al. ssflow

In this study the emission hypersurface (with $\zeta\rightarrow\eta_{s}$ the
space-time rapidity, $\eta_{s}\in[-H_{s},H_{s}]$ and $r\in[0,R]$) is

$\displaystyle\sigma(r,\phi_{x},\eta_{s})$ $\displaystyle\rightarrow$
$\displaystyle\sigma_{\mu}=[t(\eta_{s}),r\hat{e}(\phi_{x}),z(\eta_{s})].$ (30)

The differential freezeout volume is then

$\displaystyle d^{3}\sigma_{\mu}$ $\displaystyle=$
$\displaystyle\left[\frac{\partial z}{\partial\eta_{s}},0,0,\frac{\partial
t}{\partial\eta_{s}}\right]d\eta_{s}rdrd\phi_{x}$ (31)

assuming instantaneous freezeout in $r$, and with

$\displaystyle p^{\mu}d^{3}\sigma_{\mu}$ $\displaystyle=$
$\displaystyle\left[E\frac{\partial z}{\partial\eta_{s}}-p_{z}\frac{\partial
t}{\partial\eta_{s}}\right]d\eta_{s}rdrd\phi_{x}.$ (32)

The velocity field and four-momentum $u^{\mu}(x)$ and $p^{\mu}$ are as in App.
A.1.2, with “boost angles” $\eta\rightarrow\eta_{z}(z)$ and
$\rho(r)\rightarrow\eta_{t}(r)$ and with longitudinal momentum
$p_{L}\rightarrow p_{z}$. Azimuth angles are $\varphi\rightarrow\phi_{p}$ and
$\phi\rightarrow\phi_{x}$. If $z=\tau\cosh(\eta_{s})$ and
$t=\tau\sinh(\eta_{s})$, with $E=m_{t}\cosh(y_{z})$ and
$p_{z}=m_{t}\sinh(y_{z})$, then

$\displaystyle p^{\mu}d^{3}\sigma_{\mu}$ $\displaystyle\rightarrow$
$\displaystyle m_{t}\tau\cosh(y_{z}-\eta_{s})d\eta_{s}rdrd\phi_{x}.$ (33)

####  A.2.2 Tomasik et al. tomasik

In this study the single-particle momentum spectrum is defined in terms of a
“source function” $S(x,p)$

$\displaystyle E\frac{d^{3}N}{dp^{3}}$ $\displaystyle=$ $\displaystyle\int
d^{4}xS(x,p).$ (34)

The source function is defined by

$\displaystyle S(x,p)d^{4}x$ $\displaystyle\propto$ $\displaystyle
m_{t}\tau\cosh(y_{z}-\eta_{s})e^{-u^{\mu}(x)p_{\mu}/T}$
$\displaystyle\times~{}G(r)H(\eta_{s})d\eta_{s}rdrd\phi_{x}$

with $y\rightarrow y_{z}$, $\eta\rightarrow\eta_{s}$ and
$\varphi\rightarrow\phi_{x}$. The function $H(\eta_{s})$ is a Gaussian with
width $\Delta\eta_{s}$. The function $G(r)$ is Gaussian with width $R_{G}$ or
flat over $r\in[0,R_{B}]$. A unit-normal Gaussian on $\tau$ is assumed
integrated to yield mean value $\tau$. The emission hypersurface is defined by

$\displaystyle\sigma_{\mu}$ $\displaystyle=$
$\displaystyle[\tau\cosh(\eta_{s}),r\hat{e}(\phi_{x}),\tau\sinh(\eta_{s})]$
(36)

with (assuming $d\tau=0$)

$\displaystyle d^{3}\sigma_{\mu}$ $\displaystyle=$
$\displaystyle\tau\left[\cosh(\eta_{s}),0,0,\sinh(\eta_{s})\right]d\eta_{s}rdrd\phi_{x}.$
(37)

Four-vectors $u_{\mu}$ and $p_{\mu}$ are as in App. A.1.2, with flow angles
$\eta_{l}\rightarrow\eta_{z}(z)$ and $\eta_{t}(r)$ as above in Eq. (26).
Factors $m_{t}\tau\cosh(y_{z}-\eta_{s})$ in Eq. (A.2.2) result if
$d^{3}\sigma_{\mu}$ defined above is combined with $p_{\mu}$ as in App. A.1.2.
Note that here $\eta_{t}(r)=\tanh(\beta_{t})\propto r$ instead of
$\beta_{t}\propto r$ as in Ref. ssflow .

####  A.2.3 Florkowski and Broniowski Ref. wojciechs

In this study the emission hypersurface (with $\alpha_{\|}\rightarrow\eta_{s}$
introduced as the space-time rapidity) is

$\displaystyle\sigma_{\mu}$ $\displaystyle=$
$\displaystyle[\tau(\zeta)\cosh(\eta_{s}),\rho(\zeta)\hat{e}(\phi_{x}),\tau(\zeta)\sinh(\eta_{s})],~{}$
(38)

where $\zeta\in[0,1]$ is here simply a parameter that correlates proper time
$\tau$ and radius $\rho$, not space-time rapidity $\eta_{s}$. The associated
volume element (with $\phi\rightarrow\phi_{x}$) is

$\displaystyle d^{3}\sigma^{\mu}$ $\displaystyle=$
$\displaystyle\left[\frac{\partial\rho}{\partial\zeta}\cosh(\eta_{s}),\frac{d\tau}{d\zeta}\hat{e}(\phi_{x}),\frac{d\rho}{d\zeta}\sinh(\eta_{s})\right]$
$\displaystyle\times\tau(\zeta)d\zeta\rho(\zeta)d\eta_{s}d\phi_{x}.$

The velocity field $u_{\mu}(x)$ [with
$\alpha_{\perp}(\zeta)\rightarrow\eta_{t}(\zeta)$] is consistent with Eq. (27)
except $\alpha_{\|}\rightarrow\eta_{z}$ (boost angle), not $\eta_{s}$ (space-
time rapidity), and an error in factorization has been corrected. Particle
momentum $p_{\mu}$ is consistent with App. A.1.2 if $m_{\perp}\rightarrow
m_{t}$, $y\rightarrow y_{z}$ and $\varphi\rightarrow\phi_{p}$. In that case

$\displaystyle p^{\mu}d^{3}\sigma_{u}$ $\displaystyle=$
$\displaystyle\left[m_{t}\cosh(y_{z}-\eta_{s})\frac{\partial\rho}{\partial\zeta}-p_{t}\cos(\phi_{p}-\phi_{x})\frac{d\tau}{d\zeta}\right]$
(40) $\displaystyle\times\rho(\zeta)\tau(\zeta)d\zeta d\eta_{s}d\phi_{x}.$

If one assumes $d\tau/d\zeta\rightarrow 0$ as in Ref. wojciechs  and
$\rho\rightarrow r$ then

$\displaystyle d^{3}\sigma_{\mu}$ $\displaystyle\rightarrow$
$\displaystyle\tau\left[\cosh(\eta_{s}),0,0,\sinh(\eta_{s})\right]d\eta_{s}rdrd\phi_{x}$
(41)

which is consistent with Eq. (37), and

$\displaystyle p^{\mu}d^{3}\sigma_{u}$ $\displaystyle\rightarrow$
$\displaystyle
m_{t}\tau(\zeta)\cosh(y_{z}-\eta_{s})d\eta_{s}r(\zeta)drd\phi_{x}.~{}~{}$ (42)

To clarify, the distinction between space-time rapidity and longitudinal boost
angle for Ref. wojciechs  is as follows: $\alpha_{\|}$ is initially introduced
as space-time rapidity ($\rightarrow\eta_{s}$). But $\alpha_{\perp}$ is then
initially introduced as transverse boost angle ($\rightarrow\eta_{t}$) with
transverse flow $v_{r}\rightarrow\beta_{t}=\tanh{\eta_{t}}$, and $\alpha_{\|}$
is later reintroduced as longitudinal boost angle $\rightarrow\eta_{z}$) with
longitudinal flow $v_{z}\rightarrow\beta_{z}=\tanh(\eta_{z})$. The initial
assumption is that $\eta_{s}\approx\eta_{z}$ with $\alpha_{\|}$ filling both
roles. In this appendix space-time $\eta_{s}$ and velocity-field $\eta_{z}$
are maintained distinct unless indicated otherwise.

####  A.2.4 Rath et al. cleymans

In this study $u_{\mu}$ and $p_{\mu}$ are as in App. A.1.2 with
$\eta\rightarrow\eta_{z}$ (longitudinal flow angle),
$\rho\rightarrow\eta_{t}$, $y\rightarrow y_{z}$, $\phi\rightarrow\phi_{p}$ and
$\phi_{r}\rightarrow\phi_{x}$. $d^{3}\sigma_{\mu}$ is as in Eq. (37) except
$\eta\rightarrow\eta_{s}$ (space-time rapidity) and with a minus sign
corrected. With those definitions

$\displaystyle p^{\mu}d^{3}\sigma_{\mu}$ $\displaystyle=$ $\displaystyle
m_{t}\tau[\cosh(y_{z})\cosh(\eta_{s})-\sinh(y_{z})\sinh(\eta_{s})]$
$\displaystyle\times d\eta_{s}rdrd\phi_{x}$ $\displaystyle=$ $\displaystyle
m_{t}\tau\cosh(y_{z}-\eta_{s})d\eta_{s}rdrd\phi_{x}$

and

$\displaystyle u^{\mu}p_{\mu}$ $\displaystyle=$
$\displaystyle\gamma_{t}[m_{t}\cosh(y_{z}-\eta_{z})-\beta_{t}p_{t}\cos(\phi_{x}-\phi_{p})].~{}~{}~{}~{}~{}~{}$
(44)

The “Bjorken correlation in rapidity” is specified as
$(y\approx\eta)\rightarrow(y_{z}\approx\eta_{s})$ rather than the conventional
$\eta_{z}\approx\eta_{s}$. That assignment of $y_{z}$ is incorrect [see Eq.
(A.2.4)].

In summary, the above examples illustrate the differences and similarities
among several approaches to the BW model over nearly thirty years, from 1993
to 2020.

###  A.3 Single-particle spectra and blast-wave model

The several definitions of relevant quantities in the previous subsections are
combined to evaluate Eq. (23) with

$\displaystyle E\frac{d^{3}N}{dp^{3}}$ $\displaystyle=$
$\displaystyle\frac{d^{3}N}{m_{t}dm_{t}dy_{z}d\phi_{p}}.$ (45)

For each case in App. A.2 Eq. (23) is evaluated with the relevant quantities
and consistent symbols adopted herein. Distinction between longitudinal flow
angle $\eta_{z}$ and space-time rapidity $\eta_{s}$ is maintained unless
stated otherwise. But just as $\beta_{t}(r)=\tanh^{-1}(\eta_{t})\sim r$ is
often assumed $\beta_{z}(z)=\tanh^{-1}(\eta_{z})\sim z$ may also be assumed in
the form of Bjorken expansion $\eta_{z}(z)\approx\eta_{s}=\ln[(t+z)/\tau]$.

####  A.3.1 Schnedermann et al. ssflow

The formulation in Eq. (14) of Ref. ssflow  is modified by
$\zeta\rightarrow\eta_{s}$, ${\cal Z}\rightarrow H_{s}$, $y\rightarrow y_{z}$,
$\rho\rightarrow\eta_{t}$, $\eta\rightarrow\eta_{z}$ and invoking the
expression in Eq. (33), in which case

$\displaystyle E\frac{d^{3}N}{dp^{3}}$ $\displaystyle=$
$\displaystyle\frac{g}{(2\pi)^{3}}m_{t}\tau\int_{-H_{s}}^{H_{s}}d\eta_{s}\cosh(y_{z}-\eta_{s})$
$\displaystyle\times\int_{0}^{R}rdr\exp\left[-m_{t}\cosh(\eta_{t})\cosh(y_{z}-\eta_{z})/T\right]$
$\displaystyle\times\int_{0}^{2\pi}d\phi_{x}\exp[p_{t}\sinh(\eta_{t})\cos(\phi_{x}-\phi_{p})/T].$

For the integrations, integrals over $\phi_{p}$ and $\phi_{x}$ give factor
$(2\pi)^{2}I_{0}\left[p_{t}\sinh(\eta_{t})/T\right]$. Then
$\cosh(y_{z}-\eta_{s})\rightarrow\cosh(y_{z}-\eta_{z})\cosh(\eta_{z}-\eta_{s})$,
omitting the odd term $\sinh(y_{z}-\eta_{z})\sinh(\eta_{z}-\eta_{s})$ that
will not contribute to the $y_{z}$ integral. Both sides are integrated over
longitudinal rapidity $y_{z}$ leading to factor
$2K_{1}[m_{t}\cosh(\eta_{t})/T]$. The integral over $\eta_{s}$ then yields a
constant no matter what the relation between $\eta_{s}$ and $\eta_{z}$,
although equality (Bjorken expansion) is usually assumed in which case the
integral over $\eta_{s}\rightarrow 2H_{s}$ ($=2Z_{t}$ in Ref. ssflow ). The
result, assuming $\eta_{t}(r)\sim r$, is

$\displaystyle\frac{dN}{m_{t}dm_{t}}$ $\displaystyle\approx$
$\displaystyle\frac{g}{\pi}m_{t}\tau 2H_{s}\times$