UvA-DARE (Digital Academic Repository) Two-particle azimuthal correlations in photonuclear ultraperipheral Pb-Pb collisions at 5.02 TeV with ATLAS

Two-particle long-range azimuthal correlations are measured in photonuclear collisions using 1 . 7 nb − 1 of 5.02 TeV Pb + Pb collision data collected by the ATLAS experiment at the CERN Large Hadron Collider. Candidate events are selected using a dedicated high-multiplicity photonuclear event trigger, a combination of information from the zero-degree calorimeters and forward calorimeters, and from pseudorapidity gaps constructed using calorimeter energy clusters and charged-particle tracks. Distributions of event properties are compared between data and Monte Carlo simulations of photonuclear processes. Two-particle correlation functions are formed using charged-particle tracks in the selected events, and a template-ﬁtting method is employed to subtract the nonﬂow contribution to the correlation. Signiﬁcant nonzero values of the second-and third-order ﬂow coefﬁcients are observed and presented as a function of charged-particle multiplicity and transverse momentum. The results are compared with ﬂow coefﬁcients obtained in proton-proton and proton-lead collisions in similar multiplicity ranges, and with theoretical expectations. The unique initial conditions present in this measurement provide a new way to probe the origin of the collective signatures previously observed only in hadronic collisions.


Introduction
In ultrarelativistic collisions of lead nuclei at the Large Hadron Collider (LHC), the typical processes studied are those for which the nuclei have an impact parameter less than twice their radius ( 2  ).Such lead-lead (Pb+Pb) collisions are understood to create a quark-gluon plasma and result in a large number of particles in the final state which participate in collective motion as a consequence of the plasma evolution [1][2][3].In addition to the particles produced in Pb+Pb collisions, those produced in high-energy proton-proton ( ) and proton-lead (+Pb) collisions also exhibit a collective behavior which manifests as an event-wide azimuthal variation persisting broadly in pseudorapidity, initially observed as a 'ridge' [4][5][6][7].This behavior is characterized in terms of nonzero single-particle azimuthal anisotropies, given by  th -order Fourier coefficients, and here referred to as flow coefficients   .Nonzero   values have also been observed in significantly lower-energy (19.6-200GeV) +Au collisions at the Relativistic Heavy Ion Collider [8].A natural question is whether such signatures persist in even smaller collision systems [9] and, if so, how this may influence the interpretation of these signatures in   or +Pb collisions.
In the prevailing paradigm, these observed anisotropies arise from the creation of a miniature region of quark-gluon plasma [10,11], in which a hydrodynamic-like expansion of the system converts spatial nonuniformities in the initial state of the system into momentum-space anisotropies of the final-state particles.However, momentum correlations already present in the initial state of the collision may also persist into the final state [12].The relative importance of these two explanations can be tested in collision systems where one or both of the 'beams' has a significantly simpler initial state.Recently, studies were performed in archived  +  − collision data at √  = 91 GeV from the ALEPH detector [13] and in archived   collision data at √  = 316 GeV from the ZEUS detector [14] with a deep-inelastic-scattering selection are selected, and are distinguished from peripheral hadronic Pb+Pb events and other background events, using the topology of the distribution of particles in the event as measured in the zero-degree, forward, and barrel calorimeters, as well as the tracking systems through the reconstruction of pseudorapidity gaps [32,33].Properties of the selected events are compared with the expectations from Monte Carlo (MC) simulations of photonuclear processes.
Two-particle correlations as a function of relative separation in azimuth (Δ) and pseudorapidity 1 (Δ) are constructed for different selections of event charged-particle multiplicity and charged-particle kinematics.
The non-flow contributions to the two-particle correlations (for example, jet correlations which give rise to a Δ correlation structure) are suppressed by studying the correlations at large Δ, and the residual non-flow contribution is subtracted via the template-fitting method used extensively in prior ATLAS measurements [4,5,34].In the template method, the correlation in high-multiplicity events is described as a combination of the correlation in lower-multiplicity events plus a component modulated by cos(2Δ) (and higher order) Fourier terms.The template method makes particular assumptions about how the non-flow component evolves with multiplicity.Although there are differences between the system explored in this measurement and those in previous two-particle correlation measurements, the sensitivity to the assumptions of the template method can be tested within the standard approach.A test of the template method in simulated photonuclear events which do not include flow or initial-state correlation mechanisms is performed in Section 7. The resulting magnitudes of the two-particle correlations are interpreted as arising from the product of global  2 and  3 values for individual particles, and are reported as a function of the reconstructed charged-particle multiplicity ( rec ch ) and transverse momentum ( T ).The results are compared with other small collision systems at the LHC, and theoretical expectations from initial-and final-state physics mechanisms are discussed.

ATLAS detector and data sample
The ATLAS detector [35] covers nearly the entire solid angle around the collision point.The detector systems used for the measurements presented in this paper include the inner tracking detector, the electromagnetic (EM) and hadronic calorimeters, the zero-degree calorimeter (ZDC), and the trigger and data acquisition systems.The detector halves at positive and negative -values are referred to as the A and C sides, respectively.
The inner-detector system is immersed in a 2 T axial magnetic field and provides charged-particle tracking in the range || < 2.5.The high-granularity silicon pixel detector covers the vertex region and typically provides four measurements per track.The innermost layer, the insertable B-layer [36][37][38], has been operating as a part of the silicon pixel detector since 2015.The pixel detector is followed by the silicon microstrip tracker (SCT) which usually provides four measurement points per track.These silicon detectors are complemented by the transition radiation tracker, which enables radially extended track reconstruction up to || = 2.0.
The calorimeter system covers the pseudorapidity range || < 4.9.Within the region || < 3.2, electromagnetic calorimetry is provided by barrel and endcap high-granularity lead/liquid-argon (LAr) calorimeters, with an additional thin LAr presampler covering || < 1.8.Hadronic calorimetry is provided by a steel/scintillator-tile calorimeter, segmented into three barrel structures within || < 1.7, and two copper/LAr hadronic endcap calorimeters.The solid angle coverage is extended with the forward calorimeter (FCal), composed of copper/LAr and tungsten/LAr modules optimized for EM and hadronic measurements, respectively.The FCal detectors cover the regions 3.2 < || < 4.9.
The minimum-bias trigger scintillator (MBTS) detects charged particles over 2.07 < || < 3.86 using two hodoscopes of 12 counters positioned at  = ±3.6 m.The ZDC consists of EM and hadronic sections and plays a key role in identifying UPC events in heavy-ion collisions by primarily detecting neutrons resulting from the breakup of one or both nuclei.The ZDC modules are located at  = ±140 m from the IP.They measure neutral particles at pseudorapidities || ≥ 8.3 and consist of layers of alternating quartz rods and tungsten plates.
A two-level trigger system [39] is used to select events.The first-level trigger (L1) is implemented in hardware and uses a subset of the detector information to restrict the accepted rate to at most 100 kHz.This is followed by a software-based high-level trigger (HLT) stage that reduces the accepted event rate to 1-4 kHz depending on the data-taking conditions during 2018 Pb+Pb operations.
The measurements presented in this paper were performed using the √  NN = 5.02 TeV Pb+Pb dataset collected with a variety of triggers in 2018, with a total integrated luminosity of 1.7 nb −1 .Photonuclear candidate events were first selected by the trigger by requiring one ZDC side (referred to as the Pb-going side) to have a minimum amount of energy at L1,  > 1 TeV, consistent with the presence of one or more neutrons.The other side (referred to as the photon-going side) was required to have an energy below a maximum-energy cutoff,  < 1 TeV, consistent with no neutrons.The cut value of  = 1 TeV is several multiples of the energy resolution away from the single-neutron peak at  = 2.5 TeV [24].Thus, the selected topology is referred to as "0nXn" in the figures in this paper.Events were also required to satisfy an upper bound of 200 GeV on the total transverse energy deposited in the calorimeter at L1, for further rejection of hadronic Pb+Pb events.
After these requirements at L1, events had to pass either a minimum-bias (MB) trigger or one of several high-multiplicity triggers (HMT) with further requirements at the HLT stage.The MB trigger was defined by requiring at least one reconstructed online track with  T > 0.4 GeV.The HMTs were defined by requiring a larger number of tracks originating from the same vertex.The primary vertex, defined as the one with the largest scalar sum of  2 T values of associated tracks, was used.An HMT with a threshold of 15 tracks was defined by also requiring a hit in the MBTS at L1. Two HMTs with thresholds of 25 and 35 tracks were defined by requiring a minimum amount of energy in the calorimeter at L1. Finally, all HMTs had an additional HLT requirement of less than 5 GeV of transverse energy in the photon-going FCal, which rejected a large fraction of the 0nXn peripheral Pb+Pb background events.The 35-track HMT sampled the full luminosity of 1.7 nb −1 during data-taking.The 25-track HMT, 15-track HMT, and the MB trigger were configured with increasingly larger prescale factors, sampling 1.6 nb −1 , 0.13 nb −1 , and 1.0 b −1 of data respectively.

Monte Carlo simulation
Monte Carlo (MC) simulations of the relevant physics processes were used to understand the performance of the detector, and provide distributions to be compared with the data.For all the generators and configurations listed below, the simulated events were passed through a G 4 simulation [40,41] of the detector and reconstructed under the same conditions as the data.
A sample of one million peripheral Pb+Pb events was generated using H v1.383 [42] with impact parameters in the range 10 <  < 20 fm.These events were used to determine the charged-particle tracking efficiencies for the measurement, and to model the event properties of low-multiplicity, hadronic Pb+Pb events.
Several generators were used to simulate photonuclear interactions.For each of the cases below, events were generated with different minimum requirements on  rec ch to provide good statistical coverage over the  rec ch range accessed in data.First, the distribution of photon flux for 208 Pb beams at the LHC was calculated using STAR [43].The flux distribution was passed to a multipurpose generator based on the Dual Parton Model (DPM) and referred to as DPMJET-III [44,45], which simulated direct and resolved photon-lead (+Pb) interactions at the generator level, followed by a full G 4 simulation of the ATLAS detector.Thirteen million +Pb events were simulated in this way.Additionally, the flux from STAR was used to simulate two million + events using DPMJET-III, where the energy of the proton was set equal to the energy per nucleon in the Pb beams.Finally, for an alternative description of these processes, P 8.240 [46] configured with the NNPDF23LO parton distribution functions [47] and A14 set of tuned parameters [48] was used to generate twelve million + events.The photon flux in P was reweighted at the event level to match that calculated by STAR , and the simulation was configured to include both the direct and resolved photon interactions.

Reconstruction and event selection
Given the low particle multiplicities of UPC events, the charged-particle track and calorimeter energycluster reconstruction procedures follow those optimized for   data-taking [49,50].Reconstructed charged-particle tracks are used in the analysis if they satisfy quality criteria as outlined in Ref. [51].These include criteria for the minimum number of hits in the pixel and SCT detectors.Tracks are further required to have  T > 0.4 GeV, || < 2.5, and distances of closest approach to the primary reconstructed vertex in the longitudinal and transverse directions of less than 1.5 mm each.Clusters in the range || < 4.9 are constructed from topologically connected groups of calorimeter cells [52].They are required to have  T > 0.2 GeV and to meet the significance criteria for the measured energy as outlined in Ref. [32] to suppress the contribution from electronic noise fluctuations.
Events are required to have a reconstructed primary vertex within || < 90 mm.Events compatible with multiple in-time interactions are tagged by the presence of multiple reconstructed vertices.These comprise less than 1% of the sample and are rejected.Events must pass a version of the ZDC-based energy selection described in Section 2 which uses the same threshold but an improved calibration for the ZDC energies determined after data-taking.Events are further characterized by their charged-particle multiplicity,  rec ch , which is defined as in previous ATLAS measurements of correlations in small systems [4,5,7] as the total number of reconstructed tracks, without efficiency correction, with  T > 0.4 GeV and || < 2.5.For the correlation analysis described in Section 5, events with a given  rec ch were included in the analysis if they were selected by any of the HMTs which was more than 80% efficient for events with this  rec ch value.For the study of event properties in Section 4.2, each  rec ch range was only populated by the highest-sampled-luminosity HMT which was more than 99% efficient over the entire  rec ch range.In either case, events with  rec ch < 15 were populated only by the MB trigger.Reconstructed pseudorapidity gap quantities, constructed using the tracks and clusters in each event, are used to distinguish between different physics processes such as photonuclear collisions, low-activity (peripheral) hadronic Pb+Pb collisions, and  →  processes.The requirement of a rapidity gap above a minimum value in the photon-going direction can efficiently remove peripheral Pb+Pb events.Rather than the traditional pseudorapidity gap quantity [32], which determines the pseudorapidity difference between the edge of the detector and closest particle, an alternative 'sum-of-gaps' definition is used, which adds together contiguous gaps separated by particle production concentrated in a narrow pseudorapidity regions.This alternative definition is used to retain a large selection efficiency for resolved photon events which may break up a large gap with a hadronic fragment localized in pseudorapidity (see Figure 1).The quantities Σ  Δ and Σ  Δ correspond to the sum-of-gaps calculated in the photon-going and nucleus-going halves of the detector, respectively.They are calculated by first sorting the tracks and clusters in .The differences in  between adjacent particles, Δ, are included in the sum if they are larger than 0.5.The value of 0.5 was observed in simulation to retain good efficiency for resolved photon events.The position at  = 0 is treated as if it were a particle which ends the sum-of-gaps calculation.Thus, Σ  Δ and Σ  Δ range from 0 to 4.9.
In addition to the sum-of-gaps quantities defined above, the traditional edge gap quantity, constructed using clusters, Δ cluster edge , is defined as the pseudorapidity difference between the edge of the detector on the photon-going side, at  = +4.9, and the first reconstructed cluster.In the mixed-event construction discussed in Section 5, events are classified by their Δ cluster edge as a way to characterize their overall topology without using charged-particle tracks.ch , for the photonuclear event selection from multiple triggers, without prescale correction.Right: Distribution of event charged-particle multiplicity,  rec ch , for the photonuclear event selection (black) with corrections for trigger efficiency and prescale factors, and an alternative selection intended to select hadronic Pb+Pb events (red, see text).

Event properties
Figure 2 (left) shows the distribution of Σ  Δ and  rec ch values for events recorded with the MB trigger.The Σ  Δ can have a value of zero or fall into the range 0.5 ≤ Σ  Δ ≤ 4.9, with the lower value determined by the minimum gap size included in the sum.Events resulting from hadronic Pb+Pb collisions are at small Σ  Δ and have a broad multiplicity distribution, whereas events from photon-induced processes are at large Σ  Δ and have an  rec ch distribution that falls steeply.Figure 2 (right) shows the distribution of Σ  Δ and Σ  Δ values.Most hadronic Pb+Pb collisions have a broad pseudorapidity distribution of particles and hence have both small Σ  Δ and small Σ  Δ.In contrast, the photonuclear events have large Σ  Δ but smaller values of Σ  Δ.There is no significant yield of events with a large sum-of-gaps on both sides, which could signal the presence of photon-photon hadronic process backgrounds.The signal events for photonuclear collisions are defined by Σ  Δ > 2.5.No event selection is made on Σ  Δ.
Figure 3 summarizes the multiplicity distribution observed in data after the event selection described above.The  rec ch distribution for events passing the photonuclear event selection is shown in the left panel of Figure 3.The sawtooth pattern arises from the inclusion of the HMTs at different  rec ch thresholds as described above.In each  rec ch range, the selected photonuclear events have a steeply falling multiplicity distribution.The  rec ch distribution for photonuclear events, fully corrected for the different luminosities sampled by the triggers and the  rec ch -dependent trigger efficiency, is shown using black circles in the right panel of Figure 3.It is compared with the distribution from events with Σ  Δ < 1, which mostly selects Pb+Pb hadronic events, which are backgrounds in this analysis.The selected photonuclear events have a significantly more steeply falling multiplicity distribution.
Figure 4 compares the multiplicity and sum-of-gap distributions in data and simulation.The left panel of Figure 4 shows the total  rec ch distribution, corrected for the trigger efficiency and the different integrated luminosities sampled by the HMTs.This is compared with three distributions from MC generators, with  the same gap-based selection requirements as the data.The ZDC selection is applied to data but not to the MC samples as the generators do not model the nuclear-breakup processes relevant for forward neutron spectators.Including a ZDC requirement in MC events may thus impact the distributions if particle production in the detector is correlated with the nuclear fragmentation in the forward region.Despite the limitations in the modeling of forward neutrons, a generator-level check requiring neutrons in the ZDC acceptance was performed and found not to impact the level of agreement between data and DPMJET-III.
The distributions in P and DPMJET-III + are normalized to have the same integral as the data over the full  rec ch range.The models show good agreement with the data at low  rec ch , but systematically predict too low a relative yield at higher  rec ch .DPMJET-III +Pb does not describe the full distribution, but has been normalized to have the same integral as the data in  rec ch > 35 to highlight its good agreement in this region over many orders of magnitude.With this normalization, DPMJET-III +Pb systematically predicts too low a relative yield for  rec ch < 15.This comparison either suggests the presence of other, non-photonuclear, processes in data at such low  rec ch < 10 values, or points to the need for improved modeling of this region in the simulation.
The right panel of Figure 4 shows the reconstructed Σ  Δ distributions in data and simulation for selected events with  rec ch > 10, without the Σ  Δ > 2.5 requirement.Structures in the distributions correspond to transitions between detector subsystems and the change in the detector response as a function of .At large Σ  Δ values ≥ 2.5, the shape of the distribution in data is qualitatively similar to that in DPMJET-III +Pb and Pythia + simulation.However, the distributions in the simulated photonuclear events decrease at smaller Σ  Δ values, while the distribution in data rises.At low Σ  Δ, the shape in data is qualitatively similar to that in peripheral H Pb+Pb events.This comparison suggests that the trigger-selected events contain a mixture of peripheral Pb+Pb events and genuine photonuclear events, with the latter dominant at Σ  Δ > 2.5.The possible impact of residual peripheral Pb+Pb events in the set of selected events is discussed in Section 6.
Figure 5 compares the charged-particle pseudorapidity distribution,  ch /, in data and simulation.The left panel shows the  ch / in data, for charged particles with 0.4 <  T < 5 GeV, for multiple  rec ch selections in photonuclear events.The distributions are corrected for tracking efficiency on a per-track basis, which ranges from 0.7-0.9depending on track  and  T .To compare the relative shapes between  rec ch selections, the distributions are each normalized to have an integral of one.In all cases, the pseudorapidity distributions are strongly asymmetric, peaking at  = −2.5 (the nucleus-going direction) and then monotonically decreasing until  = +2.5 (the photon-going direction).This overall shape is even more asymmetric than that observed in, for example, +Pb collisions at the LHC [53].The  ch / shape remains similar over a wide range of multiplicities ( rec ch > 10).The selection with  rec ch < 5 results in an even more strongly asymmetric  ch / shape, suggesting that it may include photonuclear events at significantly lower energies or other processes.Additionally, as discussed above, the region  rec ch < 10 is not well described by simulation, and may contain a significant contribution from background events.On the basis of these observations, the two-particle correlation analysis is performed using events with  rec ch ≥ 10.Additionally, an upper bound of  rec ch = 60 is used to remove the region with very few events.The right panel in Figure 5 compares the pseudorapidity distribution of particles with 0.4 <  T < 5 GeV in data for  ch rec ≥ 10 with the four sets of simulated events, which have the same reconstruction-level selections as the data.While all three generators considered here systematically fail to predict a steep enough slope for the  ch / distribution, the distribution for DPMJET-III +Pb events is the most similar to that observed in data.In general, no single generator is able to quantitatively describe in detail all of the features observed in data.The results here may be used to benchmark future developments in the modeling of photonuclear collisions.The two-particle correlation analysis described in Section 5, due its data-driven nature, does not rely on a detailed description of photonuclear events in simulation and therefore its results are not sensitive to this mismodeling.
Since the DPMJET-III +Pb simulation and the data agree qualitatively in the  rec ch and Σ  Δ distributions, the simulated events provide an estimate of what ranges of photon energy in the laboratory frame,   , and photon-nucleon center-of-mass energy,  N , are selected for analysis by choosing various  rec ch or Σ  Δ values.These ranges are shown in Figure 6 as a function of  rec ch .In DPMJET-III, photonuclear events with a larger  rec ch have a larger average   and  N .At  rec ch = 10, the mean values are   = 26.8GeV and  N = 519 GeV, rising to   = 123 GeV and  N = 1.11TeV at  rec ch = 60.However, even a narrow range of  rec ch values selects events with a broad distribution of   and  N .At lower values of  rec ch < 10, the mean photon energy decreases rapidly and the range of photon energies increases, supporting the data-driven minimum value of  rec ch determined above.

Two-particle correlations
This section describes the two-particle correlation and template fit procedures used to extract the azimuthal anisotropies,   .These procedures have been used in previous analyses of two-particle correlations in small systems [4,5,34,54,55], and are summarized here.
Correlations between pairs of charged particles are reported as a function of the pair's relative separation in pseudorapidity, Δ =   −   , and azimuth, Δ =   −   .The labels  and  denote kinematic selections on the first and second particle, which are generally different.If two particles each meet selections  and , the pair is tabulated twice, once for each possible assignment.By convention, the distributions in Δ and Δ,  (Δ, Δ), are averaged over events by normalizing them by the total yield of particles in selection ,   , To correct for the charged-particle reconstruction efficiency, entries in the distribution are weighted by the product of the inverse efficiencies evaluated for the kinematics of the two particles.These efficiencies are derived from the H simulation described in Section 3. In this analysis, all charged-particle tracks within || < 2.5 are considered, except the tracks which lie in the region  > 4.9 − Σ  Δ (with positive  defined as the photon-going direction), if the sum-of-gaps extends into the tracker.There are two reasons to exclude these tracks.First, in simulations of photonuclear events, such particles are typically fragments of resolved photons, and not the bulk particles produced in the +Pb interaction which are being measured for possible collective signatures.Second, since these tracks are typically near the edge of the inner detector's acceptance, the |Δ| restriction described below would result in a small number of tracks being paired many times, and thus making an disproportionately large contribution to the correlation function.
The one-dimensional Δ normalized pair yields,  (Δ), are constructed by integrating over a restricted |Δ| range, typically 2.0 to 5.0, designed to suppress the contributions to the correlation from non-flow effects (such as jets, particle decays, and resonances), as: In principle, due to finite coverage and detection efficiency introduced by a real detector,  may contain acceptance effects that do not reflect physical particle pair correlations.To account for this, a mixed-event technique may be applied [56], in which particles  and  are taken from different events that have similar overall characteristics.The two events are required to have a ZDC signature on the same side (i.e.same directions of the photon and nucleus), and a similar charged-particle multiplicity, primary vertex position, and cluster edge gap Δ cluster edge .The mixed-event distribution  2  mixed /(ΔΔ) contains the effects of the pair acceptance but not the physical correlations.
The acceptance effects can be corrected for by dividing the same-event pair distribution by the mixed-event ones.The two-dimensional (2-D) correlation function corrected for acceptance effects is where   is the total yield of particles with selection ,  mixed is the yield in mixed events, and  mixed pair is the total integral of the mixed-event distribution  2  mixed /(ΔΔ).Without the mixing correction, the 2-D  (Δ, Δ) distributions have an artificial triangular shape in Δ, reflecting the convolution of the single-particle acceptances   ,   < 2.5.This effect is removed by the mixing correction, and the Δ dependence of structures in the 2-D  (Δ, Δ) distributions become more evident.
Correlation functions are constructed for various selections on  rec ch and   T , with the  T of particle  always in the range 0.4 <   T < 2 GeV.All events are given an equal statistical weight, and no corrections for the trigger efficiency or differences in sampled luminosity (for  rec ch ranges that span contributions from multiple triggers) are performed.Although the mixed-event correction is necessary to properly construct the two-dimensional correlation functions presented below, it was found to ultimately have only a minor effect on the projected, one-dimensional  (Δ) distributions, producing compatible results but with increased statistical uncertainty.Thus the mixing correction is not applied for the nominal  2 and  3 results.Instead, it is used in the determination of systematic uncertainties associated with the pair acceptance of the detectors as described in Section 6.
Examples of two-dimensional  (Δ, Δ) correlation functions are presented in Figure 7 and Figure 8.In Figure 7, particles  and  are required to have 0.4 <   T < 2 GeV , and two example  rec ch selections are shown.In Figure 8, correlation functions are presented for 20 <  rec ch < 60, with two example  T selections for particle  shown.The two-dimensional correlation functions have features which are broadly similar to those observed in   collisions.There is a localized 'near-side' peak at (Δ, Δ) ≈ (0, 0) from correlations between jet fragments, and an extended 'away-side' ridge at Δ ≈  which extends over a large Δ range from correlations between fragments of azimuthally opposite jets.A quantitative analysis and non-flow subtraction are necessary to discern if there are additional ridge structures, as detailed below.

Non-flow subtraction
To remove the contribution to the correlation function from non-flow effects that persist as long-range contributions on the away side (Δ ∼ ), a template fit procedure is employed.To perform the fit procedure, samples corresponding to events with low multiplicity (LM) and high multiplicity (HM) are selected.The shape of the non-flow contribution is assumed in the template procedure to be the same in the LM and HM samples.In this analysis, the LM events are chosen to have 15 ≤  rec ch ≤ 20.The  (Δ) in HM events is parameterized as the sum of an azimuthally modulated pedestal (which expresses the azimuthal anisotropy) and a non-flow component, as follows: Above,  LM (Δ) is the correlation function in LM events, which is parameterized as a truncated Fourier series up to the fourth order (with coefficients   for  = 0, . . ., 4).The values of  and , the three  , values, and the five parameters   , describing the LM reference, are free parameters in the fit, but  and  are constrained such that the integrals of both sides of Eq. ( 1) are the same.Modulation terms up to fourth order ( 2,2 ,  3,3 , and  4,4 ) are considered in the fit, in order to best describe the HM data.By fitting  LM (Δ) and  HM (Δ) simultaneously, the extracted uncertainty in , , and  , correctly accounts for the statistical uncertainty of both the LM and HM samples.An example of the simultaneous fit of the HM selection and LM reference is shown in Figure 9. Examples of the template fit in additional HM and   T selections are shown in Figure 10.In the bottom panels of Figs. 9 and 10, the -values are defined following the procedure described below.
The template fit is performed by minimizing the standard  2 between the data points and the functional form.However, the data points within the correlation functions contain nontrivial point-to-point correlations, since a single particle  may be used in combination with multiple particles of type .The minimum of the  2 statistic, when calculated in the traditional way, is found at the appropriate values of the fit parameters.However, the -value and the uncertainty in the parameter values, if also determined in the standard way, would be inaccurate.In order to properly account for these correlations and determine the parameter value uncertainties, a bootstrapping procedure was applied.Pseudo-experiments were generated by giving a random Poisson weight (with a mean of one) to individual events.These pseudo-experiments are independent sets of events drawn from the observed set of events with the same statistics and correlations as the real data.The parameters values extracted in the pseudo-experiments agreed with those in data and followed Gaussian distributions.The standard deviations of the parameters values extracted from the pseudo-experiments were used to set the statistical uncertainties in the final results, which are 10-30% larger than if they were estimated in the naïve way without taking the point-to-point correlations into account.The -values from the bootstrapping experiment are presented in the template fit plots that follow.
The full set of  2,2 and  3,3 values are shown in Figure 11.The filled points show the non-flow subtracted results according to the template fit procedure in Eq. ( 1).The open points show the  , values that would be obtained if instead a direct Fourier decomposition of  HM were performed, without any attempt to account for the non-flow contribution.For most of the selections considered here, the non-flow subtraction has a significant effect on the extracted  2,2 and  3,3 values.The resulting  , values are positive in all selections, with one exception: in the    T -dependent results with a single HM selection, the  2,2 value for 3 <  T < 5 GeV is negative.Additionally, the  2,2 value for 2 <  T < 3 GeV is significantly lower than that for 1.2 <  T < 2 GeV.In these selections, the  3,3 values also rise significantly.The template fits to these selections are shown in Figure 12, and are discussed further below.

Factorization test
In the flow paradigm, a two-particle azimuthal modulation characterized by a  , value arises from the product of nonzero azimuthal anisotropies,   , for each particle.These are related via  , (   T ,   T ) = 1.1 1.12 1.14 and  are selected from the identical particle  T range.Thus, a single-particle flow coefficient   (   T ) may be determined from two-particle  , values through for a given selection on reference particle .To test whether the  , values in data are compatible with this picture, a factorization test can be performed in which   values for particle  are compared for different particle  selections.The results of this test for the  2 values as a function of  rec ch are shown in Figure 13.The test demonstrates that while the  2,2 values for different   T selections may be different, the  2 values obtained for particle  as a function of  rec ch are independent of the selection for particle .This is a necessary condition of the paradigm that the observed  2,2 values arise from the product of event-wide, single-particle  2 azimuthal modulations.
Figure 14 shows the results of a factorization test as a function of particle   T .Although this test is more statistically limited than that for the  rec ch -dependent test above, the results for   T < 2 GeV are compatible with the factorization assumption.Since non-flow effects are not necessarily expected to factorize, this validates an important assumption of the template procedure.For the selection with   T > 3 GeV, the  2,2 values obtained in Figure 14 are significantly negative, and are accompanied by large  3,3 values.A negative  2,2 value violates the expected factorization for particles with independent overall modulation, and thus cannot be interpreted as arising from hydrodynamic flow.
As can be seen from the growing Fourier  2,2 values in Figure 11 and the large away-side peaks in Figure 12, the non-flow contribution in these high-   T selections is significantly larger than the magnitude of any possible modulation.Thus, the extraction of any genuine cos(Δ) modulation is highly sensitive to assumptions in the template method about the scaling of the non-flow contribution between multiplicity classes.If violated, this could manifest as a reduced  2,2 and increased  3,3 .This assumption is difficult to test given the limited data sample size for LM variations and factorization tests.Therefore, within the present uncertainties and techniques, the  2 cannot be reliably extracted at high  T .Finally, the statistical precision of the data does not allow similar tests of the  3 values.In this case, their presented values assume that factorization holds.In Section 7 the  2 and  3 values are presented for all  T for completeness, although it should be noted that factorization has been demonstrated only for  2 in the range 0.4 ≤  T ≤ 2.0 GeV.

Physics backgrounds
The measured   values are potentially biased by the presence of background events, primarily from peripheral hadronic Pb+Pb events, in the selected event sample.The purity of photonuclear events is estimated by fitting templates derived from simulated photonuclear and peripheral Pb+Pb events to distributions measured in data, and determining the fraction of selected events which are compatible with photonuclear interactions.
The purity was first estimated by fitting the Σ  Δ distribution in data to the combination of distributions from DPMJET-III +Pb events and peripheral H Pb+Pb events over the range 0.5 < Σ  Δ < 4.9, with the relative contributions from each as free parameters.An example of these two contributions and the distribution in data is shown in the right panel of Figure 4.This fit was performed separately in each of the  rec ch selections used in the measurement.The purity for the Σ  Δ > 2.5 selection is estimated to be greater than 98% in all  rec ch selections used in the analysis.Thus, for the nominal results, no correction for a possible impurity in the selected events is performed.

Systematic uncertainties
The primary sources of systematic uncertainty are grouped by category and discussed below.The systematic uncertainties are evaluated by varying the procedure applied to the data, and are not directly derived from the simulation samples described in Section 4.2.Thus, the uncertainties do not depend on an accurate modeling of the photonuclear processes in simulation.The uncertainties from each source are ultimately combined in quadrature to obtain the total systematic uncertainty.For some individual sources below where only one variation of the analysis procedure was performed, the resulting uncertainty was symmetrized by adding an uncertainty of equal magnitude in the other direction.However, in general the final systematic  1.
Pair-pseudorapidity difference.The two-particle correlation function is measured for pairs with 2.0 < |Δ| < 5.0 as a way to strongly reduce the non-flow contribution, which varies much more strongly with Δ than the genuine long-range azimuthal modulation.To test the sensitivity of the analysis to a potential residual contribution, the lower bound on Δ is varied to 1.8 and to 2.2, and the impact of the resulting variations on the data is included as an uncertainty.
Low-multiplicity reference.Photonuclear events with 15 ≤  rec ch ≤ 20 are used as the LM reference in the template fit to the other  rec ch selections.In the template method, the shape of the non-flow contribution is assumed to be the same in the LM and HM selections.To test the possible variation of the non-flow contribution with multiplicity, the analysis may be repeated using an alternative reference for LM events.Additionally, a different multiplicity range may select a different sample of photonuclear events (i.e.distribution of photon energies, or mix of resolved and direct processes [18]), to which the results may be sensitive.An alternative selection, 10 ≤  rec ch < 15, composed from the MB trigger data, is considered for the LM reference.The analysis is repeated and the difference in the results is included in the uncertainties.
Δ binning.The sensitivity to the particular choice of Δ binning was evaluated by repeating the fitting procedure after halving the number of bins.The resulting difference is included in the total uncertainty.Acceptance effects in correlation function.For the nominal analysis, the Δ-dependent correlation functions are formed by projecting the 2-D correlation functions in the region 2 < |Δ| < 5. To test the possible impact of acceptance effects that may result in a distortion of the Δ distribution, the correlation functions are corrected by the projection of the mixed-event 2-D correlation function onto Δ, as described in Section 5.The differences in the results are assigned as a systematic uncertainty.
Exclusion of tracks from Σ  Δ.The default analysis excludes reconstructed tracks that are within the sum-of-gaps range from the two-particle correlations, i.e. those with  > 4.9 − Σ  Δ.As a systematic variation, these tracks are included in constructing the correlation functions.The resulting differences are included as a systematic uncertainty.
Photonuclear event purity.The estimate of the purity of photonuclear events in data relies on the description of +Pb event properties in simulation.To test the sensitivity to this description, the fit described in Section 5.3 was repeated using P + events as the model for photonuclear events and, alternatively, by fitting to the  rec ch distribution, separately in different slices of Σ  Δ.In the latter case, the purity decreases to approximately 90% in events with  rec ch = 15-20 and to 80% in those with  rec ch = 10-15, which are the ranges used to define the nominal and alternative LM selections for the two-particle correlation analysis.
In the default analysis, no correction is made for a possible presence of backgrounds, primarily peripheral Pb+Pb events, in the selected sample of high-multiplicity +Pb events.As a systematic variation, the Δ distributions in data are corrected by subtracting the contribution from these background events, scaled to match the purity in each  rec ch range.Since the potential impact of the purity is much larger for events with low  rec ch than those with high  rec ch , the effect on the analysis is evaluated by accounting for the possible presence of backgrounds only in the LM reference.Additionally, since the purity was observed to be lower in the alternative LM reference events,  rec ch = 10-15, than in the nominal LM selection, the difference in the results using the alternative LM reference, with and without the purity correction, was taken as the uncertainty for this source.The model Δ distributions for the background events are taken from   collision data at √  = 5 TeV [57], with the same reconstruction, definition of  rec ch , and gap-based event selection, as this analysis.The resulting differences are included as a systematic uncertainty.
Other uncertainty sources and checks.Several other potential sources were investigated, and were found to have a negligible effect on the results and thus are not included in the total uncertainty.These include: comparing the results using only events where the photon is headed towards the A side of the detector with those where it is headed towards the C side; the number of Fourier terms used in parameterizing the LM reference distribution; the number of additional Fourier terms beyond  = 4 used in the template fit; the impact of a possible misalignment of the tracker on the charged-particle measurement; the impact of the correction for the finite tracking efficiency; and the impact of correcting for the trigger efficiency and different luminosities of the triggered event samples as a function of  rec ch on a per-event basis.While the Σ  Δ > 2.5 requirement is treated as an operational event definition which is not a source of uncertainty, the results using alternative definitions Σ  Δ > 2.0 and Σ  Δ > 3.0 were not significantly different from the nominal results.This demonstrates that the sample of photonuclear events in data is robust against small changes in their selection criteria and against detector effects in the measurement of Σ  Δ.
Finally, the template method carries an implicit assumption that the modulation in the HM selection is the same as that in the LM selection.If this is not the case, the extracted   values in the HM selection will be biased to be systematically closer to those in the LM selection.In recent measurements of azimuthal anisotropies in small systems [58], a 'one-step' correction was employed using two adjacent LM reference bins to account for the possibly changing   magnitudes with multiplicity.Since the measured   values were observed to have no  rec ch dependence within uncertainties, this alternative method was not utilized here.
As can be seen in Table 1, the systematic uncertainty is not dominated by any single source, and typically those from the pair-pseudorapidity difference, the choice of LM reference, and the event purity are co-dominant.The systematic uncertainties tend to be larger for the  T -dependent results than for the  rec ch -dependent results, and are typically larger for the  3 results than for the  2 results.For the  2 and  3 results as a function of  rec ch , the statistical uncertainties are typically larger than the systematic uncertainties.However, for the results as a function of  ch T , the systematic uncertainties are generally larger.This reflects, for example, the greater sensitivity to the choice of LM reference, which must be divided into different  ch T selections.: Flow coefficients  2 and  3 for charged particles with 0.5 <  T < 2.0 GeV in photonuclear events, reported as a function of charged-particle multiplicity  rec ch .The vertical error bars and colored boxes represent the statistical and total systematic uncertainties, respectively.The photonuclear data points are positioned at the average  rec ch value in each interval.The measurements in photonuclear events (solid symbols) are compared with those in   collisions at 13 TeV and +Pb collisions at 5.02 TeV [5] (open symbols), integrated over 0.5 <  T < 5.0 GeV.

Results and discussion
This section presents the measurement of the second-and third-order flow coefficients in photonuclear collisions, in both cases as a function of  rec ch and charged-particle  T .These results are compared with previous measurements of  2 and  3 for inclusive charged hadrons by ATLAS in inelastic, minimum-bias 13 TeV   collisions and in 5.02 TeV +Pb collisions [5].
Figure 15 presents the  2 and  3 values as a function of the event charged-particle multiplicity  rec ch .Significant, nonzero values for  2 and  3 are observed and they are compatible with no  rec ch dependence within uncertainties.The results are compared with previous measurements in   and +Pb collisions, presented here as a function of  rec ch for a slightly different  T selection, 0.5 <  T < 5 GeV.Given the shape of the  T distribution in high-multiplicity   collisions [59] and the  T -dependence of the  2 values, the  2 in   collisions for 0.4 <  T < 2 GeV particles (the selection used in this measurement) may be at most 10% lower than the values presented for 0.5 <  T < 5 GeV.In photonuclear events, the  2 values as a function of  rec ch are systematically below those in   and +Pb events, even after accounting for the different  T ranges in the datasets.The  3 values are compatible with those in   and +Pb events within significantly larger statistical and systematic uncertainties.
The  2 and  3 results as a function of charged-particle  T are presented in Figure 16 for events with 20 <  rec ch ≤ 60.The  2 results have central values similar to those plotted as a function of  rec ch in Figure 15 but have larger systematic and statistical uncertainties, due to the more limited size of the available data sample for the  T -dependent LM references for these measurements.For  2 , the results trend towards negative values for  T ≥ 2.0 GeV, although with very large systematic uncertainties.It is notable that there is a statistically significant confirmation of factorization, as detailed in Section 5.2, only for  2 with the selection 0.4 ≤  T ≤ 2.0 GeV.Hence the values of  2 and  3 for  T ≥ 2 GeV have an additional caveat, arising from the fact that the magnitude of the flow signal relative to that of the non-flow contribution is significantly smaller than at lower  T .In particular, the trend towards negative  2 values and rising  3 values suggests that the factorization assumption could be violated.
Figure 17 shows the same data as Figure 16, but zoomed in on the vertical axis to allow a better comparison with the analogous  T -dependent values in the   and +Pb measurements described above, with the selection  rec ch ≥ 60.In the region 0.4 <  T < 2 GeV the central values of the  2 are smaller than those in   and +Pb collisions, similar to that observed in the  rec ch -dependent results in Figure 15.However, due to the larger uncertainties in the  T -dependent case, the  2 values for photonuclear and   collisions are compatible within the uncertainties of the former in the range  T < 2 GeV.The  3 values are compatible between systems within large uncertainties.
There are currently no published theoretical predictions for flow coefficients in photonuclear collisions within a hydrodynamic or parton transport framework.In such frameworks, the elliptic and triangular flow coefficients scale with the initial geometry eccentricities,  2 and  3 respectively, and the charged-particle multiplicity  ch /.In the vector-meson dominance picture, photon-hadron interactions arise through fluctuations of the photon into hadronic states with the same quantum numbers as vector mesons, which have a nontrivial initial transverse geometry.This geometry is determined by the spectrum of these fluctuations, and while models of this spectrum exist [60], they have not yet been adapted to provide quantitative models.In the absence of a complete model, the magnitude of the eccentricities can be estimated by noting that fluctuations of the photon into light vector-meson states such as the  give the largest contribution to the cross section.The initial geometries for +Pb collisions can be computed with a Monte Carlo Glauber calculation [61] which treats the  meson as having two constituent quarks.The resulting mean values of the second-and third-order spatial eccentricities,  2 and  3 , are nearly identical to those in the +Pb case.Also, when comparing +Pb and photonuclear events with the same  rec ch , in fact the relevant  ch / is larger in the photonuclear events since the particles are distributed over a smaller pseudorapidity region.Thus, one might naïvely expect the flow coefficients to be similar in photonuclear events and +Pb collisions.However, in order to compare any such calculation with data, a full modeling of the photon fluctuations in the selected +Pb events needs to be carried out.In addition, correctly accounting for the boosted kinematics and limited acceptance using a fully three-dimensional simulation may be important.
An alternative interpretation of two-particle correlations in small collision systems involves interactions at the earliest time between gluon fields in the Color Glass Condensate (CGC) framework [62].Recently such calculations have described heavy-flavor hadron and quarkonia azimuthal anisotropies in +Pb collisions [63,64], although calculations in the CGC framework fail to describe other aspects of the data, such as the charged-hadron flow coefficients in +Pb at the LHC and small-systems collisions at RHIC [65,66].The authors have extended these calculations to consider a color dipole interacting with a Pb nucleus either at a future Electron Ion Collider or in photonuclear collisions at the LHC [31].The CGC calculation for photonuclear collisions is shown in Figure 17 and is in reasonable agreement with the  2 data within uncertainties.In these calculations, the Pb nucleus is described with a saturation scale  2  = 5 GeV 2 and typical parton transverse momentum Δ = 0.5 GeV, as used in calculations of  2 for heavy-flavor mesons and quarkonia [63,64].However, in the calculation for the photonuclear case, the parameter   = 25 GeV −2 ,  ch ≤ 60, reported as a function of particle  T .The vertical error bars and colored boxes represent the statistical and total systematic uncertainties, respectively.The photonuclear data points are positioned at the average  T value in each interval.The data are compared with the analogous measurements in   collisions at 13 TeV and +Pb collisions at 5.02 TeV for  rec ch ≥ 60 [5].The  2 data are also compared with a CGC-based theory calculation from Ref. [31].These photonuclear data are the same as in Figure 16 but with different y-axes ranges to allow comparison with additional data and theoretical predictions.which controls the transverse area of the interaction and thus the number of color domains from the Pb nucleus taking part in the interaction, is significantly larger than in the previous calculations [63,64].A larger transverse area generally contains more color domains and results in a smaller value of  2 .
Finally, it is interesting to compare the data with results from DPMJET-III, which has neither CGC-related correlations nor hydrodynamic or parton scattering effects.This comparison sets a baseline and provides a check of the possible impact of issues such as the LM template choice (which accounts for non-flow) selecting events with photon energies that differ from the HM sample.The photon energy is expected to be a function of  rec ch , which is true within DPMJET-III as shown in Figure 6.If the non-flow contribution changes with photon energy, it may be that the non-flow correlation shape would not be the same in the LM and HM selections, which is an assumption of the non-flow subtraction procedure.In order to test if this effect might mimic the  2 and  3 signals observed, the DPMJET-III output has been run through the full analysis procedure, including the LM and HM selection and non-flow template fit.The resulting  , values are shown in Figure 18 for both data and DPMJET-III.The left panel shows the values that would result from a Fourier fit to the HM correlation function, without accounting for non-flow, while the right panel shows the template-subtracted results.In DPMJET-III, the template-subtracted  , values are both negative, in contradistinction to the opposite sign for these values observed in data.Since they give negative values, which are impossible in a flow factorization interpretation, one cannot extract  2 and  3 values to compare with the data.Nevertheless, this demonstrates that, at least within DPMJET-III, there is no trivial effect in the non-flow subtraction issue which could result in a false signal.

Conclusion
This paper reports a measurement of long-range two-particle correlations in high-energy photonuclear collisions.Events are selected from 1.7 nb −1 of Pb+Pb collision data at 5.02 TeV recorded by the ATLAS detector at the LHC in 2018 using a combination of criteria for the energy deposited exclusively in one side of the zero-degree calorimeter, and requirements on the sum-of-gaps in the photon-going direction (Σ  Δ > 2.5).Azimuthal correlation functions are reported for pairs of charged particles separated by 2 < |Δ| < 5 units of pseudorapidity.A template-fitting method is used to subtract the non-flow contribution.The single-particle flow coefficients  2 and  3 are reported as a function of charged-particle multiplicity and single-particle  T .Significant, nonzero  2 and  3 values are observed in photonuclear events, indicating that particles produced in these events participate in azimuthally dependent, collective motion.The  2 values are smaller than those reported in   and +Pb collisions at similar particle multiplicities.The results presented here provide new information about long-range collective behavior in an exotic collision system with different initial conditions than those in ordinary hadronic collisions.As such, they provide a new testing ground for different pictures of the physical origin of these correlations, and a starting point for understanding the systems which will be studied at the future Electron Ion Collider.

Figure 2 :
Figure 2: Left: Correlation of  rec ch and Σ  Δ for events selected by the MB trigger before application of gap-based event selection.Right: Correlation of Σ  Δ and Σ  Δ for events selected by the MB trigger, with  rec ch > 10.

Figure 3 :
Figure3: Left: Distribution of event charged-particle multiplicity,  rec ch , for the photonuclear event selection from multiple triggers, without prescale correction.Right: Distribution of event charged-particle multiplicity,  rec ch , for the photonuclear event selection (black) with corrections for trigger efficiency and prescale factors, and an alternative selection intended to select hadronic Pb+Pb events (red, see text).

Figure 4 :
Figure 4: Left:  rec ch distribution in data, corrected for trigger and reconstruction efficiency and normalized per event (black points), compared with that in DPMJET-III +Pb (dot-dashed green histogram), DPMJET-III + (dotted red histogram), and P + (dashed blue histogram).The bottom panel shows the ratios of the MC distributions to the data distributions.Right: Σ  Δ distribution in data for  rec ch ≥ 10 (black points), normalized per event, and compared with that in DPMJET-III +Pb (dot-dashed green histogram), P + (dashed blue histogram), peripheral H Pb+Pb (solid magenta histogram), and DPMJET-III + (dotted red histogram).

Figure 5 :
Figure 5: Left: Charged-particle pseudorapidity distribution,  ch /, in selected  rec ch ranges.The distributions are normalized to the same integral and are shown in arbitrary units.Here, positive and negative  denote the photon-going and nucleus-going directions, respectively.Right:  ch / distribution in data for  rec ch > 10 (black points), normalized per event, and compared with that in DPMJET-III +Pb (dot-dashed green histogram), P + (dashed blue histogram), peripheral H Pb+Pb (solid magenta histogram), and DPMJET-III + (dotted red histogram) with the same reconstruction-level selection as the data.All distributions have been normalized to have the same value as DPMJET-III +Pb at  = 0.

Figure 6 :
Figure 6: Correlation between photon energy   and event multiplicity  rec ch in simulated DPMJET-III +Pb events.The right-hand axis shows the center-of-mass energy  for the photon-nucleon collision system.The markers and vertical errors bars show the mean and standard deviation, respectively, of   as a function of  rec ch .

Figure 7 :Figure 8 :
Figure7: Two-dimensional normalized particle pair distributions in photonuclear events, corrected for acceptance effects with the mixed-event distribution, and presented as a function of Δ and Δ.The peak at (Δ, Δ) = (0, 0) is truncated to better show the structure of the correlation function.Each panel represents a different charged-particle multiplicity  rec ch range for the selection 0.4 < (   T ,   T ) < 2.0 GeV.

Figure 9 : 2 = 3 =
Figure 9: An example of the template-fitting procedure for a selected  T range.The left plot displays the LM data with open markers and the simultaneous fit in the green dotted line.The lower panel displays the pull distribution.In the top panel of the right plot, the solid red line shows the total fit to the HM data in black markers.The dashed green line shows the scaled LM plus pedestal, while the dashed blue and dotted magenta lines indicate the two flow contributions to the fit,  ridge 2 =  [1 + 2 2,2 cos(2Δ)] and  ridge 3 =  [1 + 2 3,3 cos(3Δ)], shifted upwards by  LM (0) for visibility.The middle-right panel shows the pull distribution for the template fit in the top panel.The bottom-right panel shows the same set of data and fit components, where the scaled LM distribution has been subtracted to better isolate the modulation.

Figure 10 : 2 = 3 =vFigure 11 :
Figure 10: Selected template fit results for two   T intervals (top row) and two  rec ch intervals (bottom row).In the top panel of the right plot, the solid red line shows the total fit to the HM data in black markers.The dashed green line shows the scaled LM plus pedestal, while the dashed blue and dotted magenta lines indicate the two flow contributions to the fit,  ridge 2 =  [1 + 2 2,2 cos(2Δ)] and  ridge 3=  [1 + 2 3,3 cos(3Δ)], shifted upwards by  LM (0) for visibility.The middle panels show the pull distribution for the template fits in the top panel.The -values from the bootstrapping experiment are also shown in the middle panels.The bottom panels show the same set of data and fit components, where the scaled LM distribution has been subtracted from each to better isolate the modulation.

Figure 12 : 2 =Figure 13 :
Figure 12: Selected template fit results for the highest   T intervals.In the top panel of the right plot, the solid red line shows the total fit to the HM data in black markers.The dashed green line shows the scaled LM plus pedestal, while the dashed blue and dotted magenta lines indicate the two flow contributions to the fit,  ridge 2 =  [1 + 2 2,2 cos(2Δ)]

Figure 14 :
Figure 14: Summary of the factorization test for  2 values as a function of   T .Comparison of  2,2 (left) and derived  2 values (right) using the same particle  selection, but with different selections on particle .Only statistical uncertainties are shown.The points are displaced horizontally for visibility.

3 v 15
Figure15: Flow coefficients  2 and  3 for charged particles with 0.5 <  T < 2.0 GeV in photonuclear events, reported as a function of charged-particle multiplicity  rec ch .The vertical error bars and colored boxes represent the statistical and total systematic uncertainties, respectively.The photonuclear data points are positioned at the average  rec ch value in each interval.The measurements in photonuclear events (solid symbols) are compared with those in   collisions at 13 TeV and +Pb collisions at 5.02 TeV[5]  (open symbols), integrated over 0.5 <  T < 5.0 GeV.

Figure 17 :
Figure17: Charged-particle flow coefficients  2 (left) and  3 (right) in photonuclear events with 20 <  rec ch ≤ 60, reported as a function of particle  T .The vertical error bars and colored boxes represent the statistical and total systematic uncertainties, respectively.The photonuclear data points are positioned at the average  T value in each interval.The data are compared with the analogous measurements in   collisions at 13 TeV and +Pb collisions at 5.02 TeV for  rec ch ≥ 60[5].The  2 data are also compared with a CGC-based theory calculation from Ref.[31].These photonuclear data are the same as in Figure16but with different y-axes ranges to allow comparison with additional data and theoretical predictions.

vFigure 18 :
Figure 18: Comparison of results for raw Fourier coefficients  2,2 and  3,3 (left, without non-flow subtraction) and for non-flow subtracted coefficients  2,2 and  3,3 (right, with non-flow subtraction using the template method), shown in data (open points) and in DPMJET-III (filled points).The results in data and DPMJET-III are presented as a function of  rec ch and  truth ch , respectively.

Table 1 :
Summary of sources of systematic uncertainty, and the typical values of their absolute magnitude.The last row shows the typical total uncertainties from the quadrature sum of all uncertainties.