Correlated long-range mixed-harmonic fluctuations measured in pp, p plus Pb and low-multiplicity Pb plus Pb collisions with the ATLAS detector

Physics Letters B 789 (2019) 444–471

Contents lists available atScienceDirect

Physics

Letters

B

www.elsevier.com/locate/physletb

Correlated

long-range

mixed-harmonic

fluctuations

measured

in

pp,

p+Pb

and

low-multiplicity

Pb+Pb

collisions

with

the

ATLAS

detector

.

The

ATLAS

Collaboration



a

r

t

i

c

l

e

i

n

f

o

a

b

s

t

r

a

c

t

Articlehistory:

Received6July2018

Receivedinrevisedform7October2018 Accepted13November2018

Availableonline2January2019 Editor: D.F.Geesaman

Correlations oftwoflowharmonics vnand vm viathree- andfour-particlecumulantsaremeasuredin

13TeVpp,5.02TeVp+Pb,and2.76TeVperipheralPb+PbcollisionswiththeATLASdetectorattheLHC.

The goalistounderstandthe multi-particlenature ofthelong-rangecollectivephenomenoninthese

collisionsystems.Thelargenon-flowbackgroundfromdijetproductionpresentinthestandardcumulant

methodissuppressedusingamethodofsubeventcumulantsinvolvingtwo,threeand foursubevents

separatedinpseudorapidity.Theresultsshowanegativecorrelationbetweenv2and v3 andapositive

correlationbetweenv2andv4forallcollisionsystemsandoverthefullmultiplicityrange.However,the

magnitudes ofthecorrelationsare foundtodependontheeventmultiplicity,thechoiceoftransverse

momentumrangeand collisionsystem. Therelativecorrelation strength,obtainedbynormalisationof

the cumulants withthe v2

nfrom atwo-particle correlationanalysis, issimilar inthe threecollision

systemsanddependsweaklyontheevent multiplicityand transversemomentum.Theseresultsbased

onthesubevent methodsprovidestrongevidence ofasimilar long-rangemulti-particle collectivityin

pp,p+PbandperipheralPb+Pbcollisions.

©2019TheAuthor.PublishedbyElsevierB.V.ThisisanopenaccessarticleundertheCCBYlicense

(http://creativecommons.org/licenses/by/4.0/).FundedbySCOAP3_.

1. Introduction

One of the goals in the studies of azimuthal correlations in high-energy nuclear collisions at the Relativistic Heavy Ion Col-lider(RHIC)andtheLargeHadronCollider(LHC)istounderstand the multi-parton dynamics of QCD in the strongly coupled non-perturbative regime [1]. Measurements of azimuthal correlations in small collision systems, such as pp, p+A or d+A collisions, have revealed the ridge phenomenon [2–6]: enhanced produc-tion of particle pairs at small azimuthal angle separation,

φ

, extended over a wide range of pseudorapidity separation,



η

. The azimuthal structure has been related to harmonic modula-tion of particle densities, characterised by a Fourier expansion, dN

/

d

φ

∝

1

+

2



∞_n=1vncos n

(φ

− n

)

, where vn and



n repre-sent the magnitude and the event-plane angle of the nth-order flowharmonic.Theyarealsoconvenientlyrepresentedbytheflow vector: Vn

=

vneinn.The vn are known to depend onthe colli-sionsystem,buthaveweakdependenceoncollisionenergies [6,7]. Theridgereflectsmulti-partondynamicsearlyinthecollisionand hasgeneratedsignificantinterestinthehigh-energyphysics com-munity. A key question is whetherthe long-range multi-particle collectivityreflects initial momentumcorrelation fromgluon

sat- _{E-mailaddress:atlas.publications@cern.ch.}

uration effects [8], ora final-statehydrodynamic response to the initialtransversecollisiongeometry [9].

Further insight into the ridge phenomenon is obtained via a multi-particle correlation technique,known ascumulants, involv-ing three or more particles [10–12]. The multi-particle cumu-lants probe the event-by-event fluctuation of a single flow har-monic vn, as well as the correlated fluctuations between two flowharmonics, vn andvm.Theseevent-by-eventfluctuationsare often represented by probability density distributions p

(

vn

)

and p

(

vn

,

vm

)

, respectively. For instance, the four-particle cumulants cn{4

}

=



v4_n



−

2



v_n2



2 constrain the width of p

(

vn

)

[10], while thefour-particlesymmetriccumulantsscn,m{4

}

=



v2_nv2_m



−



v2_n

 

v2_m



quantifythelowest-ordercorrelationbetweenvn andvm[12].The three-particle asymmetric cumulants such asacn{3

}

=



V2_nV∗_2n



=



v2

nv2ncos 2n

(

n

− 

2n

)



[5,13] aresensitivetocorrelations involv-ingboththeflowmagnitudevn andflowphase



n.

Oneofthechallengesinthestudyofazimuthalcorrelationsin smallcollisionsystemsishowto distinguishthelong-rangeridge from“non-flow”correlationsinvolvingonlyafewparticles,suchas resonancedecays,jets, ordijets.Fortwo-particlecorrelations, the non-flowcontributioniscommonlysuppressedbyrequiringalarge



η

gapbetweenthe two particles ineach pair anda peripheral subtractionprocedure [3–5,7,14,15]. Formulti-particlecumulants, the non-flowcontributions canbe suppressedbyrequiring corre-lation betweenparticles fromdifferentsubevents separatedin

η

,

https://doi.org/10.1016/j.physletb.2018.11.065

0370-2693/©2019TheAuthor.PublishedbyElsevierB.V.ThisisanopenaccessarticleundertheCCBYlicense(http://creativecommons.org/licenses/by/4.0/).Fundedby SCOAP3_.

5

.

02 TeV andlow-multiplicityPb+Pbcollisionsat

√

sNN

=

2

.

76 TeV.

They are obtained using two-, three- and four-subevent cumu-lant methods and are compared with results from the standard cumulantmethod.The cumulantsare normalisedby the



v2n



ob-tainedfromatwo-particlecorrelationanalysis [7] toquantifytheir relative correlation strength. The measurements suggest that the resultsobtainedwiththe standard methodare strongly contami-natedbycorrelationsfromnon-flowsources. Theresultsobtained withthethree-subeventmethodorthefour-subeventmethod pro-videnewevidenceoflong-rangethree- orfour-particleazimuthal correlations.

TheLetterisorganisedasfollows.DetailsoftheATLASdetector, thetriggersystem, datasets,aswellaseventandtrackselections areprovided inSections 2to 4.Section 5describesthe standard andsubeventcumulantmethodsusedinthisanalysis.Theanalysis procedureandsystematicuncertaintiesaredescribedinSections6

and7,respectively.Themeasuredcumulantsarepresentedin Sec-tion8.AsummaryisgiveninSection9.

2. Detectorandtrigger

TheATLASdetector [20] providesnearlyfull solid-angle cover-agearound the collision point withtracking detectors, calorime-ters,andmuonchambers, andiswell suited formeasurementof multi-particlecorrelationsoveralarge pseudorapidityrange.1 The measurements were performed using primarily the inner detec-tor(ID), minimum-biastrigger scintillators(MBTS) and the zero-degreecalorimeters(ZDC).TheIDdetectschargedparticleswithin

|

η

|

<

2

.

5 using a combination of a silicon pixel detector, a sili-conmicrostripdetector (SCT),andastraw-tubetransitionradiation tracker, all immersed in a 2 T axial magnetic field [21]. An ad-ditionalpixellayer, the “insertable B-layer”(IBL) [22] is installed betweenthe Run-1 (2010–2013)andRun-2(2015–2018) periods. The MBTS detects charged particles within 2

.

1

 |

η

|



3

.

9 using two hodoscopes ofcounters positioned at z

= ±

3

.

6 m. The ZDC, usedonlyinp+PbandPb+Pbcollisions,arepositionedat

±

140 m from the collision point, and detect neutral particles, primarily neutronsandphotons,with

|

η

|

>

8

.

3.

The ATLAS trigger system [23,24] consistsof a first-level (L1) trigger implemented using a combination of dedicated electron-icsandprogrammablelogic,andahigh-level trigger(HLT) imple-mentedinprocessors.TheHLTreconstructscharged-particletracks

1 _{ATLAStypicallyusesaright-handedcoordinatesystemwithitsoriginat the}

nominalinteractionpoint(IP)inthe centreofthe detectorandthe z-axisalong thebeampipe.Thex-axispointsfromtheIPtothecentreoftheLHCring,andthe

y-axispointsupward.Cylindricalcoordinates(r,φ)areusedinthetransverseplane,

φbeingtheazimuthalanglearoundthebeampipe.Bydefault,thepseudorapidityis definedintermsofthepolarangleθasη= −ln tan(θ/2).However,forasymmetric

p+PborPb+p collisions,the−z directionisalwaysdefinedasthedirectionofthe Pbbeam.

andHMTtriggers.Theminimum-biastriggerrequiredeitherahit inatleastoneMBTScounter,orahitinatleastoneMBTScounter oneachside,oratleastonereconstructedtrackattheHLTseeded by a random trigger at L1. More detailed information about the triggersusedforthepp and p+Pbdataandtheirperformancecan befoundinRefs. [7,25] andRefs. [5,26],respectively.

3. DatasetsandMonteCarlosimulations

This analysisis basedon ATLAS datasets corresponding to in-tegrated luminosities of 0.9 pb−1 of pp data recorded at

√

s

=

13 TeV, 28 nb−1 of p+Pbdatarecordedat

√

sNN

=

5

.

02 TeV, and

7 μb−1 of Pb+Pb data at

√

sNN

=

2

.

76 TeV. The 2

.

76 TeV Pb+Pb

datawerecollectedin2010.The p+Pbdataweremainlycollected in2013,butalsoinclude0.3 nb−1 ofdatacollectedin2016,which increase thenumber ofeventsatmoderatemultiplicity (see Sec-tion4).Duringboth p+Pbruns,theLHCwasconfiguredtoprovide a4 TeVprotonbeamanda1.57 TeVper-nucleonPbbeam, which produced collisions at

√

sNN

=

5

.

02 TeV, with a rapidity shift of

0

.

465 of the nucleon–nucleon centre-of-mass frame towards the protonbeamdirectionrelativetotheATLASrestframe.The direc-tionofthePbbeamisalwaysdefinedtohavenegative pseudora-pidity. The 13 TeV pp data were collected during several special runsoftheLHCwithlowpile-upin2015and2016.Asummaryof thedatasetsusedinthisanalysisisshowninTable1.

Thetrackreconstructionefficiencywasdeterminedusing simu-latedMonteCarlo(MC)eventsamples (Section4).The pp events were simulated withthe Pythia8 MC event generator [27] using theA2set oftunedparameters withMSTW2008LO parton distri-bution functions [28]. The HIJING eventgenerator [29] was used to producePb+Pb and p+Pbcollisions withthesame energyand the sameboost ofthe centre-of-mass systemasinthe data.The detector response was simulated using Geant_{4 [30,31] with} de-tectorconditionsmatchingthoseduring thedata-taking.The sim-ulated events and data events are reconstructed with the same algorithms. The MC sample for Pb+Pb events in the multiplicity region of interest is very small, and so the track reconstruction efficiencyforPb+Pbwastakenfromthelargerp+Pbsample recon-structedwiththesamereconstructionalgorithm.Theefficiencyin p+Pbeventswas found to be consistentwiththe efficiencyfrom thePb+PbMCsimulation [17].

4. Eventandtrackselection

The offline eventselection forthe pp and p+Pbdata requires atleastonereconstructedvertexwithitslongitudinalposition sat-isfying

|

zvtx

|

<

100 mmrelative tothe nominalinteractionpoint.

Thevertexisrequiredtohaveatleasttwo associatedtrackswith pT

>

0

.

4 GeV. The mean number of collisions per bunch

446 The ATLAS Collaboration / Physics Letters B 789 (2019) 444–471

p+Pbdata, and0.001–0.006 forthe 2016 p+Pbdata. In orderto suppressadditional interactions in the same bunch crossing (re-ferredto aspile-up)in pp collisions,eventscontaining additional verticeswithatleast fourassociatedtracks arerejected. In p+Pb collisions,eventswithmorethanone goodvertex, definedasany vertexforwhichthescalarsumofthepToftheassociatedtracks

isgreater than5 GeV, arerejected. The remainingpile-up events are further suppressedby using the signal inthe ZDCin the di-rection of the Pb beam. This signal is calibrated to the number ofdetectedneutrons, Nn,byusingthelocationofthepeak corre-spondingtoasingleneutron.ThedistributionofNn ineventswith pile-upisbroaderthanthatfortheeventswithoutpile-up.Hence asimplerequirementontheZDCsignaldistributionisusedto fur-thersuppresseventswithpile-up,whileretainingmorethan98% of eventswithout pile-up. The impact ofresidual pile-up, atthe levelof



10−3,isstudiedbycomparingtheresultsobtainedfrom datawithdifferent

μ

values.

TheofflineeventselectionforthePb+Pbdatarequires

|

zvtx

|

<

100 mm.Theselection alsorequiresatimedifference

|

t

|

<

3 ns betweensignalsintheMBTStriggercountersoneithersideofthe interactionpoint tosuppressnon-collisionbackgrounds. A coinci-dence between the ZDC signals at forward and backward pseu-dorapidityisrequiredtoreject avarietyofbackgroundprocesses, whilemaintaininghighefficiencyforinelasticprocesses.The frac-tionofeventswithmorethanoneinteractionafterapplyingthese selectioncriteriaislessthan10−4_.

Charged-particletracksandcollisionverticesarereconstructed using algorithms optimised forimproved performance forRun-2. In order to comparedirectly with the pp and p+Pbsystems us-ing event selections based on the multiplicity of the collisions, a subset of datafrom low-multiplicityPb+Pb collisions, collected during the 2010 LHC heavy-ion run with a minimum-bias trig-ger,wasanalysedusingthesametrackreconstructionalgorithmas thatusedforp+Pbcollisions.ForthePb+Pband2013p+Pb analy-ses,tracksarerequiredtohavea pT-dependentminimumnumber

ofhitsin theSCT. The transverse (d0) andlongitudinal(z0 sin

θ

)

impactparametersofthetrackrelativeto thevertexare required tobelessthan1.5 mm.Additionalrequirements

|

d0

|/

σ

d0

<

3 and

|

z0sin

θ

|/

σ

z0

<

3 are imposed, where

σ

d0 and

σ

z0 are the

un-certainties of the transverse and longitudinal impact parameter values,respectively. A more detaileddescription of the track se-lectionforthe2010Pb+Pbdataand2013p+Pbdatacanbefound inRefs. [5,17].

Forall thedatatakensincethestartofRun-2,thetrack selec-tioncriteriamakeuseoftheIBL,asdescribedinRefs. [14,25].For the pp and2016 p+Pbanalyses,thetracks arerequiredtosatisfy

|

dBL₀

|

<

1

.

5 mm and

|

z0sin

θ

|

<

1

.

5 mm, where dBL₀ is the

trans-verseimpactparameterofthetrackrelativetothebeamline (BL). Thecumulantsarecalculatedusingtrackspassingtheabove se-lectionrequirements,andhaving

|

η

|

<

2

.

5 and 0

.

3

<

pT

<

3 GeV

or 0

.

5

<

pT

<

5 GeV. These two pT ranges are chosen because

they were often used in the previous ridge measurements at the LHC [6,7,14,15,17]. However, to count the number of recon-structed charged particles for event-class definition (denoted by Nrec

ch), tracks with pT

>

0

.

4 GeV and

|

η

|

<

2

.

5 are used for

com-patibility with the requirements in the HLT selections described above. Due to different trigger requirements, most of the p+Pb events with Nrec_ch

>

150 are provided by the 2013 dataset,while the2016datasetprovidesmostoftheeventsatlowerNrec

ch.

The efficiencyof thecombined trackreconstruction and track selection requirements is estimated using MC samples recon-structed withthesame algorithmsandselection requirements as indata.Efficiencies,

(

η

,

pT

)

, areevaluated asafunction oftrack

η

,pTandthenumberofreconstructedcharged-particletracks,but

averagedoverthefullrangeinazimuth.Theefficiencies are

simi-larforeventswiththesamemultiplicity.Forallcollisionsystems, the efficiency increases by about 4% as track pT increases from

0.3 GeVto 0.6GeV.Above 0.6GeV,the efficiencyisindependent of pT and reaches 86% (72%) for Run-1 pp and p+Pb, and 83%

(70%)forPb+PbandRun-2 p+Pbcollisions,at

η

≈

0 (

|

η

|

>

2).The efficiency is independent ofthe event multiplicity for Nrec_ch

>

40. Forlower-multiplicityeventstheefficiencyissmallerbyupto3% duetobroaderd0 andz0sin

θ

distributions [17].

The fraction of falsely reconstructed charged-particle tracks is alsoestimatedandfoundtobenegligiblysmallinalldatasets.This fractiondecreaseswithincreasingtrackpT,andevenatthelowest

transversemomentaof0.3 GeVitisbelow1%ofthetotalnumber oftracks.Therefore,thereisnocorrectionforthepresenceofsuch tracksintheanalysis.

In the simulated events, the reconstruction efficiency reduces the measured charged-particle multiplicity relative to the gener-atedmultiplicityforprimarychargedparticles.Acorrectionfactor b is used to correct Nrec_ch to obtain the efficiency-corrected aver-age number of charged particles per event,



Nch

 =

b



Nrec

ch



. The valueofthecorrectionfactorisobtainedfromtheMCsamples de-scribed above, and is found to be nearly independent of Nrec_ch in therangeusedinthisanalysis, Nrec

ch

<

400.Itsvalueandthe

asso-ciated uncertainties are b

=

1

.

29

±

0.05 forthe Pb+Pband2013 p+Pbcollisions andb

=

1

.

18

±

0.05 forRun-2 p+Pb and pp colli-sions [32].Bothscn,m{4

}

andac2

{

3

}

arethenstudiedasafunction

of



Nch



.

5. Cumulantmethod

The multi-particlecumulantmethod [10] hastheadvantageof directly reducing non-flow correlations from jets and dijets. The mathematical framework for the standard cumulant is based on theQ-cumulantsdiscussedinRefs. [11,12,33].Itwasextended re-cently to the caseof subevent cumulants in Refs. [13,16]. These methodsarebrieflysummarisedbelow.

5.1. Cumulantsinthestandardmethod

The standardcumulantmethodcalculates k-particleazimuthal correlations,

{

k

}

, in one event using a complex number nota-tion [11,12]:

{

2

}

n

 =



ein(φ1−φ2)



_,

{

₃

}

n

 =



ein(φ1+φ2−2φ3)



_,



{

4

}

n,m



=



ein(φ1−φ2)+im(φ3−φ4)



,

(1)

where “



” denotes a single-event average over all pairs, triplets orquadruplets,respectively.The averagesfromEq. (1) canbe ex-pressed intermsofper-particle normalisedflowvectors

q

n;l with l

=

1

,

2

...

ineachevent [11]: q_n;l

≡



j wl_jeinφj





j wl_j

,

(2)

wherethesumrunsoveralltracksintheeventandwjisaweight assigned to the jth track. This weight is constructed to correct for both detectornon-uniformity andtracking inefficiency as ex-plainedinSection6.

The multi-particle asymmetric and symmetric cumulants are obtainedfrom

{

k

}

as:

acn

{

3

} = {

3

}

n

 ,

scn,m

{

4

} =



{

4

}

n,m



 − {

2

}

n

 {

2

}

m

 ,

(3) where“



”representsa weighted averageof

{

k

}

over anevent ensemble with similar Nrec

Inthestandardcumulantmethoddescribedabove,allk-particle multipletsinvolvedin

{

k

}n



and



{

k

}n

,m



areselectedusingtracks inthe entire ID acceptanceof

|

η

|

<

η

max

=

2

.

5. To suppress

fur-therthenon-flowcorrelationsthattypicallyinvolveafewparticles withina localisedregion in

η

,the tracksare dividedintoseveral subevents, each covering a unique

η

interval. The multi-particle correlations are then constructed by only correlating tracks be-tweendifferentsubevents.

Inthe two-subevent cumulantmethod,the tracks are divided into two subevents, labelled by a and b, according to

−

η

max

<

η

a

<

0 and0

≤

η

b

<

η

max.Theper-eventk-particleazimuthal

cor-relationsareevaluatedas:

{

2

}

n



a|b

=



ein(φa1−φb2)



,

{

3

}

n



2a|b

=



ein(φa1+φa2−2φ3b)



,



{

4

}

n,m



2a|2b

=



ein(φa1−φ2b)+im(φa3−φ4b)



,

(6)

where the superscript or subscript a (b) indicates tracks chosen fromthe subeventa (b).Herethe three- andfour-particle cumu-lantsaredefinedas:

ac2an|b

{

3

} = {

3

}

n



2a|b

,

sc2an,m|2b

{

4

} =



{

4

}

n,m





2a|2b

− {

2

}

n



a|b

{

2

}

m



a|b

.

Thetwo-subeventmethodsuppresses correlationswithina single jet(intra-jetcorrelations),sinceparticles fromone jetusually fall inonesubevent.

Inthe three-subevent cumulant method,tracks in each event are divided into three subevents a, b and c, each covering one third of the

η

range,

−

η

max

<

η

a

<

−

η

max

/

3,

|

η

b|

≤

η

max

/

3 and

η

max

/

3

<

η

c

<

η

max.Themulti-particleazimuthalcorrelationsand

cumulantsarethenevaluatedas:

{

3

}

n



a,b|c

=



ein(φa1+φb2−2φ3c)



,



{

4

}

n,m



a,b|2c

=



ein(φ1a−φ2c)+im(φb3−φc4)



,

(7) and aca_n,b|c

{

3

} = {

3

}

n



a,b|c

,

scan,,bm|2c

{

4

} = {

4

}

n,m



a,b|2c

− {

2

}

n



a|c

{

2

}

m



b|c

.

(8) Since a dijet event usually produces particles in at most two subevents, the three-subevent method efficiently suppresses the non-flow contribution from inter-jet correlations associated with dijets. To maximise the statistical precision, the

η

range for subevent a is swapped with that for subevent b or c, and the resultsareaveragedtoobtainthefinalvalues.

precision, the

η

rangesfor thefour subevents are swapped with eachother,andtheresultsareaveragedtoobtainthefinalvalues. 5.3. Normalisedcumulants

Althoughthecumulantsreflectthenatureofthecorrelation be-tween vn and vm,their magnitudesalsodependonthesquare of singleflow harmonicsv2

n andv2m,seeEq. (4).Thedependenceon the single flow harmonics can be scaled out via the normalised cumulants [34,35]: nsc2,3

{

4

} =

sc2,3

{

4

}

v2

{

2

}

2v3

{

2

}

2

=



v2₂v2₃





v2₂

 

v2₃

 −

1

,

(11) nsc2,4

{

4

} =

sc2,4

{

4

}

v2

{

2

}

2v4

{

2

}

2

=



v2 2v24





v2 2

 

v2 4

 −

1

,

(12) nac2

{

3

} =

ac2

{

3

}

2v2

{

2

}

4

+

c2

{

4

}



c4

{

2

}

=



v2₂v4cos 4

(

2

− 

4

)





v4₂

 

v2₄



,

(13) wherethe vn

{

2

}

2

=



v2_n



areflowharmonicsobtainedusinga two-particlecorrelationmethodbasedonaperipheralsubtraction tech-nique [7,14], andc2

{

4

}

=



v4 2



−

2



v2 2



2

are four-particlecumulant resultsfromRefs. [17,18]. Thisdefinitionfornac2

{

3

}

ismotivated

byRef. [36].

6. Analysisprocedure

Themeasurement ofthescn,m{4

}

andac2

{

3

}

followsthesame

analysis procedure as for the four-particle cumulants cn{4

}

in Ref. [18].Themulti-particlecumulantsarecalculatedinthreesteps using charged particles with

|

η

|

<

2

.

5. In the first step,

{

2

}n

,

{

3

}n



and



{

4

}n

,m



from Eqs. (1), (6), (7) and (9) are calculated for each event from particles in one of two different pT ranges,

0

.

3

<

pT

<

3 GeVand 0

.

5

<

pT

<

5 GeV. The numbersof

recon-structedchargedparticlesinthesepT rangesaredenotedby Nsel1_ch

andNsel2

ch ,respectively.

In thesecond step, thecorrelators

{

k

}

for 0

.

3

<

pT

<

3 GeV

(0

.

5

<

pT

<

5 GeV) areaveraged overeventswiththe same Nsel1_ch

(N_chsel2)toobtain

{

k

}

,andthensc2,3

{

4

}

,sc2,4

{

4

}

andac2

{

3

}

.The

sc2,3

{

4

}

,sc2,4

{

4

}

andac2

{

3

}

values are then averaged in broader

multiplicityrangesoftheeventensemble,weightedbynumberof events,toobtainstatisticallysignificantresults.

In the third step, the sc2,3

{

4

}

, sc2,4

{

4

}

and ac2

{

3

}

values

ob-tained fora given Nsel1_ch or Nsel2_ch are mapped to



Nrec_ch



,the aver-agenumberofreconstructedchargedparticleswith pT

>

0

.

4 GeV.

448 The ATLAS Collaboration / Physics Letters B 789 (2019) 444–471

Themappingprocedureisnecessarysothat sc2,3

{

4

}

,sc2,4

{

4

}

and

ac2

{

3

}

obtainedforthetwodifferent pT rangescan becompared

usingacommonx-axisdefinedby



Nrec_ch



.The



Nrec_ch



valueisthen converted to



Nch



, the efficiency-corrected average number of

chargedparticleswithpT

>

0

.

4 GeV,asdiscussedinSection4.

In order to account for detector inefficiencies and non-uni-formity,particleweightsusedinEq. (2) aredefinedas:

w

(φ,

η

,

pT

)

=

d

(φ,

η

)/

(

η

,

pT

) .

Theadditionalweightfactord

(φ,

η

)

accountsfornon-uniformities in the azimuthal acceptance of the detector as a function of

η

. Allreconstructedchargedparticleswith pT

>

0

.

3 GeVareentered

intoa two-dimensional histogram N

(φ,

η

)

,andthe weightfactor isthenobtainedasd

(φ,

η

)

≡ 

N

(

η

)

 /

N

(φ,

η

)

,where



N

(

η

)



isthe trackdensityaveraged over

φ

in thegiven

η

bin.Thisprocedure removesmostofthe

φ

-dependentnon-uniformity inthedetector acceptance [17].

In order to calculate the normalised cumulants from Eqs. (11)–(13),theflowharmonicsvn{2

}

areobtainedfroma“template fit”oftwo-particle

φ

correlationasdescribedinRefs. [7,14].The vn{2

}

values are calculated identically to the procedure used in the previous ATLAS publications [7,14], butare furthercorrected fora bias,which exists onlyif vn{2

}

changes with Nrec_ch.The de-tailsofthecorrectionprocedurearegivenintheAppendixAand arediscussedbrieflybelow.

ThestandardprocedureofRefs. [7,14] firstconstructsa

φ

dis-tributionforpairsoftrackswith

|

η

|

>

2:theper-trigger-particle yield Y

(φ)

foragiven Nrec

ch range.Thedominatingnon-flow jet

peak at

φ

∼

π

is estimatedusing low-multiplicity events with Nrec_ch

<

20 andseparatedviaatemplatefitprocedure,andthe har-monic modulation of the remaining component is taken as the vn{2

}

2 _[7]: Y

(φ)

=

F Y

(φ)

peri

+

Gtmp



1

+

2 ∞



n=2 vn

{

2

,

tmp

}

2cos n

φ

,

where superscripts “peri” and “tmp” indicate quantities for the Nrec_ch

<

20 eventclassandquantitiesafterthetemplatefitforthe eventclassofinterest,respectively.ThescalefactorF andpedestal Gtmp are fixed by the fit, and vn{2

,

tmp

}

are calculated from a Fouriertransform. Thisprocedureimplicitlyassumesthat vn{2

}

is independentofNrec

ch,andrequiresasmallcorrectionifvn{2

}

does

changewithN_chrec(AppendixA).In p+PbandPb+Pbcollisions,this correctionin the Nrec_ch

>

100 regionamounts toa 2–6%reduction forv2

{

2

,

tmp

}

anda4–9%reductionforv3

{

2

,

tmp

}

andv4

{

2

,

tmp

}

.

The correction is smaller for v2

{

2

,

tmp

}

in pp collisions as it is

nearlyindependentofN_chrec[7].

7. Systematicuncertainties

The evaluation of the systematic uncertainties follows closely theprocedureestablishedforthefour-particlecumulantscn{4

}

and describedinRef. [18].Themainsourcesofsystematicuncertainties arerelatedto thedetectorazimuthal non-uniformity,track selec-tion,trackreconstructionefficiency,triggerefficiencyandpile-up. Duetotherelativelypoorstatisticsandlargernon-floweffects,the systematicuncertainties are typically larger in pp collisions. The systematic uncertainties are also generally larger, in percentage, forfour-particlecumulantsscn,m{4

}

thanforthethree-particle

cu-mulantsac2

{

3

}

,sincethe

|

scn,m{4

}|

valuesaremuchsmallerthan

thoseforac2

{

3

}

.Thesystematicuncertaintiesaregenerallysimilar

among the two- and three- and four-subevent methods, but are different from those for the standard method, which is strongly

influenced by non-flow correlations. The following discussion fo-cusesonthethree-subeventmethod,whichisthedefaultmethod usedtopresentthefinalresults.

The effect of detector azimuthal non-uniformity is accounted for usingthe weight factor d

(φ,

η

)

. The impact of theweighting procedure is studiedby fixing the weight to unityand repeating the analysis. The results are mostly consistent with the nominal results. The corresponding uncertainties for scn,m{4

}

vary in the

range of0–4%,0–2% and1–2% in pp, p+PbandPb+Pbcollisions, respectively.Theuncertaintiesforac2

{

3

}

varyintherangeof0–2%

in pp collisions, and 0–1% in p+PbandPb+Pb collisions, respec-tively.

The systematicuncertaintyassociated withthetrack selection is estimated by tightening the

|

d0

|

and

|

z0sin

θ

|

requirements.

They are each varied from the default requirement of less than 1.5 mm tolessthan1 mm.In p+PbandPb+Pbcollisions, the re-quirementon thesignificanceofimpactparameters,

|

d0

|/

σ

d0 and

|

z0sin

θ

|/

σ

z0 are also varied fromlessthan 3 to lessthan 2.For

eachvariation,thetrackingefficiencyisre-evaluatedandthe anal-ysisisrepeated. Forac2

{

3

}

,whichhasalargeflowsignal,the

dif-ferencesfromthenominalresultsareobservedtobelessthan2% forallcollisionsystems.Forscn,m{4

}

,forwhichthesignalissmall,

the differencesfrom the nominal results are found to be in the rangeof2–10%in pp collisions,2–7%in p+Pbcollisionsand2–4% inPb+Pbcollisions.Thedifferencesaresmallerforresultsobtained for0

.

5

<

pT

<

5 GeVthanthoseobtainedfor0

.

3

<

pT

<

3 GeV.

Previousmeasurementsindicatethattheazimuthalcorrelations (both the flow andnon-flow components) have a strong depen-denceon pT,buta relativelyweakdependenceon

η

[5,7].

There-fore, pT-dependent systematic effects in the trackreconstruction

efficiencycouldaffectthecumulantvalues.Theuncertaintyinthe trackreconstruction efficiencyismainlydueto differencesinthe detectorconditionsandmaterial descriptionbetweenthe simula-tion and the data. The efficiency uncertainty varies between 1% and4%,depending ontrack

η

andpT [7,17].Its impacton

multi-particlecumulantsisevaluatedbyrepeatingtheanalysiswiththe trackingefficiencyvariedupanddownbyitscorresponding uncer-taintyasafunctionoftrackpT.Forthestandardcumulantmethod,

which is more sensitive to jets and dijets, the evaluated uncer-tainty amountsto2–6%in pp collisionsandlessthan2% inp+Pb collisionsfor



Nch

 >

100.Forthesubeventmethods,theevaluated

uncertaintyistypicallylessthan3%formostofthe



Nch



ranges.

Most eventsin pp and p+Pb collisions are collected with the HMT triggerswith severalonline Nrec_ch thresholds. Inorder to es-timate the possible bias due to trigger inefficiency asa function of



Nch



,theoffline Nrecch requirementsarechanged suchthat the

HMT triggerefficiencyisatleast 50%or80%. Theresults are ob-tained independently for each variation. These results are found to be consistent with each other forthe subevent methods, and show some differencesforthe standard cumulantmethodin the low



Nch



region.Thenominalanalysisisperformedusingthe50%

efficiency selection and the differences betweenthe nominal re-sults andthose fromthe80% efficiency selection areincluded in thesystematicuncertainty.Thechangesforpp collisionsareinthe rangeof5–15%forsc2,3

{

4

}

,2–8%forsc2,4

{

4

}

and1–5%forac2

{

3

}

.

Theranges forp+Pbcollisionsaremuchsmallerduetothemuch sharperturn-onofthetriggerefficiencyandlargersignal:theyare estimatedtobe 1–3%forsc2,3

{

4

}

,2–4%forsc2,4

{

4

}

and1–2%for

ac2

{

3

}

.

Inthisanalysis,apile-uprejectioncriterionisappliedtoreject eventscontainingadditionalverticesin pp andp+Pbcollisions.In order to check the impact of residualpile-up, the analysisis re-peated without the pile-uprejection criterion.No differences are observed in p+Pbcollisions, asisexpected sincethe

μ

valuesin p+Pbaremodest.Forthe13 TeV pp dataset,thedifferenceswith

The vn{2

}

values used to obtain normalised cumulants from Eqs. (11)–(13) aremeasuredfollowingtheprescriptionofthe pre-viousATLAS publications [7,14],resulting invery similar system-aticuncertainties. The correction for the biasof the template fit procedure,asdescribedinSection6,reducesthesensitivitytothe choiceoftheperipheral Nrec

ch bin.The uncertaintiesofnormalised

cumulantsare obtainedby propagationof theuncertainties from theoriginalcumulantsandvn{2

}

,takingintoaccountthatthe cor-relatedsystematicuncertaintiespartiallycancelout.

8. Results

Theresults are presented intwo parts. Section 8.1presents a detailedcomparisonbetweenthe standard method andsubevent methods to demonstrate the ability of the subevent methods to suppress non-flow correlations. Section 8.2 compares the cumu-lantsamongpp, p+PbandPb+Pbcollisionstoprovideinsightinto thecommonnatureofcollectivityinthesesystems.

8.1.Comparisonbetweenstandardandsubeventmethods

The top row of Fig. 1 compares the sc2,3

{

4

}

values obtained

fromthe standard,two-, three- and four-subevent methodsfrom pp collisionsin0

.

3

<

pT

<

3 GeV(leftpanel)and0

.

5

<

pT

<

5 GeV

(right panel). The values from the standard method are positive overthefull



Nch



range,andare larger atlower



Nch



orinthe

higher pT range.Thisbehaviour suggeststhat the sc2,3

{

4

}

values

fromthe standard methodin pp collisions, including those from Ref. [19],are strongly influenced by non-flow effects inall



Nch



and pT ranges [16]. In contrast, the values from the subevent

methods are negative over the full



Nch



range, and they are

slightlymorenegative atlowest



Nch



andalsomorenegative at

higher pT.Theresultsare consistentamongthe varioussubevent

methods for 0

.

3

<

pT

<

3 GeV.For the high pT region of 0

.

5

<

pT

<

5 GeV,resultsfromthetwo-subeventmethodare

systemati-callylowerthanthosefromthethree- andfour-subeventmethods, suggestingthatthetwo-subeventmethodmaybeaffectedby neg-ative non-flow contributions. Such negative non-flow correlation hasbeenobservedinaPythia8 calculation [16].

Themiddlerowof Fig.1 showssc2,3

{

4

}

from p+Pbcollisions.

At



Nch

 >

140,the valuesare negative andconsistentamong all

fourmethods,reflectinggenuinelong-rangecollectivecorrelations. At



Nch

 <

140, the values are different between the standard

method and the subevent methods. The sc2,3

{

4

}

from the

stan-dardmethodchangessignaround



Nch

 ∼

80 andremainspositive

atlower



Nch



,reflecting thecontributionfromnon-flow

correla-tions.Incontrast,thesc2,3

{

4

}

fromvarioussubeventmethodsare

negativeandconsistentwitheachotherat



Nch

 <

140,suggesting

thattheymainlyreflectthegenuinelong-rangecorrelations.

observed between thetwo-subevent andthree- or four-subevent methods at low



Nch



, butthese differencesdecrease and

disap-pear for



Nch

 >

100. Within the statistical uncertainties of the

measurement,nodifferencesareobservedbetweenthethree- and four-subevent methods. This comparison suggests that the two-subeventmethodmaynotbesufficienttorejectnon-flow correla-tionsfromdijetsinpp collisions,andmethodswiththreeormore subeventsarerequiredtosuppressthenon-flowcontributionover themeasured



Nch



range.

The middlerowofFig. 2showssc2,4

{

4

}

from p+Pbcollisions.

Significantdifferencesareobservedbetweenthestandardmethod and the subevent methods over the full



Nch



range. However,

nodifferencesareobservedamongthevarioussubeventmethods. These results suggest that the standard method is contaminated by large contributions from non-flow correlations at low



Nch



,

and thesecontributions may not vanish even atlarge



Nch



val-ues.All subeventmethods suggest anincrease of sc2,4

{

4

}

toward

lower



Nch



for



Nch

 <

40,whichmayreflectsomeresidual

non-flowcorrelationsinthisregion.

ThebottomrowofFig.2showssc2,4

{

4

}

fromPb+Pbcollisions.

Thesc2,4

{

4

}

valuesincreasegraduallywith



Nch



forallfour

meth-ods. This increase reflects the known fact that the v2 increases

with



Nch



inPb+Pbcollisions [37].The valuesfromthestandard

method are systematically larger than those from the subevent methods, andthis difference varies slowly with



Nch



,similar to

the behaviour observed in p+Pb collisions in the high



Nch



re-gion.

The resultsfortheasymmetric cumulantac2

{

3

}

are presented

inFig. 3.The top rowshowstheresults obtainedfromthe stan-dard, two-subevent, and three-subevent methods from pp colli-sions in 0

.

3

<

pT

<

3 GeV (left panel) and 0

.

5

<

pT

<

5 GeV

(right panel).Theresultsare positiveforall methods.The results from the standard method are much larger than those from the subeventmethods, consistentwiththe expectationthat the stan-dardmethodismoreaffectedbynon-flowcorrelationsfromdijets. Significantdifferencesarealsoobservedbetweenthetwo-subevent and three-subevent methods at low



Nch



, but these differences

decrease and disappear at



Nch

 >

100. The ac2

{

3

}

values from

the three-subeventmethodshow a slightincrease for



Nch

 <

40

but are nearly constant for



Nch

 >

40. This behaviour suggests

thatinthethree-subeventmethod,thenon-flowcontributionmay play some role at



Nch

 <

40, but is negligible for



Nch

 >

40.

Therefore, the ac2

{

3

}

from the three-subevent method supports

theexistence ofathree-particlelong-rangecollective flowthat is nearly independent of



Nch



in pp collisions, consistentwith the



Nch



-independentbehaviourofv2 and v4 observedpreviouslyin

thetwo-particlecorrelationanalysis [7].

ThemiddleandbottomrowsofFig.3showac2

{

3

}

from p+Pb

and Pb+Pb collisions, respectively. The ac2

{

3

}

values from the

450 The ATLAS Collaboration / Physics Letters B 789 (2019) 444–471

Fig. 1. The symmetriccumulantsc2,3{4}asafunctionofNchfor0.3<pT<3 GeV(leftpanels)and0.5<pT<5 GeV (rightpanels)obtainedfor pp collisions(toprow),

p+Pbcollisions(middlerow)andlow-multiplicityPb+Pbcollisions(bottomrow).Ineachpanel,thesc2,3{4}isobtainedfromthestandardmethod(filledsymbol),the

two-subeventmethod(opencircles),three-subeventmethod(opensquares)andfour-subeventmethod(opendiamonds).Theerrorbarsandshadedboxesrepresentthe statisticalandsystematicuncertainties,respectively.



Nch

 ∼

200 inp+Pbcollisionsand



Nch

 ∼

80 inPb+Pbcollisions.

Inthesubevent methods,theinfluenceofnon-flow contributions isvery smallfor



Nch

 >

60 inboth collisionsystems,and

there-forethe



Nch



dependenceofac2

{

3

}

reflectsthe



Nch



dependence

ofthe v2 and v4. The ac2

{

3

}

valuesfrom thesubevent methods

increase with



Nch



, and the increase is stronger in Pb+Pb

colli-sions. This is consistent with previous observations that v2 and

v4increasewith



Nch



morestronglyinPb+Pbthaninp+Pb

colli-sions [17].

The valuesofsc2,4

{

4

}

andac2

{

3

}

,which areboth measuresof

correlations between v2 and v4, show significant differences

be-tween thestandardmethodandthesubevent methods,asshown in Figs. 2 and 3. The



Nch



dependence of these differences

de-creasesgraduallywith



Nch



,andisconsistentwithaninfluenceof

non-flowthatisexpectedtoscaleas1

/



Nch



.However,these

dif-ferencesseemtopersistfor



Nch

 >

200 inp+Pbcollisionsandfor



Nch

 >

150 inPb+Pbcollisions,whichisnotcompatiblewiththe

Fig. 2. The symmetriccumulantsc2,4{4}asafunctionofNchfor0.3<pT<3 GeV(leftpanels)and0.5<pT<5 GeV (rightpanels)obtainedfor pp collisions(top

row),p+Pbcollisions(middlerow)andlow-multiplicityPb+Pbcollisions(bottomrow).Ineachpanel,thesc2,4{4}isobtainedfromthestandardmethod(filledsymbol),

two-subeventmethod(opencircles),three-subeventmethod(opensquares)andfour-subeventmethod(opendiamonds).Theerrorbarsandshadedboxesrepresentthe statisticalandsystematicuncertainties,respectively.

large



Nch



mayarisefromlongitudinalflowdecorrelations [38,39],

whichhavebeenmeasuredbyCMS [40] andATLAS [41]. Decorre-lationeffectsarefoundto belargefor v4 andstronglycorrelated

with v2, andtherefore they are expected to reduce the sc2,4

{

4

}

andac2

{

3

}

in the subeventmethod. Therefore, the observed

dif-ferencesbetweenthe standardmethod andsubevent method re-flectthecombinedcontributionfromnon-flowcorrelations,which dominates in the low



Nch



region, and decorrelation, which is

moreimportant atlarge



Nch



(see furtherdiscussion inthe

Ap-pendixB).

The results presented above suggest that the three-subevent method is sufficient to suppress mostof the non-flow effects. It isthereforeusedasthedefaultmethodforthediscussionbelow. 8.2. Comparisonbetweencollisionsystems

Fig. 4 shows a direct comparison of cumulants for the three collisionsystems.Thethreepanelsinthetoprowshowtheresults forsc2,3

{

4

}

,sc2,4

{

4

}

andac2

{

3

}

,respectively,for0

.

3

<

pT

<

3 GeV.

These results support the existence of a negative correlation be-tween v2 and v3 and a positive correlation between v2 and v4.

452 The ATLAS Collaboration / Physics Letters B 789 (2019) 444–471

Fig. 3. The asymmetriccumulantac2{3}asafunctionofNchfor0.3<pT<3 GeV(leftpanels)and0.5<pT<5 GeV(rightpanels)obtainedforpp collisions(toprow),p+Pb

collisions(middlerow)andlow-multiplicityPb+Pbcollisions(bottomrow).Ineachpanel,theac2{3}isobtainedfromthestandardmethod(filledsymbol),two-subevent method(opencircles),andthree-subeventmethod(opensquares).Theerrorbarsandshadedboxesrepresentthestatisticalandsystematicuncertainties,respectively.

Such correlation patterns havepreviously been observed inlarge collisionsystems [42–44],butarenowconfirmedalsointhesmall collisionsystems,oncenon-floweffectsareadequatelysuppressed. Inthemultiplicityrangecoveredbythe pp collisions,



Nch

 <

150,

theresultsforsymmetriccumulantssc2,3

{

4

}

andsc2,4

{

4

}

are

sim-ilaramongthethree systems.Intherange



Nch

 >

150,

|

sc2,3

{

4

}|

andsc2,4

{

4

}

arelargerinPb+Pbthaninp+Pbcollisions.Theresults

for ac2

{

3

}

are similar among the three systems at



Nch

 <

100,

buttheydeviatefromeachother athigher



Nch



.The pp dataare

approximatelyconstantordecreaseslightlywith



Nch



,whilethe

p+PbandPb+Pb data show significant increasesas a functionof



Nch



.The bottomrowshowstheresultsforthehigher pT range

of0

.

5

<

pT

<

5 GeV,wheresimilartrendsareobserved.

Fig. 5 shows the results for normalised cumulants, nsc2,3

{

4

}

,

nsc2,4

{

4

}

and nac2

{

3

}

, compared among the three systems. The

normalised cumulants generally show a much weaker



Nch



de-pendence at



Nch

 >

100, where the statistical uncertainties are

small.Thisbehaviourimpliesthatthestrong



Nch



dependenceof

thescn,m{4

}

andac2

{

3

}

valuesreflectsthe



Nch



dependenceofthe

vn values,andthesedependencesare removedinthe normalised cumulants.Thenormalised cumulantsarealsosimilaramong dif-ferent collisionsystems atlarge



Nch



,althoughsome differences

Fig. 4. The Nchdependenceofsc2,3{4}(leftpanels),sc2,4{4}(middlepanels)andac2{3}(rightpanels)in0.3<pT<3 GeV(toprow)and0.5<pT<5 GeV(bottomrow)

obtainedfor pp collisions(solidcircles),p+Pbcollisions(opencircles)andlow-multiplicityPb+Pbcollisions(opensquares).Theerrorbarsandshadedboxesrepresentthe statisticalandsystematicuncertainties,respectively.

attherelativelevelof20–30%areobservedforsmaller



Nch



.The

onlyexception isnsc2,3

{

4

}

,whose valuesinthe pp collisions are

verydifferentfromthoseinp+PbandPb+Pbcollisions.Incontrast, the sc2,3

{

4

}

values in Fig. 4 are close among different systems.

Thissuggeststhatthe



v2₃



valuesfromthetemplatefitmethod [7] maybesignificantlyunderestimated.AspointedoutinRef. [7] and emphasised in Appendix A, the template fit method, and other methods based on peripheral subtraction in general [5,15], tend tounderestimatetheoddflow harmonics,duetothepresenceof alargeaway-side peak at

φ

∼

π

inthetwo-particle correlation function.Thecomparisonofsc2,3

{

4

}

andnsc2,3

{

4

}

amongdifferent

collision systems provides indirect evidenceof this underestima-tionof



v2₃



.

Fig.5 showsthat thenormalisedcumulantsareconsistent be-tween0

.

3

<

pT

<

3 GeVand0

.

5

<

pT

<

5 GeV.Ontheotherhand,

themagnitudesofthecumulantsinFig.4differby alargefactor betweenthetwopTranges:aboutafactorofthreeforsc2,3

{

4

}

and

sc2,4

{

4

}

,anda factoroftwoforac2

{

3

}

.Theseresultssuggest that

the pT dependenceofsc2,3

{

4

}

,sc2,4

{

4

}

andac2

{

3

}

largelyreflects

thepT dependenceofthevn atthesingle-particlelevel.

9. Discussion

Three- and four-particle cumulants involving correlations be-tweentwo harmonicsofdifferentorder vn andvm are measured in

√

s

=

13 TeV pp,

√

sNN

=

5

.

02 TeV p+Pb,andlow-multiplicity

√

sNN

=

2

.

76 TeV Pb+Pb collisions with the ATLAS detector at

theLHC,with totalintegratedluminosities of0.9 pb−1, 28 nb−1, and 7 μb−1, respectively. The correlation between vn and vm is studied using four-particle symmetric cumulants, sc2,3

{

4

}

and

sc2,4

{

4

}

,and the three-particleasymmetric cumulant ac2

{

3

}

.The

symmetriccumulants scn,m{4

}

=



v_n2v2_m



−



v2_n

 

v_m2



probe the cor-relation ofthe flow magnitudes,while the asymmetric cumulant ac2

{

3

} =



v2₂v4cos 4

(

2

− 

4

)



issensitivetocorrelationsinvolving

both the flowmagnitude vn andflow phase



n.They are calcu-lated using the standard cumulant method,as well asthe two-, three- and four-subevent methods to suppress non-flow effects. Thefinalresultsarepresentedasafunctionoftheaveragenumber ofchargedparticleswithpT

>

0

.

4 GeV,



Nch



.

Significant differences are observed between the standard method and the subevent methods over the full



Nch



range in

pp collisions, as well as over the low



Nch



range in p+Pb and

Pb+Pbcollisions. The differencesarelarger forparticles athigher pT oratsmaller



Nch



.When analysedwiththestandard method

in pp collisions,thisbehaviouriscompatiblewiththedominance ofthe non-flowcorrelationsratherthan thelong-rangecollective flow correlations. Systematic, but much smaller, differences are alsoobserved inthelow



Nch



regionbetweenthetwo-subevent

methodandthree- orfour-subeventmethods,whichindicatethat the two-subevent method may still be affected by correlations arising fromjets. Onthe other handno differencesare observed between the three-subevent and four-subevent methods, within experimentaluncertainties,suggestingthatmethodswiththreeor moresubeventsaresufficienttorejectnon-flow correlationsfrom jets. Therefore,thethree-subevent methodisusedto presentthe mainresultsinthisanalysis.

The three-subevent methodprovides a measurement of nega-tive sc2,3

{

4

}

and positive sc2,4

{

4

}

andac2

{

3

}

over nearly thefull



Nch



range and in all three collision systems. These results

454 The ATLAS Collaboration / Physics Letters B 789 (2019) 444–471

Fig. 5. The Nchdependenceofnsc2,3{4}(leftpanels),nsc2,4{4}(middlepanels)andnac2{3}(rightpanels)in0.3<pT<3 GeV(toprow)and0.5<pT<5 GeV(bottom

row)obtainedforpp collisions(solidcircles),p+Pbcollisions(opencircles)andlow-multiplicityPb+Pbcollisions(opensquares).Theerrorbarsandshadedboxesrepresent thestatisticalandsystematicuncertainties,respectively.

correlationbetweenv2andv4.Suchcorrelationpatternshave

pre-viously beenobserved inlarge collision systems [42–44], butare now confirmed in small collision systems, once non-flow effects areadequatelysuppressed.Thevaluesofsc2,3

{

4

}

andsc2,4

{

4

}

are

consistent in pp and p+Pb collisions over the same



Nch



range,

buttheir magnitudesatlarge



Nch



are muchsmaller thanthose

forPb+Pb collisions.The valuesof ac2

{

3

}

are similar atvery low



Nch



among the three systems, but are very different at large



Nch



.Ontheotherhand,afterscalingbythe



v2_n



estimatedfrom atwo-particle analysis [7,14], theresultingnormalisedcumulants nsc2,3

{

4

}

,nsc2,4

{

4

}

andnac2

{

3

}

showamuchweakerdependence

on



Nch



,and their valuesare much closerto each other among

the three systems. The magnitudes of the normalised cumulants are also similar to each other for 0

.

5

<

pT

<

5 GeV as well as

0

.

3

<

pT

<

3 GeV. This suggests that the



Nch



, pT and system

dependenceofthesc2,3

{

4

}

,sc2,4

{

4

}

andac2

{

3

}

reflectmostlythe



Nch



,pTandsystemdependenceof



v2n



,buttherelativestrengths ofthecorrelationsaresimilarforthethreecollisionsystems.

The new results obtained with the subevent cumulant tech-nique provide further evidence that the ridge is indeed a long-rangecollectivephenomenoninvolvingmanyparticlesdistributed across a broad rapidity interval. The similarity between differ-entcollisionsystemsfornsc2,3

{

4

}

,nsc2,4

{

4

}

andnac2

{

3

}

,andthe

weakdependenceoftheseobservablesonthepTrangeand



Nch



,

largelyfree fromnon-flow effects,provideanimportantinput to-wardsunderstandingthespace–timedynamicsandtheproperties of the medium created in small collision systems. These results provideinputstodistinguishbetweenmodelsbasedoninitial-state momentumcorrelationsandmodelsbasedonfinal-state hydrody-namics.

Acknowledgements

We thank CERN forthe very successful operation ofthe LHC, as well asthe supportstaff fromour institutions withoutwhom ATLAScouldnotbeoperatedefficiently.

WeacknowledgethesupportofANPCyT,Argentina;YerPhI, Ar-menia; ARC, Australia; BMWFW and FWF, Austria; ANAS, Azer-baijan; SSTC, Belarus; CNPq and FAPESP,Brazil; NSERC, NRC and CFI, Canada; CERN; CONICYT,Chile; CAS,MOSTand NSFC,China; COLCIENCIAS, Colombia; MSMT CR, MPO CR and VSC CR, Czech Republic;DNRFandDNSRC,Denmark;IN2P3-CNRS,CEA-DRF/IRFU, France; SRNSFG, Georgia; BMBF, HGF, andMPG, Germany; GSRT, Greece; RGC,Hong KongSAR, China;ISFandBenoziyoCenter, Is-rael; INFN, Italy; MEXT and JSPS, Japan; CNRST, Morocco; NWO, Netherlands; RCN,Norway; MNiSW andNCN, Poland; FCT, Portu-gal;MNE/IFA,Romania; MESofRussiaandNRCKI,Russian Feder-ation;JINR;MESTD,Serbia;MSSR,Slovakia;ARRSandMIZŠ, Slove-nia; DST/NRF, South Africa; MINECO, Spain; SRC and Wallenberg Foundation, Sweden;SERI,SNSFandCantonsofBernandGeneva, Switzerland;MOST,Taiwan;TAEK, Turkey;STFC,UnitedKingdom; DOE and NSF, United States of America. In addition, individ-ual groupsandmembers havereceived support fromBCKDF, the Canada Council, CANARIE,CRC, Compute Canada,FQRNT, andthe OntarioInnovation Trust,Canada; EPLANET,ERC, ERDF,FP7, Hori-zon 2020and Marie Skłodowska-Curie Actions, European Union; Investissements d’Avenir Labex and Idex, ANR, Région Auvergne andFondationPartagerleSavoir,France;DFGandAvHFoundation, Germany;Herakleitos,ThalesandAristeiaprogrammesco-financed by EU-ESFandtheGreek NSRF;BSF,GIFandMinerva,Israel;BRF, Norway; CERCA Programme Generalitat de Catalunya,Generalitat

Fig. 6. The valuesofvn{2,tmp}2obtainedfollowingthetemplatefitproceduregiveninEq. (14) [7] inpp collisionsforn=2 (leftpanel),n=3 (middlepanel)andn=4

(rightpanel).Ineachpanel,thevaluesarecalculatedforthreeperipheralNrecch intervals:N rec ch <20,N

rec

ch <10 and10≤N rec

ch <20.Onlystatisticaluncertaintiesareshown.

Valenciana,Spain;theRoyalSocietyandLeverhulmeTrust,United Kingdom.

The crucialcomputing support fromall WLCG partners is ac-knowledged gratefully, in particular from CERN, the ATLAS Tier-1 facilities at TRIUMF (Canada), NDGF (Denmark, Norway, Swe-den),CC-IN2P3(France),KIT/GridKA(Germany),INFN-CNAF(Italy), NL-T1(Netherlands),PIC(Spain),ASGC(Taiwan),RAL(UK)andBNL (USA),theTier-2facilitiesworldwideandlargenon-WLCGresource providers.Majorcontributorsofcomputingresources arelistedin Ref. [45].

Appendix A. Improvementtothetemplatefitprocedure

In order to separate the long-range ridge from other non-flowsources, especiallydijets,the ATLAS Collaborationdeveloped atemplatefittingproceduredescribedinRefs. [7,14].Thefirststep istoconstructa

φ

distributionofparticlepairswithlarge pseu-dorapidityseparation

|

η

|

>

2,theso-called“per-trigger”particle yield, Y

(φ)

, for a given N_chrec range. The

|

η

|

>

2 requirement suppresses the intra-jet and other short-range correlations, and in small collision systems the resulting Y

(φ)

distributions are knowntobedominatedbyaway-sidejetcorrelations [4,5,14].This away-sidenon-flowcomponentispeakedat

φ

∼

π

,andleadsto asignificantbiasintheflowcoefficientsvn,especiallyfortheodd harmonics.

Tosubtracttheaway-sidejetcorrelations,themeasuredY

(φ)

distribution in a given Nrec

ch interval is assumed to be a sum

ofa scaled “peripheral” distribution Y

(φ)

peri, obtainedfor low-multiplicity events Nrec

ch

<

20,anda constantpedestal modulated

bycos

(

n

φ)

forn

≥

2 [7,14]: Y

(φ)

=

F Y

(φ)

peri

+

Gtmp



1

+

2 ∞



n=2 vn

{

2

,

tmp

}

2cos n

φ

.

(14)

The scale factor F and pedestal Gtmp are fixed by the fit, and vn{2

,

tmp

}

arecalculated froma Fouriertransform. Onthe other hand,both Y

(φ)

andY

(φ)

peri contain a dijetcomponent and flowcomponent:

Y

(φ)

=

Y

(φ)

cent_jet

+

Gcent



1

+

2 ∞



n=2 vn

{

2

}

2cos n

φ

,

(15)

Y

(φ)

peri

=

Y

(φ)

peri_jet

+

Gperi



1

+

2 ∞



n=2 vn

{

2

,

peri

}

2cos n

φ

.

(16) Withtheassumptionthat theshapeofthedijetcomponentis in-dependent of Nrec_ch, and the magnitudes of the dijet components arerelatedbythescalefactorF :Y

(φ)

_jetcent

=

F Y

(φ)

peri_jet ,Eq. (14) canbewrittenas:

Y

(φ)

=

Y

(φ)

cent_jet

+ (

Gtmp

+

F Gperi

)

+

2 ∞



n=2



Gtmpvn

{

2

,

tmp

}

2

+

F Gperivn

{

2

,

peri

}

2



×

cos n

φ.

Comparing with Eqs. (15) and (16), one obtains Gcent

₌

_Gtmp

₊

F Gperiandthefollowingrelation:

vn

{

2

}

2

=

vn

{

2

,

tmp

}

2

−

F Gperi Gcent



vn

{

2

,

tmp

}

2

−

vn

{

2

,

peri

}

2



,

which shows that vn{2

,

tmp

}

from the template fit differs from the true vn{2

}

by a correction term that vanishes ifand only if vn{2

}

is independent of Nrec_ch. Since the true flow harmonics in the peripheral interval vn{2

,

peri

}

are unknown in principle, the correction is applied starting from the third-lowest Nrec

ch interval

(40

≤

Nrec_ch

<

60) inthis analysis, by using vn{2

,

tmp

}

of the sec-ond Nrec_ch interval (20

≤

N_chrec

<

40) as an estimate of the true flow harmonics.Sincethe non-flow contributionprimarily affects the oddharmonics, the v3

{

2

,

tmp

}

2 maybecome negative in the

firstfew Nrec_ch intervalsin pp collisions.Insuchcases,the correc-tionstartsfromthesecond Nrec_ch intervalwithpositive v3

{

2

,

tmp

}

2

(60

≤

Nrec

ch

<

80)byusingv3

{

2

,

tmp

}

fromtheprevious Nrecch

inter-val(40

≤

Nrec_ch

<

60).

One important feature of the template fit analysis is the as-sumption that the dijet component Y

(φ)

jet is independent of



Nch



. In Ref. [7], the uncertainty associated with this

assump-tion is studied by changing the default peripheral interval from Nrec_ch

<

20 to Nrec_ch

<

10 and10

≤

Nrec_ch

<

20.Itwas foundthat the vn{2

,

tmp

}

valuesarerelativelyinsensitivetothechoiceof periph-eralintervalforn

=

2 andn

=

4,butthesensitivityismuchlarger for n

=

3. This finding is reproduced in Fig. 6 for pp collisions, whichshowsthatthev3

{

2

,

tmp

}

2valuesobtainedviaEq. (14)

dif-fersubstantiallyforthedifferentNrec_ch ranges.

Inadditiontothetemplatefitwithandwithouttheabove men-tioned correction procedure, the ATLAS and CMS collaborations

456 The ATLAS Collaboration / Physics Letters B 789 (2019) 444–471

Fig. 7. The v2(leftcolumn),v3 (middlecolumn)and v4 (rightcolumn)obtainedfromtwo-particlecorrelationsin0.3<pT<3 GeVinpp (toprow),p+Pb(middlerow)

andPb+Pb(bottomrow)collisions.Ineachpanel,theyarecomparedbetweenthreemethods:directFouriertransformation(solidcircles),templatefit(opencircles)andthe improvedtemplatefit(opensquares).Theerrorbarsandshadedboxesrepresentthestatisticalandsystematicuncertainties,respectively.

alsocalculateddirectlythe vn{2

}

valuesviaaFouriertransformof theY

(φ)

distributionwithoutdijetsubtraction [7,19].The differ-encesbetweenthedirectFouriertransformandtemplatefitreflect mainlythe away-sidejet contributionsubtracted by thetemplate fitprocedure,andthereforegivea senseofthemagnitudeof un-known systematic uncertainties associated with the template fit procedure. If these differences are too large, the vn

{

2

,

tmp

}

val-uesmaybesensitivetothesystematiceffectsassociatedwiththe assumptionthattheshapeofY

(φ)

jetisindependentofNrec_ch.

Fig.7 comparesthe vn{2

}

in 0

.

3

<

pT

<

3 GeV obtainedfrom

Y

(φ)

usingthreemethods:a directFouriertransform(solid cir-cles), a template fit (open circles) and a template fit corrected for the bias (open squares), as described above. The systematic uncertainties forthe template fit results are nearly the same as those fromRef. [7]. Fig. 7showsthat the changes introduced by

the correction procedure described above are small in all cases andforallharmonics.Thevaluesoftheeven-orderharmonics, v2

and v4, are also quite similar to those obtained from the direct

Fourier transformation, reflecting the fact that the dijet correla-tions havevery little influence on the even-orderharmonics. On the other hand,significant differences are observed between the direct Fouriertransformandtemplate fitfor v3,especially inthe

pp collisions, dueto the influence ofY

(φ)

jet,a trend observed

anddiscussedpreviouslyinRefs. [7,15].Thetemplatefitprocedure isabletosubtractthedijetcorrelationsandchangethesignofv3,

butalsointroducesalargeuncertaintyassociatedwiththe proce-dure. AsdiscussedinSection8.2,thebehaviourofthesymmetric cumulantssc2,3

{

4

}

inFig.4andnormalisedcumulantsnsc2,3

{

4

}

in

Fig.5inpp collisions,suggestthatthev3valuesfromthetemplate