Review of nonflow estimation methods and uncertainties in relativistic heavy-ion collisions

Abstract

Collective anisotropic flow, where particles are correlated over the entire event, is a prominent phenomenon in relativistic heavy-ion collisions and is sensitive to the properties of the matter created in those collisions. It is often measured by two- and multi-particle correlations and is therefore contaminated by nonflow, those genuine few-body correlations unrelated to the global event-wise correlations. Many methods have been devised to estimate nonflow contamination with various degrees of successes and difficulties. Here, we review those methods pedagogically, discussing the pros and cons of each method, and give examples of ballpark estimate of nonflow contamination and associated uncertainties in relativistic heavy-ion collisions. We hope such a review of the various nonflow estimation methods in a single place would prove helpful to future researches.

1 Introduction

Collective anisotropic flow is a hallmark of heavy-ion (nucleus-nucleus) collisions. It refers to correlations among particles over the entire event – every particle is correlated with every other particle. One particular example is hydrodynamic flow, where the initial geometry anisotropy in finite impact parameter collisions is converted by interactions into momentum space anisotropy of particles [1]. All particles are thus correlated to the ellipse-shaped collision geometry. The anisotropy is nonzero even in head-on (zero impact parameter) collisions because of position fluctuations of nucleons inside the colliding nuclei giving rise to finite eccentricities [2, 3]. Because of the same reason, nonzero anisotropy can also emerge in small-system collisions, such as proton-proton (pp), proton-nucleus (pA), deuteron-nucleus (dA), and helium-nucleus collisions [4, 5].

The interactions in the system created in relativistic heavy-ion collisions at RHIC (Relativistic Heavy-Ion Collider) and the LHC (Large Hadron Collider), presumably the quark-gluon plasma (QGP), are governed by quantum chromodynamics (QCD) [6, 7, 8]. How exactly those interactions convert the initial-state geometry anisotropy into the final-state momentum anisotropy is not well settled. It is generally believed that ultra-strong interactions are required to produce the observed large anisotropy (or flow) in heavy-ion collisions, and the QGP created in those collisions is a nearly perfect fluid [9] and can be well described by hydrodynamics with little viscosity [10]. In peripheral heavy-ion collisions and small-system collisions, the interactions may not be intense enough where hydrodynamics could be applicable and the escape mechanism may be at work [11, 12, 13, 14]. It is important to keep in mind that the word “flow”, while suggestive, does not necessarily mean hydrodynamic flow.

The collision geometry–the impact parameter vector of a heavy-ion collision or generally the geometric shape of the interaction zone–is experimentally unknown. It is often reconstructed as a proxy from the final-state particle azimuthal distribution as the symmetry harmonic plane [15]. The particle azimuthal distribution can be written in the Fourier series [16]:

$$\frac{dN}{d\phi}=\frac{N}{2\pi}\left(1+2\sum_{n=1}^{\infty}v_{n}\cos n(\phi-\psi_{n})\right)\,,$$ (1)

where $\phi$ is the particle azimuthal angle and $\psi_{n}$ is that of the $n^{\rm th}$ order harmonic plane. The anisotropic harmonic flow is then given by

$$v_{n}=\langle{\cos n(\phi-\psi_{n})}\rangle\,.$$ (2)

Here $n=1,2,3,...$; $v_{1}$ is called directed flow, $v_{2}$ elliptic flow, and $v_{3}$ triangular flow. Because $\psi_{n}$ is reconstructed from particles, Eq. (2) is equivalent to two-particle correlations. Namely, it can be obtained from two-particle cumulant, or the Fourier coefficient $V_{n}$ of the two-particle correlation function,

$$\frac{dN_{\rm pair}}{d\Delta\phi}=\frac{N_{\rm pair}}{2\pi}\left(1+2\sum_{n=1}^{\infty}V_{n}\cos n\Delta\phi\right)\,,$$ (3)

where $\Delta\phi=\phi_{1}-\phi_{2}$ is the azimuthal angle difference of a particle pair. The Fourier coefficient is simply

$$V_{n}=\langle{\cos n\Delta\phi}\rangle\equiv c_{n}\,,$$ (4)

where $c_{n}$ is also used in literature [17, 18]. Under the presence of only collective flow, it follows straightforwardly that

$$V_{n}=v_{n}^{2}\{2\}\,,$$ (5)

where $v_{2}\{2\}$ denotes two-particle cumulant flow.

If measured by the two-particle cumulant method of Eqs. (3,4), the $V_{n}$ unnecessarily must be related to the collision geometry. For example, the gluon field from each incoming nucleus can be correlated and this correlation may result in a final-state anisotropy $V_{n}$ [19, 20]. Such an initial-state correlation permeates over an entire collision event, so it is flow, but not part of hydrodynamic flow.

There is, however, an important contribution to $V_{n}$ (or more generally, any variables) measured via correlations, and that is nonflow: two- and few-particle genuine correlations that have nothing to do with the collision geometry [21, 22, 23]. Examples of nonflow correlations are those between daughter particles from a resonance decay, particles from a jet shower originated from an energetic parton (quark or gluon) [24, 25], hadrons from string fragmentation [26], Hanbury-Brown Twiss (HBT) interferometry [27], and global momentum conservation [28]. Figure 1 illustrates those correlations in a two-particle angular correlation plot in $(\Delta\eta,\Delta\phi)$ where $\Delta\eta=\eta_{1}-\eta_{2}$ is the two-particle pseudorapidity difference. Most of the nonflow correlations (HBT, resonance decays, intra-jet correlations, string fragmentation) are small-angle correlations contributing to the near-side ($\Delta\phi\sim 0$) peak at $\Delta\eta\sim 0$. HBT is relevant at very small angle differences and generally considered minor. There is a long-range contribution from dijet correlations contributing to the away-side ridge at $\Delta\phi\sim\pi$ that has little $\Delta\eta$ dependence (because of the stochastic sampling of kinematics of the underlying parton-parton scattering producing the dijet). It is noteworthy that the nonflow contributions from jets are those intra- and inter-jet hadron correlations, whereas the azimuthal anisotropy of jet-axis orientations, a result of path-length dependent partonic energy loss or jet quenching [29, 30], is related to the collision geometry and is thus part of flow. Global momentum conservation contributes to the away-side ridge, but only to the first harmonic [28].

Figure 1: Illustration of nonflow correlations in two-particle $(\Delta\eta,\Delta\phi)$ angular distribution. The near-side peak at $\Delta\phi\sim 0$ is mainly short ranged ($\Delta\eta\sim 0$). There could be various near-side short-range correlations: HBT is very short ranged (sharp peak), resonance decay and intra-jet correlations are short ranged with a $\Delta\eta$ width on the order of one unit. The away-side ridge at $\Delta\phi\sim\pi$ is long ranged and comes mainly from dijet correlations and global momentum conservation. Underlying these nonflow peaks are the majority pair distribution (zero-suppressed), modulated by anisotropic flows of various harmonic orders in $\Delta\phi$, which could have weak dependencies on $\eta$ (on single particle level) and $\Delta\eta$ (on two-particle level, referred to as flow decorrelation).

All those nonflow correlations contribute to the $V_{n}$ of Eq. (4), so the factorization in Eq. (5) is no longer valid. A keen interest in relativistic heavy-ion collisions is to measure collective anisotropic flows arising from final-state interactions to probe the properties of the QGP by, for example, comparing to hydrodynamic calculations. Thus, the goal is to measure the collective anisotropic flows from global event-wise correlations related to the collision geometry. To this end, nonflow contamination in $V_{n}$ must first be subtracted with faithful systematic uncertainties [31]. Many experimental methods have been devised to estimate and/or subtract nonflow contributions. This note aims to give a pedagogical review of those methods in a single place, hopefully useful for future researches. The possible contributions from geometry-unrelated initial-state gluon correlations [19, 20] need also to be considered. This is however outside the scope of this note.

2 Nonflow estimation methods

2.1 $\Delta\eta$-gap Methods

Nonflow correlations are primarily short ranged in $\Delta\eta$. Those short-range nonflow correlations can be suppressed by applying a $\Delta\eta$ gap between particle pairs used in the correlation analysis [15]. The $\Delta\eta$ gap cannot suppress those long-range nonflow correlations.

2.1.1 Simple $\Delta\eta$-gap method

The simplest way to suppress nonflow is to impose a $\Delta\eta$ gap between the two hadrons in two-particle cumulant measurement of $V_{n}$ [15]. This method is easy and straightforward. The shortcomings are obvious:

1. the method is not clean–how much nonflow is eliminated depends on the $\Delta\eta$-gap size relative to the width of the near-side short-range correlations;

2. the away-side dijet correlations cannot be eliminated and how much they contribute to $V_{n}$ is unknown a priori;

3. $\eta$-dependent $v_{n}$ [32, 33] and longitudinal decorrelations [34, 35, 36] would yield different measurements with different $\Delta\eta$ gaps.

The last is not related to nonflow, but flow and flow fluctuations. However, since they cannot be distinguished by a $\Delta\eta$-gap analysis, the effects of flow and nonflow are mixed. All those shortcomings make the uncertainty estimation of a particular $\Delta\eta$-gap result difficult.

2.1.2 Two-subevent method

Instead of applying $\Delta\eta$-gap between the two particles, one may use the two-subevent method where one particle is taken from one subevent in a given $\eta$ region and the other from another subevent in a different $\eta$ region, and the two subevents are separated in $\eta$ with a certain $\Delta\eta$ gap [15]. This has the advantage that the cumulant calculation can be applied which involves only single-particle loops [37, 38], thus saving tremendous computing time compared to double loops of particle pairs. Obviously, all the shortcomings of the simple $\Delta\eta$-gap method described in Sect. 2.1.1 are present in the two-subevent method.

2.1.3 Three- and four-subevent methods

The away-side jet-correlations are long ranged in $\Delta\eta$, and thus the simple $\Delta\eta$-gap method and the two-subevent method cannot suppress the away-side nonflow correlations. One can use three- and four-subevent method to calculate multi-particle cumulants to suppress the away-side nonflow correlations from dijets by separating those subevents well in $\eta$ [39, 40, 41]. Figure 2 illustrates the idea behind the three-subevent method of four-particle cumulant, where dijets appear in at most two subevents and the third subevent is free of nonflow dijet correlations once the $\Delta\eta$-gaps are made wide enough. As a result, the correlations among the particles–two from one subevent and one from each of the other two subevents–are due only to collective flow. The four-subevent method for four-particle cumulant is similar, except that one particle is taken from each of the four subevents.

Figure 2: Sketch of the three-subevent method. Four-particle cumulant is formed by taking two particles from one subevent and one particle from each of the other two subevents. All nonflow correlations are suppressed including back-to-back dijet correlations, except global momentum conservation.

The multi-subevent method suppresses essentially all nonflow correlations (except that in $V_{1}$ from global momentum conservation [28]). This method also takes advantage of cumulant calculations [37, 38], greatly reducing the computing demand. The shortcomings are as same as those from the simple $\Delta\eta$-gap method described in Sect. 2.1.1, except (ii).

2.2 Low-multiplicity subtraction methods

The general idea behind this category of methods is that particle correlations in peripheral and low-multiplicity collisions are primarily due to nonflow, and those in central and high-multiplicity collisions come from collective flow with some contamination from nonflow correlations. The key questions are:

• at what low multiplicity the correlations can be considered all as nonflow, and

• how this low-multiplicity nonflow can be modeled or utilized to gauge nonflow in high-multiplicity/central events.

2.2.1 Inverse multiplicity scaling method

Resonance abundances are expected to scale approximately with the final-state multiplicity. If the scaling is exact and if resonance kinematic distributions, and thus the average correlations between decay daughters, do not vary with event multiplicity, then nonflow contribution from resonance decays to $V_{n}$ by Eq. (4) is inversely proportional to multiplicity ($N$). Take the number of particles to be $N_{\rm low}$ in low-multiplicity events and $N_{\rm high}$ in high-multiplicity events, and take the number of pairs to be $N^{2}$ (valid under Poisson statistics), then

$$N^{2}_{\rm high}V_{2}^{\rm high} = N^{2}_{\rm high}V_{2}^{\rm sub}+N^{2}_{\rm low}V_{2}^{\rm low}\cdot\frac{N_{\rm high}}{N_{\rm low}}\,,$$ (6)

$$V_{2}^{\rm sub} = V_{2}^{\rm high}-\frac{N_{\rm low}}{N_{\rm high}}V_{2}^{\rm low}\,.$$ (7)

Here, it is assumed that the number of nonflow “sources” scales with multiplicity, and those nonflow particles, while keeping the nonflow correlations among themselves, have attained flow individually in high-multiplicity events. The notation $V_{2}^{\rm sub}$ denotes low-multiplicity subtracted $V_{2}$ with the ultimate goal to be the nonflow-subtracted anisotropic flow. This method is sometimes called the $c_{0}$ method where $c_{0}\equiv N$ [17, 18].

The $1/N$ scaling is a strong assumption for nonflow contribution from jet correlations:

• Jets are produced by hard processes whose abundance increases with $N$ more strongly than linearly [29, 30]. The $1/N$ scaling would be an underestimate of this part of nonflow. However, the majority nonflow contribution comes from relatively low transverse momentum ($p_{T}$) jets or minijets, whose production may be more closely proportional to $N$.

• Jets are modified by the nuclear medium created in relativistic heavy-ion collisions. Such modifications broaden/suppress the angular correlations between jet fragments [42, 43, 44], making nonflow weaker.

• On the other hand, jets lose energy via collisional and radiative partonic energy loss [45, 46], resulting in more lower $p_{T}$ particles and particle pairs, and thus stronger nonflow correlations.

Quantitatively, these effects are $p_{T}$ dependent, and the interplay among them determines the final dependency of nonflow on collision centrality/multiplicity. Simulations by the hijing model (a jet production model without hydrodynamic flow) indicate only a modest excess of nonflow than the $1/N$ scaling by 10–20% from peripheral to central Au+Au collisions at $\sqrt{s_{\textsc{nn}}}=200$ GeV [31]. Nevertheless, this is one inferior part of the method.

The $1/N$ scaling implies not only a scaling in the abundance of nonflow sources, but also that the physics of nonflow correlations does not change from low- to high-multiplicity collisions. In addition, the method explicitly assumes that the correlations in low-multiplicity collisions are all nonflow. The flow result in high-multiplicity events extracted from this method will thus inevitably depend on what low-multiplicity events are considered as the pure-nonflow baseline. The selection of very low-multiplicity or very peripheral events seems to be the natural choice, but they could be biased towards too soft underlying nucleon-nucleon scatterings or diffractive interactions. This is a selection bias [47, 48]. The nonflow in those events may not be a good reflection of nonflow in central collisions; in other words, the lowest multiplicity collisions are unnecessarily the best baseline for nonflow subtraction. This issue is more severe in analyzing anisotropies in small-system collisions, such as pp, pA, and dA collisions. Because of multiplicity selection biases [47, 48], the low- and high-multiplicity collisions of those small systems can be vastly different in physics, including nonflow correlations. High-multiplicity pp events are likely biased towards jet production, whereas low-multiplicity events are likely biased towards softer-than-average interactions.

To get a quantitative feeling, we examine Au+Au collisions at the top energy at RHIC. The measured $v_{2}$ in peripheral 70-80% and top 5% central Au+Au collisions are $v_{2}\{2\}^{\rm peri}\approx 6.9\%$ and $v_{2}\{2\}^{\rm cent}\approx 2.4\%$, respectively [49]. Assuming the former is all nonflow, the $1/N$ scaling of nonflow would result in a nonflow of $(v_{2}^{\rm peri})^{2}\frac{N_{\rm ch}^{\rm peri}}{N_{\rm ch}^{\rm cent}}\approx 1.7\times 10^{-4}$, i.e. 30% nonflow in $(v_{2}^{\rm cent})^{2}$ of the latter. Choosing 90-100% centrality as the low-multiplicity baseline would yield a smaller nonflow estimate [50]. However, it is possible that the 90-100% centrality, which is determined by event multiplicity, is biased towards too soft interactions, thus underestimating the nonflow to be used for central collisions. On the other hand, the higher 70-80% centrality data used above likely contain flow, and thus the 30% estimate is likely an overestimate of nonflow.

Multiplicity selection bias is likely insignificant in central heavy-ion collisions–there are many underlying nucleon-nucleon interactions. It is unlikely, up to certain high-multiplicity limit, to have all those interactions to produce jets, for example. Therefore, minimum-bias (MB) pp collisions would be the best baseline at our disposal for nonflow subtraction for central heavy-ion collisions within the scope of the $1/N$ scaling method. STAR has measured accumulative correlations ($NV_{2}$) as functions of $p_{T}$ in pp and Au+Au collisions (Fig. 1 in Ref. [51]). The ratio (pp/Au+Au) can be directly related to nonflow fraction in Au+Au under the assumption of $1/N$ scaling. The peripheral 80-100% Au+Au collisions are comparable to pp, suggesting that these peripheral collisions are dominated by nonflow. In the top 5% central collisions, the measured ratio is on average 15%, or 12% after correcting for the occupancy-dependent detector efficiency of 20% from pp to central Au+Au [52]. One may consider that the mean $\langle{p_{T}}\rangle$ in central Au+Au collisions is larger than in pp, and because $NV_{2}$ is a strongly increasing function of $p_{T}$, the $N$-weighted pp/Au+Au ratio (or the ratio of $NV_{2}$ at the corresponding $\langle{p_{T}}\rangle$’s) is smaller, approximately 8% [53]. However, since the multiplicity $N$ is already included in the $NV_{2}$ measure, this extra weighting may not be justified. On the other hand, one may argue that the pp/Au+Au ratio should be taken at higher $p_{T}$ for pp than for Au+Au because of jet energy loss in the latter, which would yield a larger ratio. These considerations suggest that nonflow estimations by the $1/N$ scaling have large uncertainties because of strong assumptions made in this method.

In Eq. (6), it is implicitly assumed that all particles in central events have attained the global collective flow by the $N_{\rm high}^{2}$ in front of $V_{2}$. This is reasonable for heavy-ion collisions but may not be justified in small-system collisions. It is possible that those nonflow particles, such as jets produced early, have exited the collision zone without participating in the final-state interactions that are responsible for the generation of flow. In such a case, the $V_{2}$ estimated by Eq. (7) would be an underestimate of flow in those high-multiplicity small-system collisions.

To recap, the advantage of the $1/N$ scaling method is that it is easy and straightforward. The shortcomings are the inherent assumptions, namely,

• it assumes that all correlations are nonflow in low-multiplicity events, which begs the question how low in multiplicity is a good baseline;

• it assumes the $1/N$ scaling of nonflow correlations, implying no change in the physics of nonflow correlations from low- to high-multiplicity collisions, a strict proportionality to $N$ of the abundance of nonflow sources, and no selection biases for those low-multiplicity events;

• it assumes that those nonflow particles, while still genuinely correlated among themselves, have also attained individually the same collective flow as the rest of the collision event.

The strong assumption of ii) and the fact that how nonflow correlations vary with the collision centrality/multiplicity make it difficult to assess the robustness and the associated uncertainties of the estimated nonflow by the simple $1/N$ scaling method.

2.2.2 Near-side jet yield scaled subtraction

Nuclear effects, like jet quenching, result in modifications to jet-correlations [24], thus nonflow effects change from peripheral to central collisions. Particle production mechanisms and distributions change, for example, yielding the baryon-over-meson enhancement from peripheral to central collisions [54, 55, 56], which would modify the nonflow effects from resonance decays. Any those changes will cause nonflow to deviate from the simple $1/N$ scaling.

Both jet-correlations and resonance decays yield a near-side correlation peak at $(\Delta\eta,\Delta\phi)\sim(0,0)$. Dijets contribute to an away-side correlation at $\Delta\phi\sim\pi$ but more or less uniform in $\Delta\eta$. One may take the difference between small $\Delta\eta$ (short-range) and large $\Delta\eta$ (long-range) correlations, properly normalized, to arrive at a near-side correlated yield (primarily composed of contributions from resonance decays and intra-jet fragments). One can then compare this near-side yield in high-multiplicity events to that in low-multiplicity events. Any difference would indicate modification of nonflow correlations from low- to high-multiplicity collisions, multiplicity selection biases [47, 48], or likely both. One may take the ratio of the correlated near-side yields ($Y_{\rm high}/Y_{\rm low}$) as a scaling factor to apply on to the low-multiplicity nonflow in Eq. (7) to possibly take into account modifications to nonflow correlations [57, 48, 58]. The estimated flow in high-multiplicity events would then be

$$V_{2}=V_{2}^{\rm high}-\frac{N_{\rm low}}{N_{\rm high}}\frac{Y_{\rm high}}{Y_{\rm low}}V_{2}^{\rm low}\,.$$ (8)

Such a scaling, while an improvement to the simple $1/N$ scaling, comes with its own issues:

• the near-side nonflow correlation shape may change/broaden, the effect of which is not included by the simple scaling of the correlated yields;

• the away-side jet-correlated yield unlikely scales with the near-side one because of “trigger” biases (for example, surface bias at high $p_{T}$), and even if it scales, the scaling factor $Y_{\rm high}/Y_{\rm low}$ includes resonance decays besides near-side jet contribution and is therefore already an incorrect scaling factor, needless to say that the away-side correlation shape can be significantly modified because of large jet-quenching effects;

• the correlated yield analysis is difficult at low $p_{T}$ because of large combinatorial background.

2.2.3 Template fit

In the template fit method [59, 17], the two-particle correlation is assumed to be composed of a scaled correlation from low-multiplicity events (considered all as nonflow) and a series of Fourier harmonics (to represent flow) except the $n=1$ component. This is illustrated in Fig. 3. Namely, the following function is fit to two-particle correlations in high-multiplicity events,

$$\frac{dN_{\rm pair}^{\rm high}}{d\phi}=F\frac{dN_{\rm pair}^{\rm low}}{d\phi}+\frac{G}{2\pi}\left(1+2\sum_{n=2}^{\infty}V_{n}\cos n\Delta\phi\right)\,,$$ (9)

where $F$ is a fit parameter to scale the low-multiplicity correlations. Note that in the harmonic series, the first harmonic is excluded. In other words, the fitting constraint is to assume that the first harmonic is all nonflow; there is no flow harmonic of the first order. (Of course, without any imposed constraint, the fit would not yield any unique result because anything can be described by a Fourier series.) In effect, the template fit method assumes the nonflow component in high-multiplicity events to scale according to $V_{1}$. If $V_{1}$ is inversely proportional to multiplicity (which is approximately the case, see e.g. Ref. [58]), then the template fit method is similar to the $1/N$ scaling method in Sect. 2.2.1. Any deviation of the fit parameter $F$ from the $1/N$ scaling factor can be considered as a modification in nonflow correlations from low to high multiplicity, taken care of by a constant scaling without change in shape. This is similar to the near-side jet yield scaled subtraction described in Sect. 2.2.2, except that the scaling here is determined by the $V_{1}$ dipole component (mainly the away-side correlation, probably dominated by global momentum conservation [28]) instead of the near-side correlation.

Figure 3: Sketch of the template fit method and the dipole method. The heavy-ion pair correlations are considered to be composed of nonflow correlations scaled from pp and flow correlations. In the template fit, the extracted flow magnitudes are with respect to the pair multiplicity between the two curves, i.e., the nonflow pairs scaled from pp is first subtracted. In the dipole method, the flow magnitudes are with respect to the total pair multiplicity underneath the upper curve. The methods also apply to small system analysis, where the label “heavy-ion” is to be replaced by “high-multiplicity events” and “scaled pp” by “low-multiplicity events.”

Within fit errors, the following equities follow from Eq. (9):

$$N_{\rm pair}^{\rm high} = FN_{\rm pair}^{\rm low}+G\,,$$ (10)

$$N_{\rm pair}^{\rm high}V_{1}^{\rm high} = FN_{\rm pair}^{\rm low}V_{1}^{\rm low}\,,$$ (11)

$$N_{\rm pair}^{\rm high}V_{n}^{\rm high} = FN_{\rm pair}^{\rm low}V_{n}^{\rm low}+GV_{n}^{\rm sub}\,.$$ (12)

Simple algebra leads to

$$V_{n}^{\rm sub}=\left(V_{n}^{\rm high}-\frac{V_{1}^{\rm high}}{V_{1}^{\rm low}}V_{n}^{\rm low}\right)\left/\left(1-\frac{V_{1}^{\rm high}}{V_{1}^{\rm low}}\right)\right.\,.$$ (13)

Eq. (13) is similar to Eqs. (7,8) except the denominator. The denominator comes about because the normalization $G$ in Eq. (12) is not all the pairs in the high-multiplicity event but that excluding those nonflow particle pairs. In other words, the physics picture is that the high-multiplicity event is composed of nonflow particles (the white area in Fig. 3) and flow particles (the hatched area in Fig. 3) and the nonflow particles do not participate in the flow.

The template fit method attempts to account for nonflow changes by scaling the correlations in low-multiplicity events. It is similar in spirit to the near-side jet yield scaled subtraction method, differing in the assumption of how nonflow changes. So a similar set of assumptions and shortcomings, namely:

• it assumes that the first harmonic is all nonflow and there is no flow $V_{1}$, so nonflow scales from low to high multiplicity according to the $V_{1}$ magnitude, without change in the correlation shape;

• the away-side dijet correlations likely change shape due to jet quenching in medium-to-central heavy-ion collisions, and thus unlikely scale with $V_{1}$;

• the near-side nonflow correlations may not scale like $V_{1}$, and the correlation shape can also change from low- to high-multiplicity collisions.

• it assumes that the nonflow pairs do not participate in the final-state collective flow, and the extracted flow is relative to the reduced average pair multiplicity–the normalization $G$ in Eq. (9) or the area between the two curves in Fig. 3–after subtraction of those nonflow pairs. However, the cross pairs (one from a nonflow source and the other from the collective flow bulk) are included in the normalization $G$ and thus are considered as part of flow. This normalization issue is mostly a small effect but may be relevant for small-system flow studies.

The last assumption is unique to the template fit method. For all other nonflow estimation methods, the final extracted flows are based on the total event (pair) multiplicity, including nonflow particles–see bullet iii) in Sect. 2.2.1 (not explicitly stated for other nonflow estimation methods).

2.2.4 Dipole method

This method [17, 18] is very similar to the template fit method in Sect. 2.2.3, also assuming the first harmonic $V_{1}$ to be all nonflow, and the $V_{1}$ in low-multiplicity events is scaled to match that in high-multiplicity events, as in Eq. (11). The $V_{n}^{\rm sub}$ is, however, defined in terms of the total pair multiplicity $N_{\rm pair}^{\rm high}$, not $G$ as in Eq. (12). Namely, $N_{\rm pair}^{\rm high}V_{n}^{\rm high}=FN_{\rm pair}^{\rm low}V_{n}^{\rm low}+N_{\rm pair}^{\rm high}V_{n}^{\rm sub}$, which leads to

$$V_{n}^{\rm sub}=V_{n}^{\rm high}-\frac{V_{1}^{\rm high}}{V_{1}^{\rm low}}V_{n}^{\rm low}\,.$$ (14)

This assumes effectively that all those nonflow particles, while keeping their inter-particle nonflow correlations, have gained collective flow themselves individually. This is as same as the assumption made in Sects. 2.2.1 and 2.2.2. In the template fit method, on the other hand, the low-multiplicity correlated yield is scaled and subtracted first, and those particles are not counted as part of collective flow as aforementioned. Effectively, the dipole method is as same as the template method except that the scaled low-multiplicity correlation is first lowered to average zero. The dipole method is sometimes called “$c_{1}$ method” (where $c_{1}\equiv V_{1}$) [17, 18].

2.3 Data-driven fitting methods

As seen from Sects. 2.1 and 2.2, nonflow contamination is difficult to remove, involving strong assumptions. Another approach to estimate nonflow contamination is quite different from those described in Sects. 2.1 and 2.2, and is data driven. It relies on correlation structures observed in data, and identifies local peaks and ridges to attribute to nonflow correlations and the more or less $\Delta\eta$-independent underlying distribution to collective flow. It then involves fitting the correlation data with pre-defined functional forms in a data-driven way.

The advantage of these fitting methods is that they are data driven, having minimal assumptions about the physics of nonflow and no reliance on the evolution of nonflow over centrality/multiplicity. The downside, however, is that the fit functional forms to nonflow correlations are ad hoc and that the fitting can sometimes be tedious.

2.3.1 2D fitting in $(\Delta\eta,\Delta\phi)$

One data-driven technique is to perform 2-dimensional $(\Delta\eta,\Delta\phi)$ fits to two-particle correlations [60, 44, 61, 62, 63]. Raw two-particle correlations are predominantly of a triangular shape in $\Delta\eta$, which comes from the approximately uniform single-particle density distribution within limited $\eta$ acceptance (for example, $|\eta|<1$ in the STAR experiment [64]). The triangular acceptance can be largely corrected by the mixed-event technique, resulting in two-particle $(\Delta\eta,\Delta\phi)$ correlations, such as the one sketched in Fig. 1, that are dominated by an overall pedestal. The fine structures over the pedestal are nonflow correlations and a flow modulation along $\Delta\phi$ (with possibly a $\Delta\eta$ dependence). These features have been generally described in the introduction.

The near-side nonflow peaks can be modeled by 2D Gaussians, and the flow modulation is described by a Fourier series. The Fourier coefficients, corresponding to flow harmonics, are typically assumed to be $\Delta\eta$ independent. The two-particle correlations can then be fit with a functional form including Gaussians and Fourier harmonics. This 2D fitting method is data-driven; the nonflow correlation shapes are dictated by the data structure and modeled. The method fits the data in a given centrality or multiplicity bin and does not rely on assumptions using low-multiplicity baseline events. A recent 2D fit study indicates an approximately 40% nonflow in central isobar collisions [62, 63], in line with the central Au+Au data considering multiplicity dilution of nonflow.

It is noteworthy that the Fourier flow harmonics can be $\Delta\eta$ dependent, for example, due to flow decorrelations [34, 35, 65, 36]. In addition, the $v_{n}$, although long ranged, can be $\eta$ dependent [32]; such a dependence can result in $\Delta\eta$-dependent Fourier coefficients as well. Including these possible $\eta$-dependencies of the flow harmonics as free parameters in the fit model would not be fruitful as the correlation data would not have sufficient constraining power over these parameters (including those modeling nonflow). Usually the effects due to the possible $\eta$-dependent flow harmonics are assessed as part of the systematic uncertainties.

However, the information of the $\Delta\eta$-dependent harmonics can be obtained from other analyses and can be factored into the 2D fit function. We leave discussions on this to the next section.

2.3.2 1D fitting in $\Delta\eta$

The spirit of the 1D fitting method is similar to that of the 2D fitting–decomposing the flow and nonflow in the inclusive measurement of $V_{n}(\Delta\eta)$ as a function of $\Delta\eta$. Figure 4 illustrates a typical $V_{2}(\Delta\eta)$ dependence on $\Delta\eta$. It is generally a decreasing function of $\Delta\eta$, primarily because of the short-range nonflow contributions. For tracking detectors, track merging is sometimes important, which creates a dip at $\Delta\eta\sim 0$ as indicated by the dashed curve in Fig. 4. When track merging effects are significant, the average $V_{n}$ from full-event analysis can be even smaller than subevent analysis which automatically excludes the dip. In such a case, the difference between full-event and subevent $V_{n}$ in the $\Delta\eta$-gap methods would not be a correct assessment of systematic uncertainties.

Compared to the 2D fitting method in Sect. 2.3.1, 1D fitting is more straightforward because the $\Delta\phi$ dimension is collapsed in the $V_{n}$ measurement. It is also easier to include the possible $\Delta\eta$ dependence of flow. There are two sources for $\Delta\eta$ dependence of flow as aforementioned: $\eta$-dependent single-particle flow and flow $\Delta\eta$ decorrelations.

The single particle flow $v_{n}$ can depend on $\eta$ [32]. Such a dependence can be examined, for the case of $n=2$ for example, by four-particle cumulant measurement $v_{2}\{4\}(\eta)$ where nonflow is largely eliminated, and by $v_{2}\{{\textsc{zdc}}\}(\eta)$ with respect to the first order harmonic plane of spectator neutrons measured by Zero-Degree Calorimeters (ZDC) which is free of nonflow. Note that $v_{2}\{{\textsc{zdc}}\}$ measures approximately the average $\langle{v_{2}}\rangle$ and $v_{2}^{2}\{4\}\approx\langle{v_{2}}\rangle^{2}-\sigma^{2}$ where $\sigma$ denotes Gaussian width of $v_{2}$ fluctuations, whereas the flow component in the two-particle measurement $V_{2}$ is $v_{2}^{2}\{2\}=\langle{v_{2}}\rangle^{2}+\sigma^{2}$ [15]. As an approximation, one may obtain $v_{2}\{2\}(\eta)$ from $v_{2}\{4\}(\eta)$ and $v_{2}\{{\textsc{zdc}}\}(\eta)$. Alternatively, one may make the reasonable assumption that the effect of $v_{2}$ fluctuations is proportional to $\langle{v_{2}}\rangle$ and apply a scale factor as a free fit parameter to $v_{2}\{4\}(\eta)$ or $v_{2}\{{\textsc{zdc}}\}(\eta)$ to obtain $v_{2}\{2\}(\eta)$. The $\Delta\eta$ dependence can be obtained by folding two single-particle distributions of $v_{2}(\eta_{1})$ and $v_{2}(\eta_{2})$.

Another source of $\Delta\eta$ dependence is the flow decorrelation, $r_{n}(\Delta\eta)\approx 1-F_{n}\Delta\eta$ [65, 36]. Such decorrelations can result from fluctuations of the harmonic flow magnitude as well as harmonic planes over $\Delta\eta$ [66]. The decorrelation parameter $F_{n}$ of a few percent has been measured at RHIC [67] and the LHC [65, 36].

The $\Delta\eta$ dependence of the flow component in two-particle $V_{2}(\Delta\eta)$ is the combination of the two effects above. This is illustrated by the magenta curve in Fig. 4. The nonflow component in $V_{2}(\Delta\eta)$ is modeled by Gaussians in $\Delta\eta$. Fits are applied to $V_{2}(\Delta\eta)$ excluding the track merging $\Delta\eta\sim 0$ region, or alternatively including a negative Gaussian at $\Delta\eta=0$ to model it, treating the flow fluctuation scale factor as one parameter in the fit. The blue curve is an illustration of the fitted flow, scaled up by the fit parameter from the shape given in magenta. The idea of the fit is to have the best description of $V_{2}(\Delta\eta)$ with the known shape of flow and Gaussian-shaped nonflow.

Figure 4: Illustration of the 1D fit method. The black curve indicates a $V_{2}$ measurement, where the dip at $\Delta\eta\sim 0$ indicated by the dashed curve is caused by track merging in a typical tracking detector. The magenta curve incorporates $\Delta\eta$ dependencies from single-particle $v_{2}(\eta)$ and flow decorrelation $1-F_{n}\Delta\eta$. The blue curve indicates the flow component in $V_{2}$ scaled up from the magenta curve to include flow fluctuations, with the scaling factor treated as a fit parameter. The difference between the black curve and the blue curve is nonflow contribution to the $V_{2}$ measurement. The vertical axis is zero-suppressed and the values are only order-of-magnitude indications.

The away-side dijet correlations are always the most notorious nonflow to handle. If the away-side correlations in $\Delta\eta$ are just like flow correlations that both are relatively uniform, then 2D fitting cannot distinguish them. One may resort to 2D fitting if the away-side correlations are relatively sharp-peaked in $\Delta\phi$ such that one is reasonably confident that it is nonflow. If it is broadly distributed, then there is really no way to distinguish it from flow as Fourier components can describe any functional shape. However, it is more likely that the away-side correlations in $\Delta\eta$ differ from flow; for instance, away-side jet correlations could be more peaked at $\Delta\eta=0$ at relatively high $p_{T}$, or at low $p_{T}$ perhaps dipped at $\Delta\eta=0$. As long as the away-side correlations in $\Delta\eta$ are different from flow, then it is possible to distinguish it in 1D fitting by examining the data structure in $V_{2}(\Delta\eta)$.

2.3.3 Symmetry method in $(\eta_{1},\eta_{2})$

STAR performed another data-driven analysis by examining $V_{2}(\eta_{1},\eta_{2})$ as a function of the two particles’ pseudorapidities [68]. This method exploits the $\eta$ reflection symmetry in symmetric heavy-ion collisions by comparing two pairs, one at $(\eta_{1},\eta_{2})$ and the other at $(\eta_{1},-\eta_{2})$, and thus does not assume any particular shapes for flow as functions of $\eta$. A $\Delta\eta$-dependent component and a $\Delta\eta$-independent component are identified in $V_{2}(\eta_{1},\eta_{2})$. The former is associated with the combination of nonflow and $\Delta\eta$-dependent flow fluctuations. Linear fits are performed to data as functions of $2\eta_{2}$ and the linear fit parameters are examined as functions of $\eta_{1}$. This leads to a refined functional form as the sum of an exponential and a Gaussian in $\Delta\eta$ to describe the data well. The $\Delta\eta$-independent part is associated with flow plus $\Delta\eta$-independent flow fluctuations. They are further found to be independent of $\eta$ within the limited STAR acceptance of $|\eta|<1$. STAR has also supplemented the analysis with the four-particle cumulant as a function of the particles’ pseudorapidities, albeit large uncertainties, to separate flow fluctuations from the average flow magnitude.

This method has the least model assumption and the functional dependencies on $\Delta\eta$ and $\eta$ are determined from data. The disadvantages of the method are the required high statistics–because $V_{2}(\eta_{1},\eta_{2})$ and especially the four-particle cumulant $V_{2}$ are analyzed in 2-dimensions–and subsequently the large uncertainties, and the reduction of data to the functional forms to describe the various components of flow and nonflow requires a high level of attention. The estimated nonflow in central Au+Au collisions from this method is on the order of 20% with relatively large uncertainty [68].

3 Discussions and Summary

Nonflow includes all few-body correlations in a collision event except the global flow correlations where all particles are correlated over the entire event. Those nonflow correlations cannot possibly be fully measured, and various estimation/subtraction methods have been devised, as reviewed here in a single place. They can be broadly divided into three categories: (1) $\Delta\eta$-gap, (2) low-multiplicity subtraction, and (3) data-driven fits. All of them have assumptions, some of which are strong and some are less so; they come with different pros and cons. Because of the various degrees of assumptions, assessments of systematic uncertainties on nonflow are challenging.

1. The $\Delta\eta$-gap methods, while easy to implement, are not clean. The interpretation of results of a given $\Delta\eta$-gap analysis is subject to issues like the $\Delta\eta$-decorrelation and $\eta$-dependence of flow. The systematic uncertainties are hard to quantify–comparing to full-event (no $\Delta\eta$ gap) benchmark can, for example, be bitten by the usual track-merging detector artifact. An improvement would be to examine the results as functions of the applied $\Delta\eta$-gap size, which would constitute the data-driven 1D fitting method.

2. The low-multiplicity subtraction methods come with the large uncertainty in the assumptions of the multiplicity/centrality dependence of nonflow and the arbitrariness in the choice of low-multiplicity event class. Strictly confining within a given set of assumptions, one may arrive at a nonflow estimate with relatively small systematic uncertainty. However, this is only valid when the assumptions are correct; more likely, the estimated nonflow is systematically biased, and loosing the assumptions would yield a wide range of uncertainties.

3. The data-driven fitting methods are probably the best at our disposal with the fewest assumptions and least model dependency. However, it requires careful study of data and identification of nonflow correlation shapes. The fitting is often tedious and requires thorough attention and systematic assessment.

In general, nonflow contamination is severe in peripheral collisions and become less so towards more central collisions. However, since the elliptic flow also decreases with increasing centrality because of the more spherical collision zone, nonflow contamination in central collisions can still be appreciable, whereas it is generally the smallest in midcentral collisions. In the top 5% central Au+Au collisions at RHIC, for example, the data-driven fitting methods indicate a nonflow fraction of 20% with a typical relative systematic uncertainty of 20% [68, 69, 63, 62]. From low-multiplicity subtraction, one can get a wide range of nonflow estimate, probably 10–30%, depending on assumptions of multiplicity dependence of nonflow and what low-multiplicity events are taken as nonflow baseline (see Sect. 2.2.1). The former depends on many factors, including a faster increase in the number of nonflow sources than the multiplicity and a change in the correlation of each nonflow source with multiplicity/centrality.

Nonflow subtraction is particularly challenging in small system collisions to search for signs of collectivity in high-multiplicity events. The reasons are many-fold: the nonflow contamination is high and likely dominates the $V_{2}$ measurement, multiplicity selection biases are significant for both low- and high-multiplicity events, how nonflow varies from low- to high-multiplicity collisions is largely unknown partially because of the selection biases. The task is easier at the LHC than at RHIC because of the larger detector longitudinal acceptances, the larger multiplicities produced, and the more likelihood to have collectivity at the higher energies of LHC. At RHIC the analysis is significantly more difficult; for example, the difference between PHENIX [70, 71] and STAR [17, 18] on flow in p+Au, d+Au, and ³He+Au is not fully settled.

The difficulties in estimating/subtracting nonflow reflects the fact that nonflow cannot be thoroughly measured and the physics evolution of nonflow with collision system and centrality/multiplicity is not well understood. It is therefore important to examine various estimation/subtraction methods, when strong assumptions are involved, in order to arrive at a robust nonflow estimate with faithful systematic uncertainties, which often takes the majority effort of data analysis. This is particularly important when data are compared to theoretical models, such as hydrodynamics, where nonflow effects are not fully incorporated, to draw quantitative physics conclusions.