1 Introduction

The large azimuthal anisotropy observed for particles produced in heavy-ion collisions at RHIC [1,2,3,4] and the LHC [5,6,7,8] is one of the main signatures of the formation of strongly interacting matter called quark–gluon plasma (QGP). A standard picture of an ultrarelativistic heavy-ion collision is that the initial, asymmetric ‘almond’ shape of the colliding nuclei’s overlap region leads to the formation of pressure gradients in the QGP. These pressure gradients transform the initial shape into an azimuthal anisotropy of the final-state particle distributions through a nearly ideal hydrodynamic evolution and subsequent QGP hadronisation process [9]. The azimuthal anisotropy is customarily decomposed into Fourier components with the amplitude of the nth term denoted by \(v_n\) and known as a flow harmonic [10]. Theoretical hydrodynamical models successfully describe observed flow phenomena at low particle transverse momenta [11]. The properties of QGP were recently studied with measurements of correlations between flow harmonics of different order [12,13,14,15,16] as well as with analyses of event shapes [16,17,18,19,20]. It is expected that in lead–lead (Pb+Pb) collisions the magnitudes of the azimuthal flow harmonics [6, 7] should be correlated with the mean transverse momentum \([p_{\mathrm{T}}]\) of the particles on an event-by-event basis [21]. In this paper, that correlation is called the \(v_n\)\([p_{\mathrm{T}}] \) correlation. In proton–lead (\(p\)+Pb) collisions, the measurements of multi-particle correlations [22] show evidence of collective phenomena. The spectra of identified particles in \(p\)+Pb collisions are consistent with a presence of the radial flow [23] while the nuclear modification factor at high \(p_{\mathrm{T}}\) approaches unity [24]. Despite intensive studies, the mechanism responsible for the collective behaviour in small collision systems still remains unknown [9]. In \(p\)+Pb collisions the \(v_n\)\([p_{\mathrm{T}}] \) correlation could provide constraints on the initial geometry of the particle source, thereby reducing the overall modelling uncertainty. According to the hydrodynamical model predictions [25], in \(p\)+Pb collisions the \(v_n\)\([p_{\mathrm{T}}] \) correlation is sensitive to the distribution of energy deposition in the first stage of the collision. For a larger source a positive \(v_2\)\([p_{\mathrm{T}}] \) correlation is expected while for a compact source the negative correlation is obtained. Simultaneous measurements of \(v_n\)\([p_{\mathrm{T}}] \) correlations in small and large systems may help disentangle the role of initial conditions and subsequent dynamical QGP evolution in final-state particle distributions.

To measure the strength of the \(v_n\)\([p_{\mathrm{T}}] \) correlation, the Pearson correlation coefficient (PCC) R [25] is used where

$$\begin{aligned} R = \frac{\mathrm {cov}( v_{n}\{2\}^2, [p_{\mathrm{T}}] )}{\sqrt{ \mathrm {Var} ( v_{n}\{2\} ^2)} \sqrt{ \mathrm {Var} ([p_{\mathrm{T}}]) }}. \end{aligned}$$
(1)

The term \( v_{n}\{2\} ^2\) is the square of the nth-order flow harmonic obtained using the two-particle correlation method [26], \(\mathrm {cov}( v_{n}\{2\}^2, [p_{\mathrm{T}}] \)) is the covariance between \( v_{n}\{2\} ^2\) and \([p_{\mathrm{T}}] \), and \(\mathrm {Var} ( v_{n}\{2\} ^2)\) and \(\mathrm {Var} ([p_{\mathrm{T}}])\) are the variances of the \( v_{n}\{2\} ^2\) and \([p_{\mathrm{T}}] \) distributions, respectively. Experimentally, however, the finite event-by-event charged-particle track multiplicity results in an additional broadening of the \( v_{n}\{2\} ^2\) and \([p_{\mathrm{T}}]\) distributions due to statistical fluctuations. Thus, the values of the respective variances are increased, especially for \([p_{\mathrm{T}}]\). The magnitude of this broadening depends on the choice of kinematic region and on detector performance, making direct comparisons between experimental results and with theoretical calculations difficult. To overcome this problem, a modified correlation coefficient \(\rho \), less sensitive to the charged-particle multiplicity than R, was suggested in Ref. [25]. To reduce the auto-correlation effects and those due to the finite charged-particle multiplicity in an event, the variances of the \( v_{n}\{2\} ^2\) and \([p_{\mathrm{T}}]\) distributions are replaced by corresponding dynamical variables, which are more sensitive to intrinsic initial-state fluctuations. The variance of \( v_{n}\{2\} ^2\) is replaced by its dynamical counterpart [27]

$$\begin{aligned} \mathrm {Var}\big (v_{n}\{2\}^2\big )_{\mathrm{dyn}} = v_n\{2\}^4 - v_n\{4\}^4 = \langle \mathrm {corr} _{n}\{4\}\rangle - \langle \mathrm {corr} _{n}\{2\}\rangle ^2, \end{aligned}$$
(2)

where \(\mathrm {corr} _{n}\{2\}\) and \(\mathrm {corr} _{n}\{4\}\) are the two- and four-particle correlations [26] and where angular brackets denote that they are averaged over events. These correlations are described in detail in Sect. 4.

The variance of \([p_{\mathrm{T}}]\) is replaced by the dynamical \(p_{\mathrm{T}}\) fluctuation magnitude [28, 29] \(c_k\) defined as

$$\begin{aligned} c_k= \left\langle \frac{1}{N_\mathrm {pair}} \sum _i \sum _{j\ne i} ({p_{\mathrm{T}}}_{,i} - \langle [p_{\mathrm{T}}] \rangle ) ({p_{\mathrm{T}}}_{,j}-\langle [p_{\mathrm{T}}] \rangle ) \right\rangle \end{aligned}$$
(3)

where \( \langle [p_{\mathrm{T}}] \rangle \) is the average \([p_{\mathrm{T}}] \) over the all analysed events. The modified PCC \(\rho \) is thus defined as

$$\begin{aligned} \rho = \frac{\mathrm {cov}( v_{n}\{2\}^2, [p_{\mathrm{T}}] )}{ \sqrt{\mathrm {Var}\big (v_{n}\{2\}^2\big )_{\mathrm{dyn}} } \sqrt{c_k}}. \end{aligned}$$
(4)

It was demonstrated in Ref. [25] that the \(\rho \) coefficient calculated using realistic and finite multiplicities provides a reliable estimate of the true value of R found in the limit of infinite multiplicity, whereas the coefficient R, calculated using Eq. (1) for finite multiplicity underestimates the true value.

The ALICE experiment measured [20] that the charged-particle \(p_{\mathrm{T}}\) spectrum is correlated with the magnitude of the elliptic flow. It is measured to be harder in collisions with the higher second flow harmonics and softer in collisions where the elliptic flow is smaller. The magnitude of spectra modification is observed to increase with \(p_{\mathrm{T}}\), starting to be significant at around 1 \(\text {Ge} \text {V}\) and reaching a few percent at around 5 \(\text {Ge} \text {V}\). The modification is found to be most significant in the mid-central collisions, decreasing in the most central ones. The ALICE results suggest that the value of the correlation coefficient should be significant in mid-central and central collisions and that its magnitude and centrality dependence should be sensitive to the scale of intrinsic fluctuations of \(v_2\) and \(p_{\mathrm{T}}\). Including particles of higher \(p_{\mathrm{T}}\) in the measurement is expected to result in increased values of the \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \). The \([p_{\mathrm{T}} ]\) correlations with \(v_2\) in peripheral Pb+Pb collisions, \(v_3\) and \(v_4\) in wide centrality range as well as for the \(v_2\) in high multiplicity \(p\)+Pb are unexplored by measurements.

This paper reports on the first measurement of the \(\rho \) coefficient with the ATLAS detector in Pb+Pb and \(p\)+Pb collisions at a centre-of-mass energy per nucleon pair of 5.02 \(\text {Te} \text {V}\). The Pb+Pb data sample with a total integrated luminosity of 22 \(\mu \text {b}^{-1} \) was collected in 2015, and the \(p\)+Pb sample with \(28~\text{ nb }^{-1}\) in 2013.

This paper is organised as follows. Section 2 gives a brief description of the ATLAS detector. Details of the event selection and charged-particle reconstruction are provided in Sect. 3. Section 4 describes the analysis procedure for calculating the \(\rho \) coefficient. Systematic uncertainties are described in Sect. 5 and Appendix A. Results are presented in Sect. 6, followed by a summary in Sect. 7.

2 Experimental setup

The ATLAS experiment [30] at the LHC is a multipurpose particle detector with a forward–backward symmetric cylindrical geometry and a near \(4\pi \) solid angle coverage. The inner detector (ID) covers the pseudorapidityFootnote 1 range \(|\eta | < 2.5\) and is surrounded by a thin superconducting solenoid providing a 2 T axial magnetic field. The ID consists of silicon pixel, silicon microstrip (SCT), and straw tube tracking detectors. After the 2013 \(p\)+Pb run, an additional pixel silicon layer, the insertable B-layer [31, 31, 32], was installed prior to the 5.02 \(\text {Te} \text {V}\) Pb+Pb data-taking to attain more precise tracking. Lead/liquid-argon (LAr) sampling calorimeters provide electromagnetic (EM) energy measurements with high granularity. A steel/scintillator tile hadronic calorimeter covers the central pseudorapidity range (\(|\eta | < 1.7\)). The endcap and forward regions are instrumented with LAr calorimeters for EM and hadronic energy measurements up to \(|\eta | = 4.9\). The forward calorimeter (FCal) covers \(3.2<|\eta |< 4.9\) and is used for centrality estimation [10]. The minimum-bias trigger scintillators (MBTS) are located on each side of the detector at \(z = \pm 3.6~\hbox {m}\) and detect charged particles with \(2.07< |\eta | < 3.86\). The zero-degree calorimeter (ZDC), located in the LHC tunnel and covering \(|\eta |>8.3\), is used for triggering on collision events and pile-up event rejection. It is calibrated to resolve an individual neutron originating from the collision spectators.

A two-level trigger system selects events [33, 34]. The level-1 trigger is implemented in hardware and preselects up to \(10^5\) events per second for further decisions by the high-level trigger (HLT). The software-based HLT tuned for Pb+Pb collision data selects up to 1000 events per second for recording. This analysis primarily uses charged-particle tracks in the ID, but information from the central calorimeters and the ZDC is also used for triggering, event selection, and analysis.

3 Event and track selection

The Pb+Pb data in this analysis were selected using two mutually exclusive minimum-bias triggers. Events with semi-central and central collisions were selected if the scalar sum of transverse energy in the entire ATLAS calorimeter system exceeded 50 \(\text {Ge} \text {V}\). Peripheral events, i.e. those with large impact parameter of the colliding Pb nuclei, fail the 50 \(\text {Ge} \text {V}\) selection and were instead selected by requiring a deposition in the ZDC corresponding to at least one neutron and by requiring at least one track reconstructed in the HLT. Data in this analysis are required to come from periods when the entire detector was functioning normally. The events are required to have a reconstructed vertex within 100 mm of the nominal interaction point. The contribution from events containing more than one inelastic interaction (pile-up) is studied by exploiting correlations between the transverse energy measured in the FCal (\(\Sigma E_{\mathrm{T}}^{\mathrm{FCal}}\)) with the estimated number of neutrons in the ZDC, and with the number of tracks associated with a primary vertex [27, 35]. The distribution of \(\Sigma E_{\mathrm{T}}^{\mathrm{FCal}}\) and the distribution of the number of neutrons in events with more than one collision are broader than the corresponding distributions in events with only one collision. Pile-up events are suppressed by rejecting events with abnormally large values of either \(\Sigma E_{\mathrm{T}}^{\mathrm{FCal}}\) or the number of neutrons in the ZDC compared with the charged-particle multiplicity in the event. Approximately 0.2% of the events are rejected by these requirements.

The \(p\)+Pb data in this analysis were selected using minimum-bias triggers and high-multiplicity triggers (HMT). The minimum-bias trigger required signals in both sides of the MBTS system with a timing difference of less than 10 ns to eliminate non-collision backgrounds. The HMT required the total transverse energy in the calorimeter at level-one and the number of ID track candidates reconstructed in the HLT to be above predefined thresholds. Six combinations of thresholds were used to optimise data-taking during periods with different luminosities. Samples of events collected by these triggers are combined by applying event weights to reproduce the charged-particle multiplicity distribution of the minimum-bias trigger. Further details of the data selection are given in Refs. [22, 36]. The average pile-up probability in the \(p\)+Pb dataset is approximately 3% but can be significantly larger in high-multiplicity events. Events with more than one reconstructed vertex are removed from the sample. Similarly to the Pb+Pb dataset, to remove events where the two interaction vertices are too close to resolve as independent ones, the ZDC signal on the Pb fragmentation side is used. The distribution of the number of neutrons, which is broader in events with pile-up than that for the events without pile-up is exploited for that purpose [36]. The fraction of rejected events varies with the event activity and reaches a maximum of 10% for events with the highest multiplicities.

The analysis for both collision systems is performed in narrow bins of event activity defined by the charged-particle multiplicity \(N_{\mathrm{ch}}\) (described in Sect. 4), which estimates the collision centrality. In addition, the Pb+Pb results are presented as a function of collision centrality expressed by the average number of nucleons participating in the collision, \(N_{\mathrm{part}}\), to allow comparison with theoretical predictions [37]. The centrality is estimated from the \(\Sigma E_{\mathrm{T}}^{\mathrm{FCal}}\) distribution [6, 10] using the Glauber model [38]. The number of events passing the selection requirements is \(1.3\times 10^8\) for Pb+Pb within the 0–80% centrality interval. For the \(p\)+Pb system, about \(0.64 \times 10^8\) events enter the analysis.

The charged-particle tracks reconstructed in the ID are required to satisfy selection criteria in order to suppress the contribution of incorrectly reconstructed tracks and secondary products of particle decays. The selection criteria include the requirement that the number of hits in the pixel and SCT detectors should be greater than two and eight, respectively, for the Pb+Pb data and greater than one and six for the \(p\)+Pb data. The track impact parameters relative to the collision vertex in the transverse direction, \(|d_0|\), and longitudinal direction, \(|z_0\sin \theta |\), are required to be less than 1 mm for tracks in the Pb+Pb data sample and less than 1.5 mm in the \(p\)+Pb sample. In addition, in \(p\)+Pb collisions, the track impact parameter significances must satisfy \(|d_0/\sigma _{d_0} |<\) 3 and \(|z_0 \sin \theta /\sigma _{z} |<\) 3, where \(\sigma _{d_0}\) and \(\sigma _{z}\) are the uncertainties in \(d_0\) and \(z_0 \sin \theta \) determined from the covariance matrix of the track fit. The different selection criteria for Pb+Pb and \(p\)+Pb optimise the performance of the track reconstruction in differing running conditions.

Corrections needed due to track reconstruction effects are evaluated using 4\(\times 10^6\) Pb+Pb and \(10^7\) \(p\)+Pb minimum-bias Monte Carlo (MC) events generated by the HIJING v1.38b [39] event generator. After the generation, an azimuthal flow is implemented using the afterburner technique [40], and the \(p_{\mathrm{T}}\) spectrum is reweighted to match the data. Generated events were simulated in the detector by the Geant 4-based [41] ATLAS detector simulation programs [42] and reconstructed using the same procedures and detector conditions as the data. Track reconstruction corrections are applied to each selected track using weights to account for the tracking efficiency \(\epsilon \) and the fake-track fraction f. The efficiency is defined as the fraction of primary MC charged particles that are matched to reconstructed tracks, and f is the fraction of tracks that are not matched to primary MC particles or are produced from random combinations of hits in the ID. A similar analysis procedure is described in Refs. [10, 16]. The fake-track fraction and tracking efficiency are determined as functions of the track \(p_{\mathrm{T}}\) and \(\eta \) and of the track multiplicity in the event. Tracks included in the analysis are weighted with the factor \((1-f)/\epsilon \). An additional multiplicative weight evaluated from data is applied to the data to correct for detector non-uniformity in the azimuthal angle. These weights are obtained by requiring the tracks to be distributed uniformly in azimuth in all pseudorapidity slices of width 0.1.

In the Pb+Pb data, the contribution of fake tracks is largest in central collisions at the lowest analysed track \(p_{\mathrm{T}}\) of 0.5 \(\text {Ge} \text {V}\) and at the largest \(|\eta |\), reaching up to 20%. The fake-track rate is below 1% for tracks with \(p_{\mathrm{T}}\) above 2 \(\text {Ge} \text {V}\) and \(|\eta |<1.5\). The tracking efficiency depends weakly on centrality, and in the most central events it is about 3% less than in more peripheral events. The efficiency increases with the track \(p_{\mathrm{T}}\) from about 50% at the lowest analysed \(p_{\mathrm{T}}\) to 70% above 2 \(\text {Ge} \text {V}\). It is highest at mid-rapidity and drops by about 15% for \(|\eta |>1\). For \(p\)+Pb collisions, with \(p_{\mathrm{T}}\) increasing from 0.3 to 1 \(\text {Ge} \text {V}\) the efficiency increases from about 75% (60%) to 83% (70%) at \(\eta \approx \) 0 (\(|\eta | > 2\)). The \(p\)+Pb tracking efficiency is independent of the event’s multiplicity for \(N_{\mathrm{ch}} \ge 10\), i.e. in the multiplicity range used in the analysis. The fake rate in \(p\)+Pb collisions is very low, below 1% (3%) at \(\eta \approx \) 0 (\(|\eta | > 2\)).

4 Correlation coefficient \(\rho \)

In each event, charged-particle tracks are grouped into three regions of subevents based on their pseudorapidity: region A with \(-2.5<\eta <-0.75\), central region B with \(|\eta |<0.5\) and region C with \(0.75< \eta < 2.5\). The \(v_{n}^2\) for the n = 2–4 harmonics are calculated by correlating charged-particle tracks from subevents A and C, which are separated in pseudorapidity to suppress non-flow contributions. Tracks in central region B are used to obtain the mean value of the charged-particle transverse momentum in the event, \([p_{\mathrm{T}} ]\), defined as

$$\begin{aligned}{}[p_{\mathrm{T}}] = \frac{1}{\sum _{b}{w_b}} \sum _{b } w_b p_{\mathrm {T}b} \end{aligned}$$

where the summation is over tracks in region B, labelled by index b. The variable \(c_k\) (Eq. (3)) is also calculated using tracks from region B. Here, and in following formulas, the weights w include the fake-track fraction, efficiency, and azimuthal non-uniformity corrections, as discussed in Sect. 3.

The covariance term from the numerator of Eq. (4) is defined as

$$\begin{aligned}&\mathrm {cov} ( v_{n}\{2\}^2, [p_{\mathrm{T}}] ) \nonumber \\&\quad = \mathrm {Re}\left( \left\langle \frac{1}{\sum _{a,c}{w_a w_c}} \sum _{a, c } w_a w_c \mathrm {e}^{\mathrm{i}n\phi _a - \mathrm {i}n\phi _c} ([p_{\mathrm{T}}]- \langle [p_{\mathrm{T}}] \rangle )\right\rangle \right) , \end{aligned}$$
(5)

where \(\phi \) is the azimuthal angle and indices a and c span the tracks in regions A and C, respectively.

The two- and four-particle correlations used to define the dynamical variance in Eq. (2), which enters the denominator of Eq. (4), are calculated as in Ref. [26]

$$\begin{aligned} \langle \mathrm {corr} _{n}\{2\} \rangle= & {} \mathrm {Re} \left( \left\langle \frac{1}{\sum _{a,c}{w_a w_c}} \sum _{a, c } w_a w_c \mathrm {e}^{\mathrm{i}n\phi _a - \mathrm {i}n\phi _c} \right\rangle \right) \nonumber \\= & {} \mathrm {Re} \left( \langle q_{n,a} q_{n,c}^{*} \rangle \right) \end{aligned}$$
(6)

where the \(q_a\) and \(q_c\) are the complex flow vectors of subevent A and subevent C, respectively, and the asterisk denotes the complex conjugate. The flow vectors are

$$\begin{aligned} q_{n,a}=\frac{1}{\sum _{a}{w_a}} \sum _{a}{w_a \mathrm {e}^{\mathrm{i}n\phi _a}}\quad \mathrm {and} \quad q_{n,c}=\frac{1}{\sum _{c}{w_c}} \sum _{c}{w_c \mathrm {e}^{\mathrm{i}n\phi _c}}. \end{aligned}$$

The four-particle correlation is obtained from the expression

$$\begin{aligned} \langle \mathrm {corr} _{n}\{4\} \rangle = \mathrm {Re}\left( \left\langle \frac{(Q_{n,a}^2-Q_{2n,a})(Q_{n,c}^2-Q_{2n,c})^{*}}{S_a S_c} \right\rangle \right) , \end{aligned}$$
(7)

where for subevent A

$$\begin{aligned}&Q_{n,a}=\sum _{a}{w_a \mathrm {e}^{\mathrm{i}n\phi _a}}, \qquad Q_{2n,a}= \sum _{a}{w_a^2 \mathrm {e}^{\mathrm{i}2n\phi _a}}, \\&S_a= \left( \sum _{a}{w_a} \right) ^2- \sum _{a}{w_a}^2, \end{aligned}$$

and similarly for subevent C. Equation (7) represents the sum \(\sum {\mathrm {e}^{\mathrm{i}n(\phi _1^a+\phi _2^a - \phi _3^c - \phi _4^c)}}\) over all particles from subevents A and C normalised by the number of quadruplets without auto-correlations in each subevent.

The second factor in the denominator of Eq. (4), the mean \(p_{\mathrm{T}}\) fluctuation in the event class \(c_k\), is defined by Eq. (3) and in this analysis it is calculated as

$$\begin{aligned} c_k= & {} \left\langle \frac{1}{(\sum _b w_b)^2 - \sum _b w_b^2}\right. \\&\qquad \left. \sum _b \sum _{ b'\ne b} w_b (p_{\mathrm {T},b} - \langle [p_{\mathrm{T}}] \rangle ) w_{b'} (p_{\mathrm {T},b'}-\langle [p_{\mathrm{T}}] \rangle ) \right\rangle . \end{aligned}$$

The summation indices b and \(b'\) run over all charged particles in region B.

The correlation coefficient expressed by Eq. (4) is evaluated for the range \(0.5<p_{\mathrm{T}} < 2~\text {Ge} \text {V}\) in Pb+Pb collisions and \(0.3<p_{\mathrm{T}} < 2~\text {Ge} \text {V}\) in \(p\)+Pb collisions. These intervals, called ‘main’, contain a large number of soft particles and constitute the main result of the analysis which can be compared with hydrodynamical models. For each system, two additional \(p_{\mathrm{T}}\) ranges are considered: \(0.5<p_{\mathrm{T}} < 5~\text {Ge} \text {V}\) and \(1<p_{\mathrm{T}} < 2~\text {Ge} \text {V}\) in the analysis of Pb+Pb collisions, and \(0.3<p_{\mathrm{T}} < 5~\text {Ge} \text {V}\) and \(0.5<p_{\mathrm{T}} < 2~\text {Ge} \text {V}\) in \(p\)+Pb collisions. These ranges facilitate the study of the sensitivity of \( \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]) \) to the high \(p_{\mathrm{T}}\) part of the particle spectrum and to the lower charged-particle multiplicity from the higher minimum \(p_{\mathrm{T}}\) value. The charged-particle \(p_{\mathrm{T}}\) range \(0.5<p_{\mathrm{T}} < 2~\text {Ge} \text {V}\) is common to both systems and can be used to compare the \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \) results from Pb+Pb and \(p\)+Pb collisions.

The quantities of interest, i.e. \(\mathrm {cov}( v_{n}\{2\}^2, [p_{\mathrm{T}}] ) \), \(\mathrm {Var}\big (v_{n}\{2\}^2\big )_{\mathrm{dyn}} \), \(c_k\), and \( \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}] \)), are determined in bins of reconstructed track multiplicity \( M_{\mathrm{AC}}\) measured in the combination of regions A and C. This is done to avoid a negative correlation between the multiplicity in subevents A+C and B that occurs if the analysis is binned in multiplicity in the entire ID. Narrow \( M_{\mathrm{AC}}\) bins are also chosen due to the sensitivity to multiplicity fluctuations of the multi-particle correlations that are used to obtain the \(\mathrm {Var}\big (v_{n}\{2\}^2\big )_{\mathrm{dyn}}\)  [27]. The events are grouped in fine bins with a width of ten in \( M_{\mathrm{AC}}\) for \(0.5< p_{\mathrm{T}} < 5~\text {Ge} \text {V}\) in the Pb+Pb analysis and \(0.3< p_{\mathrm{T}} < 5~\text {Ge} \text {V}\) in the \(p\)+Pb analysis. It was cross-checked that the variables of interest obtained with a finer binning in \( M_{\mathrm{AC}}\) are consistent with the measurement with the nominal binning.

To enable comparisons with the theoretical predictions and with future experimental results, measurements obtained in \( M_{\mathrm{AC}}\) are presented as a function of the ATLAS ID multiplicity \(N_{\mathrm{ch}}\) of \(0.5< p_{\mathrm{T}} < 5~\text {Ge} \text {V}\) and \(|\eta |<2.5\). They are projected from the \( M_{\mathrm{AC}}\) values taking into account tracking efficiency and fake-track production as described in the previous section. A similar analysis procedure is described in Ref. [22]. For the \(N_{\mathrm{part}}\) dependencies in the Pb+Pb system, the results measured in \( M_{\mathrm{AC}}\) multiplicity intervals are averaged, with weights equal to the probabilities to find any given \( M_{\mathrm{AC}}\) value in the centrality intervals.

The formulation of the modified PCC \( \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]) \) requires that there should be at least two tracks in each region (A, B, and C). Further, \(\mathrm {Var}\big (v_{n}\{2\}^2\big )_{\mathrm{dyn}}\) calculated according to Eq. (6) can be negative at low multiplicities due to statistical fluctuations, which renders Eq. (4) invalid because of the \(\surd \mathrm {Var}\big (v_{n}\{2\}^2\big )_{\mathrm{dyn}} \) term. For each \( M_{\mathrm{AC}}\) bin, \(p_{\mathrm{T}}\) interval, and harmonic, a criterion is applied that \(\mathrm {Var}\big (v_{n}\{2\}^2\big )_{\mathrm{dyn}}\) needs to be positive at a level of at least one standard deviation of its statistical uncertainty. Results presented as a function of \(N_{\mathrm{ch}}\) are produced only for those \( M_{\mathrm{AC}}\) intervals. For the \(N_{\mathrm{part}}\) dependencies in the Pb+Pb system, it is additionally required for each centrality interval that the fraction of rejected events due to this criterion does not exceed 1%.

5 Systematic uncertainties

The systematic uncertainty is estimated by varying individual aspects of the analysis. The systematic uncertainties for the main \(p_{\mathrm{T}}\) interval are discussed for each collision system. Systematic uncertainties for the other \(p_{\mathrm{T}}\) intervals behave consistently with the ones for the main \(p_{\mathrm{T}}\) interval. Since the modified PCC \( \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]) \) is a ratio of quantities which are calculated using tracks, many variations largely cancel out and the resulting systematic uncertainties are small. To suppress the statistical fluctuations and to get more robust estimation of systematic uncertainties, they are averaged over several, wide ranges of the charged-particle multiplicity. For each uncertainty source and for each measurement point, the maximum variation from the baseline measurement is used. The total resulting uncertainty is the sum of the individual contributions combined in quadrature. The following sources of systematic uncertainties are considered.

Track selection The tracking performance has a relatively small impact on \( v_{n}\{2\} \), but it directly affects the \([p_{\mathrm{T}}]\) and \(c_k\)via the admixture of the fake tracks, especially at low \(p_{\mathrm{T}}\). To assess the impact on \( \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]) \), the measurement is repeated with tracks selected with looser and tighter track quality criteria, thus increasing and decreasing the fake-track rate, respectively. The weights used in the evaluation of measured quantities take the modified selection into account. The loose track selection in the Pb+Pb analysis relaxes requirements on the number of pixel and SCT hits to at least one and six, respectively. Additionally, the requirements on the transverse and longitudinal impact parameters of the track are relaxed to 1.5 mm. The tighter selection in the Pb+Pb analysis tightens the requirement on the transverse and longitudinal impact parameters of the track to 0.5 mm. For the \(p\)+Pb analysis, the loose selection relaxes the requirements on the transverse and longitudinal impact parameters of the track to 2 mm and on the impact parameter significances to less than 4. In the tight selection, the impact parameter values and their significances must be less than 1 mm and 2, respectively. For each of the two track selections the absolute difference is calculated with respect to the baseline measurement: \(| \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]) ^{\mathrm{base}}- \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]) ^{\mathrm{loose}}|\) or \(| \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]) ^{\mathrm{base}}- \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]) ^{\mathrm{tight}}|\). The largest difference is taken as a systematic uncertainty.

Detector material Since the tracks that are used in the calculation of \( \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]) \) are weighted by the inverse of the tracking efficiency, a bias in its estimation due to inaccurate modelling of the material in the detector may change the balance between low- and high-\(p_{\mathrm{T}}\) tracks in the sums. Based on simulations, the estimated uncertainty in the detector description is obtained [43, 44]. The resulting \(p_{\mathrm{T}}\)- and \(\eta \)-dependent uncertainties in the track efficiency of up to 4% are used to determine the systematic uncertainty.

Tracking azimuthal uniformity In this analysis, the weighting factors w correct for any non-uniformity in the azimuthal angle distribution of reconstructed tracks. The weights are obtained from the data by requiring azimuthal uniformity over the two-dimensional distribution of reconstructed tracks in the \(\eta \)\(\phi \) plane. The effect of that correction on the result is conservatively estimated by comparing the baseline measurement and the measurement obtained without applying this weight. The uncertainty is small, and it envelopes potential effects of imperfections in the weighting factors determination, including their dependence on the transverse momentum, collision centrality, run-by-run differences, on dead module maps or the vertex position.

Fig. 1
figure 1

The systematic uncertainty of \( \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}] \)) as a function of \(N_{\mathrm{ch}}\) measured with tracks from main \(p_{\mathrm{T}}\) intervals for each collision system for the a second, b third, and c fourth harmonics in Pb+Pb collisions, and for d \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}] \)) in \(p\)+Pb collisions. The total uncertainty is also shown

Residual pile-up events The selection criteria discussed in Sect. 3 suppress the fraction of pile-up events accepted for analysis to almost zero in central Pb+Pb collisions. To estimate the systematic uncertainty related to pile-up, the measurement is conservatively repeated without this event rejection, resulting in at most a 1% difference in the most central Pb+Pb events for the \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}] \)) coefficient. The \(p\)+Pb data were taken with higher pile-up than the Pb+Pb data. To estimate the impact of contamination by residual pile-up events, \(p\)+Pb results were obtained with only the vertex criteria applied. The variation covers the estimated residual pile-up fraction in events of the highest track multiplicity [36].

Centrality selection The minimum-bias trigger is fully efficient for the 0–85% centrality interval. However, the total fraction of inelastic Pb+Pb events selected is known only to 1% accuracy due to trigger inefficiency and possible sample contamination in more peripheral interactions. The centrality is estimated using the \(\Sigma E_{\mathrm{T}}^{\mathrm{FCal}}\) distribution [6, 10] and the Glauber model [38] to obtain the mapping from the observed \(\Sigma E_{\mathrm{T}}^{\mathrm{FCal}}\) to the number of nucleons participating in the collision, \(N_{\mathrm{part}}\). The modified PCC uncertainty is evaluated by repeating the analysis with the altered centrality selections on the \(\Sigma E_{\mathrm{T}}^{\mathrm{FCal}}\) distribution, which results in \(\pm 1\%\) uncertainty in the total fraction of inelastic Pb+Pb events. The centrality selection contributes mainly to uncertainties for peripheral collisions.

Figure 1 shows the magnitude of the systematic uncertainties \(\updelta \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]) \) for \(n = 2-4\) in Pb+Pb collisions as a function of \(N_{\mathrm{ch}}\). In Pb+Pb collisions, the systematic uncertainty of the measured correlation coefficients across different order harmonics and centralities is not dominated by a single source. One of the largest uncertainties comes from restoring the azimuthal uniformity, and dominates for the second order harmonic in the most central collisions and for the third and fourth order harmonics almost over the full centrality range. A sizeable contribution to the uncertainty for all three harmonics is due to the track selection. The impact of the detector material is rather small except for a significant contribution for the forth order harmonic in the most central events. The residual pile-up in Pb+Pb collisions gives a negligible contribution. Figure 1d shows systematic uncertainties for \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \) coefficients in \(p\)+Pb collisions for the main interval of \(0.3< p_{\mathrm{T}} < 2~\text {Ge} \text {V}\) as a function of event activity. In \(p\)+Pb interactions the largest uncertainty in the most active collisions (\(N_{\mathrm{ch}} > 150\)) originates from pile-up. The track selection is a source of sizeable uncertainty for this collision system, while the azimuthal uniformity correction procedure and the detector material have a small impact.

Details on the contributions to systematic uncertainties from different sources of \(c_k\), \(\mathrm {Var}\big (v_{n}\{2\}^2\big )_{\mathrm{dyn}}\) and \(\mathrm {cov}( v_{n}\{2\}^2, [p_{\mathrm{T}}] \)) are included in the Appendix.

6 Results

6.1 The constituents of the modified PCC

The constituents of the modified PCC, \(c_k\) , \(\mathrm {Var}\big (v_{n}\{2\}^2\big )_{\mathrm{dyn}} \) and \(\mathrm {cov}( v_{n}\{2\}^2, [p_{\mathrm{T}}] ) \) and are combined, using Eq. (4), to obtain \(\rho \). Figure 2 shows the dynamical \(p_{\mathrm{T}}\) fluctuation coefficient \(c_k\) as a function of charged-particle multiplicity in Pb+Pb and \(p\)+Pb collision systems for tracks in three different \(p_{\mathrm{T}}\) intervals. A strong decrease of \(c_k\) with increasing \(N_{\mathrm{ch}}\) is observed in all measured results. A similar decrease was seen for \(c_k\) in Au+Au and Pb+Pb data at lower centre-of-mass energies [28, 29], evaluated for lower \(p_{\mathrm{T}}\) range, \(0.15<p_{\mathrm{T}} < 2~\text {Ge} \text {V}\), not accessible with the ATLAS detector. For the same \(N_{\mathrm{ch}}\), the \(c_k\) values differ by an order of magnitude for different \(p_{\mathrm{T}}\) ranges of tracks used in the analysis. For the intervals with the same lower \(p_{\mathrm{T}}\) limit, the \(c_k\) values are higher for the interval with the larger upper \(p_{\mathrm{T}}\) limit.

Fig. 2
figure 2

The variable \(c_k\) for three \(p_{\mathrm{T}}\) ranges as a function of the charged-particle multiplicity \(N_{\mathrm{ch}}\) of a Pb+Pb and b \(p\)+Pb collisions. The statistical and systematic uncertainties are shown as vertical error bars (smaller than symbols) and boxes, respectively

Figure 3 shows \(\mathrm {Var}\big (v_{n}\{2\}^2\big )_{\mathrm{dyn}} \) for \(n = 2-4\) as function of \(N_{\mathrm{ch}}\) for Pb+Pb collisions. For low multiplicities, \(\mathrm {Var}\big (v_{n}\{2\}^2\big )_{\mathrm{dyn}} \) increases with increasing \(N_{\mathrm{ch}} \), reaching a maximum at \(N_{\mathrm{ch}}\) of approximately 500 (1000) for \(n = 2\) (\(n = 3\)), respectively. At higher \(N_{\mathrm{ch}}\) values the variances decrease with multiplicity. The dynamical variance for \(n=4\), measured for \(N_{\mathrm{ch}} \gtrapprox 500\), decreases with increasing \(N_{\mathrm{ch}}\). The ordering \(\mathrm {Var}\big (v_{2}\{2\}^2\big )_{\mathrm{dyn}}> \mathrm {Var}\big (v_{3}\{2\}^2\big )_{\mathrm{dyn}} > \mathrm {Var}\big (v_{4}\{2\}^2\big )_{\mathrm{dyn}} \) and the multiplicity dependence of \(\mathrm {Var}\big (v_{n}\{2\}^2\big )_{\mathrm{dyn}} \) are similar to the ordering and centrality dependence of \( v_{n}\{2\} \) measured by ATLAS [10]. Also shown in Fig. 3 is \(\mathrm {Var}\big (v_{2}\{2\}^2\big )_{\mathrm{dyn}} \) for \(p\)+Pb collisions as a function of \(N_{\mathrm{ch}}\). The dependence is monotonic, similarly to \( v_{2}\{2\} \) [45]. In both collision systems and for all harmonics, the same ordering of \(\mathrm {Var}\big (v_{n}\{2\}^2\big )_{\mathrm{dyn}} \) depending on the \(p_{\mathrm{T}}\) interval is observed. The largest variances are observed for the \(p_{\mathrm{T}}\) intervals with an increased lower limit. This is expected as the \( v_{n}\{2\} \) value increases strongly with \(p_{\mathrm{T}}\) below 3 \(\text {Ge} \text {V}\) [10]. Additionally, the interval in which the upper limit on \(p_{\mathrm{T}}\) is set to 5 \(\text {Ge} \text {V}\) integrates the region with the highest values of \( v_{n}\{2\} \) (which occur around 3 \(\text {Ge} \text {V}\)) and thus the values of the variance are expected to be larger than that for the main \(p_{\mathrm{T}}\) range.

Fig. 3
figure 3

The variance \(\mathrm {Var}\big (v_{n}\{2\}^2\big )_{\mathrm{dyn}} \) for \(n = 2-4\) for ac Pb+Pb collisions and \(\mathrm {Var}\big (v_{2}\{2\}^2\big )_{\mathrm{dyn}} \) for d \(p\)+Pb collisions for the three \(p_{\mathrm{T}}\) intervals as a function of charged-particle multiplicity \(N_{\mathrm{ch}}\). The statistical and systematic uncertainties are shown as vertical error bars and boxes, respectively

In Fig. 4, the covariances \(\mathrm {cov}( v_{n}\{2\}^2, [p_{\mathrm{T}}] ) \) are shown for the 2nd-, 3rd-, and 4th-order harmonics in Pb+Pb collisions and for the second-order harmonics in \(p\)+Pb collisions. They are presented as a function of \(N_{\mathrm{ch}}\) for three \(p_{\mathrm{T}}\) intervals. Significant positive correlations between \( v_{n}\{2\} \) and \([p_{\mathrm{T}}]\) are observed in the Pb+Pb events. The measured covariances depend on the charged-particle multiplicity and the \(p_{\mathrm{T}}\) range of the charged particles. In Pb+Pb collisions, a strong dependence on the multiplicity is observed for \(n=2\) and 4. The \(\mathrm {cov}( v_{3}\{2\}^2, [p_{\mathrm{T}}] ) \) depends only weakly on \(N_{\mathrm{ch}}\). A negative \(\mathrm {cov}( v_{2}\{2\}^2, [p_{\mathrm{T}}] ) \) is measured at multiplicities \(N_{\mathrm{ch}} < 200\) and a negative \(\mathrm {cov}( v_{3}\{2\}^2, [p_{\mathrm{T}}] ) \) for 1 \(<p_{\mathrm{T}}\ < 2~\text {Ge} \text {V}\) below \(N_{\mathrm{ch}} < 1800\). The covariances \(\mathrm {cov}( v_{2}\{2\}^2, [p_{\mathrm{T}}] \)) in \(p\)+Pb events are negative in the entire measured \(N_{\mathrm{ch}}\) range and show weak \(N_{\mathrm{ch}}\) dependence. Unlike in Pb+Pb events, the \(\mathrm {cov}( v_{2}\{2\}^2, [p_{\mathrm{T}}] \)) in \(p\)+Pb events have similar magnitudes for different \(p_{\mathrm{T}}\) intervals.

Fig. 4
figure 4

The covariance \(\mathrm {cov}( v_{n}\{2\}^2, [p_{\mathrm{T}}] ) \) for \(n = 2-4\) in ac Pb+Pb collisions and \(\mathrm {cov}( v_{2}\{2\}^2, [p_{\mathrm{T}}] ) \) in d \(p\)+Pb collisions for three \(p_{\mathrm{T}}\) ranges as a function of the charged-particle multiplicity \(N_{\mathrm{ch}}\). The statistical and systematic uncertainties are shown as vertical error bars and boxes, respectively

6.2 The modified PCC

Fig. 5
figure 5

The PCC \( \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]))\) for \(n = 2-4\) in ac Pb+Pb collisions and d \(p\)+Pb collisions as a function of the charged-particle multiplicity \(N_{\mathrm{ch}}\) for three \(p_{\mathrm{T}}\) ranges. The statistical and systematic uncertainties are shown as vertical error bars and boxes, respectively

The modified PCC \( \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]) \) for \(n = 2-4\) in Pb+Pb collisions and for \(n=2\) in \(p\)+Pb collisions is shown in Fig. 5. In Pb+Pb collisions, the behaviour of \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}] \)) is similar for all \(p_{\mathrm{T}}\) intervals. It starts at negative values for \(N_{\mathrm{ch}} < 200\) and rapidly increases with multiplicity up to \(\sim 1500\) particles where the increase slows down and reaches the maximum at \(N_{\mathrm{ch}} \approx 4500\) of 0.24–0.3, depending on the \(p_{\mathrm{T}}\) interval. At even higher \(N_{\mathrm{ch}}\), the \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}] \)) value decreases rapidly. The significant correlation observed for mid-central events suggests a connection between anisotropic and radial [46] flows which might be attributed to stronger hydrodynamic response (larger pressure gradients) to the large initial-state eccentricities [47]. The modified PCC multiplicity dependence could reflect a balance between stronger radial flow observed in central collision and the larger initial eccentricity seen in peripheral interactions. The decrease observed in central collisions, for \(N_{\mathrm{ch}} \gtrsim 5000\), might be related to the increased role of initial-state fluctuations in anisotropic flow [27]. However, a complete understanding of this effect would require a more precise modelling of heavy ion collisions. The correlation coefficients calculated with the upper \(p_{\mathrm{T}}\) limit of 2 \(\text {Ge} \text {V}\) are 10–20% smaller than the values obtained with a \(p_{\mathrm{T}}\) limit of 5 \(\text {Ge} \text {V}\). The correlation coefficient \( \rho (v_{3}\{2\}^2, [p_{\mathrm{T}}]) \) is evaluated in Pb+Pb collisions for the same three \(p_{\mathrm{T}}\) ranges. The magnitudes measured for \( \rho (v_{3}\{2\}^2, [p_{\mathrm{T}}] \)) are significantly smaller than those measured for \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}] \)) and similar to the magnitudes of \( \rho (v_{4}\{2\}^2, [p_{\mathrm{T}}] \)). All three curves increase with \(N_{\mathrm{ch}}\) in the range of \( 1000< N_{\mathrm{ch}} < 5000 \). At low values of \(N_{\mathrm{ch}} \), a flattening of the trend can be noticed. In the most central collisions, a breakdown of the rise is seen, similarly to the \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \). Above \(N_{\mathrm{ch}} \sim 1500\), the curves for the two intervals with the same maximum \(p_{\mathrm{T}}\) are consistent with each other and are below the curve for the interval which uses tracks with \(p_{\mathrm{T}}\) up to 5 \(\text {Ge} \text {V}\). The largest values of \( \rho (v_{4}\{2\}^2, [p_{\mathrm{T}}]) \) are observed at \(N_{\mathrm{ch}} \approx 1000\). For high \(N_{\mathrm{ch}}\), \( \rho (v_{4}\{2\}^2, [p_{\mathrm{T}}]) \) decreases with \(N_{\mathrm{ch}}\) up to about \(N_{\mathrm{ch}} \approx 4000\) and rises slowly at higher values. The trends obtained for \(p_{\mathrm{T}}\) intervals with the same minimum value are consistent above \(N_{\mathrm{ch}} \sim 1500\) as is the case for \( \rho (v_{3}\{2\}^2, [p_{\mathrm{T}}]) \). The decrease for \(N_{\mathrm{ch}} < 4000\) might be due to a contribution to \(v_{4}\) from a non-linear term containing \(v_{2}^2\), decreasing with increasing centrality [13]. However, a theoretical modelling of the initial state and its subsequent evolution would be required to support this interpretation. Similarly to the \( \rho (v_{3}\{2\}^2, [p_{\mathrm{T}}] \)), the \( \rho (v_{4}\{2\}^2, [p_{\mathrm{T}}] \)) correlations measured with the larger upper \(p_{\mathrm{T}}\) limit have larger magnitudes. The results for the larger upper \(p_{\mathrm{T}}\) limit show the sensitivity of the \( \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]) \) coefficients to the high \(p_{\mathrm{T}}\) part of the particle spectrum contaminated with non-flow correlations from jets. On the other hand, the correlations measured for the intervals with fixed upper \(p_{\mathrm{T}}\) limit (2 \(\text {Ge} \text {V}\)) and varied lower \(p_{\mathrm{T}}\) limits are similar, demonstrating insensitivity of the modified PCC coefficients to a significant change of the event charged-particle multiplicity as expected [25]. The fourth-order correlations are weaker than those for the second-order flow harmonic and for \(N_{\mathrm{ch}} >4000\) are comparable to \( \rho (v_{3}\{2\}^2, [p_{\mathrm{T}}]) \). The results for all harmonics indicate a change in the trend in events with high \(N_{\mathrm{ch}}\) around 4500, which suggests a change in the nature of the correlations in those events [47].

In \(p\)+Pb collisions, \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \) exhibits much weaker \(N_{\mathrm{ch}}\) dependence than that in Pb+Pb collisions. For the main \(p_{\mathrm{T}}\) interval, the modified PCC assumes a negative value of approximately \(-0.1\) and is almost constant within uncertainties. Values for different lower \(p_{\mathrm{T}}\) limits are similar, and the \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \) magnitudes for the larger upper \(p_{\mathrm{T}}\) limit are smaller. The magnitude (and sign) of the modified PCC in \(p\)+Pb collisions is expected to be related to the distribution of the energy deposition in the initial state, as predicted by the hydrodynamic model [25]. In hydrodynamics, in \(p\)+Pb collision, for small sources a higher initial pressure gradients and smaller eccentricities are expected to be generated. This mechanism could lead to the negative correlation of the final state observables, this is the mean transverse momentum and higher order flow harmonics. Thus, the negative value of the modified PCC for \( v_{2}\{2\} \) in \(p\)+Pb and peripheral Pb+Pb that is measured should provide valuable constraints for models describing the collectivity in small systems.

6.3 Comparison of \(p\)+Pb and Pb+Pb results

Fig. 6
figure 6

Comparison of a \(c_k\), b \(\mathrm {Var}\big (v_{2}\{2\}^2\big )_{\mathrm{dyn}} \), c \(\mathrm {cov}( v_{2}\{2\}^2, [p_{\mathrm{T}}] ) \), and the d \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \) for the range \(0.5< p_{\mathrm{T}} < 2~\text {Ge} \text {V}\) as a function of the charged-particle multiplicity \(N_{\mathrm{ch}}\). The statistical and systematic uncertainties are shown as vertical error bars and boxes, respectively

Figure 6 shows a comparison of \(p\)+Pb and Pb+Pb results shown in Figs. 2, 3, 4, 5 for the common \(p_{\mathrm{T}}\) interval of \(0.5< p_{\mathrm{T}} < 2~\text {Ge} \text {V}\). The values of the \(c_k\) (Fig. 6a) are similar for \(p\)+Pb and Pb+Pb collisions in this \(p_{\mathrm{T}}\) interval, while the behaviour of the dynamical variance \(\mathrm {Var}\big (v_{2}\{2\}^2\big )_{\mathrm{dyn}} \) (Fig. 6b) is very different due to the different initial eccentricities in the overlap regions in Pb+Pb and \(p\)+Pb collisions. Only a small rise with the multiplicity is observed for \(p\)+Pb collisions, which is in agreement with a slow increase of \( v_{2}\{2\} \) with growing event activity [22, 36, 45]. For \(N_{\mathrm{ch}} \approx 50\), the dynamical variances are comparable between Pb+Pb and \(p\)+Pb collisions. The \(N_{\mathrm{ch}}\) dependence of \(\mathrm {cov}( v_{2}\{2\}^2, [p_{\mathrm{T}}] \)) is significantly different for Pb+Pb and \(p\)+Pb collisions. A steady rise from negative to positive values with \(N_{\mathrm{ch}}\) is observed for peripheral Pb+Pb collisions, and approximately constant values are obtained for \(p\)+Pb collisions. The \(N_{\mathrm{ch}}\) dependence of \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \) is different for the two collision systems. Much weaker \(N_{\mathrm{ch}}\) dependence of modified PCC is observed in \(p\)+Pb collisions compared to Pb+Pb collisions. For \(N_{\mathrm{ch}}\ <100 \) the values of \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \) are consistent between Pb+Pb and \(p\)+Pb collisions. The negative \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \) coefficients for the small systems in \(p\)+Pb and Pb+Pb collisions may suggest a more compact source model [25]. The comparison of the systems underlines the importance of the initial stage in the correlations described by the \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \) coefficient. The theoretical predictions for midcentral and central Pb+Pb collisions suggests that for a large system an increase of the mean transverse momentum indicates a stronger transverse flow and a stronger collective response to the initial geometry of the source, characterized by the positive value of the modified PCC.

6.4 Comparison to theoretical predictions

To compare the Pb+Pb results with a theoretical prediction in Ref. [25], the \( \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]) \) coefficients for \(0.5< p_{\mathrm{T}} < 2~\text {Ge} \text {V}\) are obtained as a function of centrality intervals expressed by \(N_{\mathrm{part}} \) using the procedure described in Sect. 4. Figure 7 shows the \(N_{\mathrm{part}} \) dependence of \( \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]) \) for \(n = 2-4\) in Pb+Pb collisions. It resembles the trends observed in Fig. 5, which show the modified PCC as a function of \(N_{\mathrm{ch}} \), a measure of event activity. The theoretical predictions of the \( \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]) \) coefficient are based on a model in which the initial conditions were generated with nucleon positions by a MC Glauber model [48]. These initial conditions are then evolved using the pressure-driven 3+1D hydrodynamical simulations with viscous effects followed by the statistical particle emission to match multiplicities observed experimentally [37]. The modified Pearson correlation coefficient is then extracted from the final-state particles. The predictions for all harmonics are consistent with the data within the large model uncertainties except for the most central collisions where the predictions underestimate the measured \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \) and for the semi-peripheral collisions, for \(N_{\mathrm{part}} \sim 130\), where the predictions overestimate the \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \) and underestimate \( \rho (v_{4}\{2\}^2, [p_{\mathrm{T}}]) \).

Fig. 7
figure 7

The PCC \(\rho ( v_{n}\{2\}^2, [p_{\mathrm{T}}] )\) for a \(n = 2\), b \(n = 3\), and c \(n = 4\) in Pb+Pb collisions as a function of \(N_{\mathrm{part}} \) for three \(p_{\mathrm{T}}\) ranges. The statistical and systematic uncertainties are shown as vertical error bars and boxes, respectively. A comparison with model predictions [37] is also shown with a line added to guide the eye

7 Summary

The first measurement of the modified PCC \( \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]) \), which quantifies the correlation between the flow harmonics and the mean transverse momentum, is performed by ATLAS experiment at the LHC. The measurement uses \(22~\mu \text {b}^{-1} \) of Pb+Pb data and \(28~\text{ nb }^{-1}\) of \(p\)+Pb data at the same centre-of-mass energy per nucleon pair of 5.02 \(\text {Te} \text {V}\).

The correlation coefficient for several charged-particle \(p_{\mathrm{T}}\) ranges is measured as a function of the number of charged particles \(N_{\mathrm{ch}}\) and, in Pb+Pb collisions, the average number of nucleons participating in the collision, \(N_{\mathrm{part}}\). For the \(2^{\mathrm{nd}}\)-, \(3^{\mathrm{rd}}\)-, and \(4^{\mathrm{th}}\)-order harmonics, the measured quantities exhibit a dependence on the choice of charged-particle \(p_{\mathrm{T}}\) range. Measurements with an upper limit of 5 \(\text {Ge} \text {V}\) on \(p_{\mathrm{T}}\) indicate a stronger correlation than those with an upper limit of 2 \(\text {Ge} \text {V}\). For mid-central and central collisions, when varying the lower \(p_{\mathrm{T}}\) limit, consistent values of \( \rho (v_{3}\{2\}^2, [p_{\mathrm{T}}]) \) and \( \rho (v_{4}\{2\}^2, [p_{\mathrm{T}}]) \) coefficients are obtained, whereas for the \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \) coefficient a difference of 10–20% is seen. As a function of event activity, for Pb+Pb collisions, a strong positive correlation \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \) is observed in mid-central and central collisions while negative values are measured for peripheral events. The correlation \( \rho (v_{3}\{2\}^2, [p_{\mathrm{T}}]) \) is found to be weaker, yet non-zero. The values of \( \rho (v_{4}\{2\}^2, [p_{\mathrm{T}}]) \) are also positive in the studied centrality range. Non-monotonic behaviour is observed in central Pb+Pb collisions. That trend observed for \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \) in Pb+Pb collisions is in line with expectations drawn from the ALICE results [20]. In \(p\)+Pb collisions, the value of \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \) is negative and approximately independent of \(N_{\mathrm{ch}}\).

The modified PCC is a valuable tool for studying the dynamics of heavy-ion collisions. It provides a reliable estimate of the magnitude of correlations calculated using finite multiplicities. In comparison with existing results, it allows quantitative comparisons between the experimental data and theoretical models. The precise measurements of this observable, presented in this paper, provide useful insights into the interplay of the azimuthal anisotropies (azimuthal flow) and the mean event \(p_{\mathrm{T}}\) (radial flow), providing input for a better understanding of QGP dynamics and for constraining the theoretical models. The obtained \( \rho (v_{n}\{2\}^2, [p_{\mathrm{T}}]) \) coefficients for \(0.5< p_{\mathrm{T}} < 2~\text {Ge} \text {V}\) were compared with a theoretical prediction based on the pressure-driven 3+1D hydrodynamical simulations with viscous effects. The predictions for all harmonics are consistent with the data within the large model uncertainties. The only exception are the most central collisions, where the predictions underestimate the measured \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \) and the semi-peripheral collisions, where the predictions overestimate the \( \rho (v_{2}\{2\}^2, [p_{\mathrm{T}}]) \) and underestimate \( \rho (v_{4}\{2\}^2, [p_{\mathrm{T}}]) \). Sizeable positive correlations observed for non-peripheral Pb+Pb collisions support a qualitatively expected scenario in which the azimuthal flow originates from the pressure gradients.

In small system collisions the magnitude of the transverse flow is expected to be very sensitive to the size of the initial source in the hydrodynamic model. In particular, in the compact source scenario in \(p\)+Pb collisions, the smaller source sizes are expected to yield larger transverse flow and smaller initial eccentricities. The negative sign of the modified PCC measured in \(p\)+Pb collisions seems to support the compact source scenario, and indicates the role of the initial conditions in these systems.