UvA-DARE (

Azimuthal anisotropies of muons from charm and bottom hadron decays are measured in Pb+Pb collisions at √ s NN = 5 . 02 TeV. The data were collected with the ATLAS detector at the Large Hadron Collider in 2015 and 2018 with integrated luminosities of 0 . 5 nb − 1 and 1 . 4 nb − 1 , respectively. The kinematic selection for heavy-ﬂavor muons requires transverse momentum 4 < p T < 30 GeV and pseudorapidity | η | < 2 . 0. The dominant sources of muons in this p T range are semi-leptonic decays of charm and bottom hadrons. These heavy-ﬂavor muons are separated from light-hadron decay muons and punch-through hadrons using the momentum imbalance between the measurements in the tracking detector and in the muon spectrometers. Azimuthal anisotropies, quantiﬁed by ﬂow coeﬃcients, are measured via the event- plane method for inclusive heavy-ﬂavor muons as a function of the muon p T and in intervals of Pb+Pb collision centrality. Heavy-ﬂavor muons are separated into contributions from charm and bottom hadron decays using the muon transverse impact parameter with respect to the event primary vertex. Non-zero elliptic ( v 2 ) and triangular ( v 3 ) ﬂow coeﬃcients are extracted for charm and bottom muons, with the charm muon coeﬃcients larger than those for bottom muons for all Pb+Pb collision centralities. The results indicate substantial modiﬁcation to the charm and bottom quark angular distributions through interactions in the quark-gluon plasma produced in these Pb+Pb collisions, with smaller modiﬁcations for the bottom quarks as expected theoretically due to their larger mass. © 2020 The Author(s). Published by Elsevier B.V. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/). Funded by SCOAP 3 .


Introduction
The paradigm for the time evolution of heavy-ion collisions at the Relativistic Heavy Ion Collider (RHIC) and the Large Hadron Collider (LHC) involves the formation and hydrodynamic expansion of a region of hot and dense quark-gluon plasma (QGP) with a small ratio of the shear viscosity to entropy density. In this paradigm, the QGP is considered to be a nearly perfect fluid [1,2]. Initial geometric inhomogeneities of the QGP are translated into momentum anisotropies of the final-state hadrons via large pressure gradients. Extensive measurements of light-hadron azimuthal anisotropies have been performed, in which the single-particle azimuthal distributions are expressed in terms of a Fourier expansion: dN dφ ∝ 1 + 2 ∞ n=1 v n cos(n(φ − n )), (1) E-mail address: atlas .publications @cern .ch.
where the event-plane angle, n , specifies the orientation of the initial density profile in the transverse plane [3], and Fourier coefficients, v n , quantify the magnitude of the modulation with respect to the event-plane angle. The second-and third-order v n coefficients are referred to as elliptic (v 2 ) and triangular (v 3 ) flow coefficients, respectively, with the term 'flow' invoking the hydrodynamic paradigm.
Heavy-flavour (charm and bottom) quarks have masses much larger than the temperature of the QGP (m c,b > T ), with maximum temperatures at early times ranging between 300 and 500 MeV [4]. Thus, thermal production of heavy quarks during the QGP phase is highly suppressed. Instead, heavy quarks are typically produced at the earliest times via high-momentum-transfer collisions between incoming partons. Once created, the heavy quarks persist throughout the dynamical time evolution of the QGP and thus act as sensitive probes of the hot and dense medium. Owing to their larger masses, radiative energy loss of heavy quarks in the QGP is suppressed relative to that of light quarks [5]. However, it was still postulated that charm quarks interact strongly enough to flow with the QGP [6]. Experimental data at RHIC and then at the LHC reveals that heavy-quark hadrons, as well as their decay leptons, have transverse momentum (p T ) distributions that are strongly modified by the QGP relative to observations in proton-proton (pp) collisions [7][8][9][10][11][12]. Charm hadrons [13,14] and heavy-quark hadron decay leptons [7,15] are also observed to have significant azimuthal anisotropies, suggesting that they participate in the overall collective flow of the medium. For recent reviews of heavy-flavour measurements in heavy-ion collisions, see Refs. [16][17][18].
For p T 4-6 GeV, it was proposed that heavy quarks can be described via a Langevin approach with drag and diffusion terms [19]. Modified p T distributions and azimuthal anisotropies of D mesons have been used to constrain heavy-quark transport coefficients [20,21]. Other models of heavy-flavor kinematics in the QGP, including a Boltzmann approach, have also been explored [22][23][24][25][26]. At higher momenta p T 5 − 10 GeV, heavy-quark energy loss is thought to dominate, with collisional and induced radiative processes both contributing [27]. At intermediate p T hadronization effects can be important as azimuthal anisotropies for the deconfined heavy-quark is transferred to the heavy-flavor hadron [28]. There are numerous theoretical predictions for the azimuthal anisotropies of bottom quarks, e.g. in Refs. [29][30][31]; however, only limited experimental data are currently available. Precision experimental data for p T distributions and azimuthal anisotropies is crucial as this over-constrains the calculations that depend on the heavy-quark to QGP coupling as well as the QGP space-time evolution.
The flow coefficients v 2 , v 3 , and v 4 of inclusive heavy-flavour muon production, which includes both muons from charm hadron decays ("charm muons") and muons from bottom hadron decays ("bottom muons"), have been measured by the ATLAS experiment [7] and ALICE experiment [32] in Pb+Pb collisions at √ s NN = 2.76 TeV. The measurement indicates significant elliptic flow for heavy-flavour muons with 4 < p T < 10 GeV. Recently, the heavyflavour muon v 2 has also been measured in high-multiplicity √ s = 13 TeV pp collisions [33]. Unlike the earlier Pb+Pb measurement, the pp measurement examined charm muons and bottom muons separately, finding a non-zero v 2 for charm muons while the v 2 for bottom muons is consistent with zero within uncertainties.
In the measurement presented in this paper, the procedure of the previous √ s NN = 2.76 TeV Pb+Pb analysis using the eventplane method is followed [7], and is extended to extract separate flow coefficients for charm and bottom muons. These measurements extends the previously published ones to the higher √ s NN = 5.02 TeV Pb+Pb beam energy, using a larger event sample provided by the 2015 and 2018 combined data sets. The larger data sample enables measurements over a larger momentum range 4 < p T < 30 GeV for inclusive heavy-flavour muons. It also allows the separation of the inclusive heavy-flavour muons into charm and bottom contributions. Results for charm and bottom muon elliptic v 2 and triangular v 3 flow coefficients are presented as a function of muon p T for various ranges of overlap between the colliding nuclei, referred to as "centrality".

ATLAS detector
The ATLAS detector [34][35][36] at the LHC is a multipurpose particle detector with a forward-backward symmetric cylindrical geometry and a near 4π coverage in solid angle. 1 It consists of an inner tracking detector surrounded by a thin superconducting solenoid 1 ATLAS uses a right-handed coordinate system with its origin at the nominal interaction point (IP) in the centre of the detector and the z-axis along the beam pipe. The x-axis points from the IP to the centre of the LHC ring, and the y-axis points upwards. Cylindrical coordinates (r, φ) are used in the transverse plane, φ being the azimuthal angle around the z-axis. The pseudorapidity is defined in terms of the polar angle θ as η = − ln tan(θ/2). Angular distance is measured in units of providing a 2 T axial magnetic field, electromagnetic and hadron calorimeters, and a muon spectrometer. and consist of layers of alternating quartz rods and tungsten plates. A two-level trigger system [37] is used to select events. The firstlevel trigger (L1) is implemented in hardware and uses a subset of the detector information to reduce the accepted rate to below 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 Pb+Pb collisions.

Event selection
Data used in this analysis were recorded with the ATLAS detector in 2015 and 2018 from Pb+Pb collisions at √ s NN = 5.02 TeV with integrated luminosities of 0.5 nb −1 and 1.4 nb −1 , respectively. Events were selected online using a set of muon triggers that require a muon at the HLT stage with p T larger than 3, 4, 6, or 8 GeV [37]. The muon trigger selecting p T > 8 GeV sampled the full luminosity in both 2015 and 2018, while triggers with lower p T thresholds were prescaled to reduce the overall data rate. Thus the higher-threshold triggers are utilised at a given muon p T to sample a larger fraction of the full luminosity. The resulting sampled luminosities are 0.3 nb −1 , 0.6 nb −1 , and 1.9 nb −1 for muons with 4 < p T < 6 GeV, 6 < p T < 8 GeV and p T > 8 GeV, respec- Muons with 4 < p T < 30 GeV and |η| < 2.0 reconstructed in both the ID and the MS are selected and required to pass 'medium' selection requirements, detailed in Ref. [40], without any requirement on the number of TRT hits, due to the high occupancy in heavy-ion running. Candidate muons are required to be matched with a muon reconstructed by the event trigger. Each muon is assigned a weight given by the inverse of the reconstruction and trigger efficiency for the specific muon kinematics. The muon reconstruction and trigger efficiencies are determined using the J /ψ → μ + μ − tag-and-probe method as detailed in Ref. [40]. The muon reconstruction efficiency is factorised as the product of ID and MS reconstruction efficiencies. The ID reconstruction efficiency is obtained from Pb+Pb data directly, while the MS reconstruction efficiency is obtained both from Pb+Pb data and by overlaying Pb+Pb minimum-bias events on simulated J /ψ produced by Pythia 8 [41] with the CTEQ6L1 [42] parton distribution functions, using the same tag-and-probe method. The events in overlay simulations are then given weights such that the E F Cal T distribution matches the muon-triggered Pb+Pb data distribution. The MS reconstruction efficiency obtained from simulation is used as the central value for MS reconstruction, with additional data-to-MC scale factors applied to account for residual differences between data and overlay simulation. The same J /ψ → μ + μ − tag-and-probe method is used to determine the muon trigger efficiency. The central value of the single-muon trigger efficiency correction is obtained from simulations without overlaying Pb+Pb minimum-bias events. Additional correction factors of order 1-10% are obtained from data and applied to the selection to correct for trigger detector performance differences between data and simulation, as well as the centrality dependence of the muon trigger efficiency in Pb+Pb data.

Analysis procedure
The analysis follows the event-plane method for measuring flow coefficients as used in previous ATLAS measurements [7,39] and is briefly summarised here. Each Pb+Pb event has a geometric orientation of the impact parameter vector, and the event can also have a tilt relative to that due to fluctuations in the geometry of the resulting QGP. In any particular Pb+Pb collision, one can estimate the orientation, represented by the FCal-determined n thorder event-plane angle n . The azimuthal distribution of transverse energy deposited in the forward and backward rapidity FCal is used to determine the event plane. A comparison of event planes, as measured separately in the forward-rapidity FCal and the backward-rapidity FCal, enables a determination of the eventplane resolution Res{n n } as detailed in Ref. [39]. In each Pb+Pb centrality interval and in each muon p T selection, the muons are divided into a finite number of intervals in φ − n , where φ is the azimuthal angle of the muon. As different harmonic orders are orthogonal to each other, the Fourier decomposition of the angular distribution (introduced in Eq. (1)) at a given order n can be expressed as Similarly to previous ATLAS publications [7,33,43], different sources of muons are separated using two variables. The first is the momentum imbalance, ρ = (p ID − p MS )/p ID , where p ID is the muon momentum measured in the ID, and p MS is that measured in the MS corrected for the energy loss inside the calorimeter.
Real muons have a ρ distribution peaked around zero while the π/K background has a broader ρ distribution that is shifted toward higher values. The different shapes of the ρ distribution for the π/K background and other muons enable the isolation of the π/K background muons. The analysis is repeated using the transverse momentum imbalance, as opposed to the total momentum imbalance ρ, and no difference is observed. The second variable is the transverse impact parameter, d 0 , relative to the event's primary vertex [44]. Charm and bottom muons have different d 0 distributions due to the different decay lengths of charm and bottom hadrons.
A two-step fit in ρ and d 0 is performed in data, using ρ and d 0 line-shape templates for different sources of muons obtained from simulations. First, the yields of inclusive heavy-flavour and π/K background muons are extracted from a fit to the ρ distribution. The relative yields of light/onia background muons and inclusive heavy-flavour muons are fixed to the fractions obtained from Pythia 8 simulations, which are approximately 4% on average. Then, with the extracted π/K background yields fixed, a fit to the d 0 distribution is performed to determine the relative fraction of charm and bottom muons within the yield of inclusive heavyflavour muons.
The muon ID momentum resolution in Pythia 8 simulations overlaid with minimum-bias Pb+Pb events is found to be worse than the resolution in Pb+Pb data. Thus, the ρ templates are obtained from simulation without Pb+Pb event overlay. The singlemuon ID and MS momentum responses in the Pythia 8 simulation are shifted and smeared in order to match those in Pb+Pb data.  on muon p T or event centrality but a moderate dependence on parent charm and bottom hadron p T due to the strong correlation between decay length and particle velocity. Additional reweighting is applied to the charm and bottom muon signal samples to match the input charm and bottom hadron p T spectra to those measured in Pb+Pb collisions by ALICE [11] and CMS [8,47]. The fits are performed independently in differential intervals of muon p T , centrality, n|φ − n | and two intervals of muon η. The two muon η intervals (|η| < 1 and 1 < |η| < 2) are fitted independently to minimise residual data/MC difference in the barrel and endcap muon detectors separately, and then combined to obtain charm, bottom, and inclusive heavy-flavour muon yields in the given p T , centrality, and n|φ − n | intervals as reported in the results. Fluctuations in the simulation templates are incorporated in the fitting procedure. Examples of selected fits in ρ and d 0 based on simulation templates are shown in Fig. 1 for two different muon

Systematic uncertainties
Systematic uncertainties are presented for different categories covering each step of the analysis procedure: 1) muon efficiency; 2) ρ fit; 3) d 0 fit; 4) light/onia background; 5) other background; 6) ρ-d 0 correlation; 7) event-plane resolution; and 8) jet bias. Table 1 summarises the contributions from different sources of sys-tematic uncertainty to the final flow-coefficient results. Systematic uncertainties from all sources are summed in quadrature to determine the total uncertainty.
The systematic uncertainties from the MS reconstruction efficiency and muon trigger efficiency corrections are dominated by the uncertainty in determining the data-to-MC scale factor. The scale factor uncertainties are evaluated following the procedure from previous ATLAS measurements [40] including variations in the tag-and-probe efficiency extraction method, object-matching resolution, and purity of the sample. The systematic uncertainty in the muon trigger efficiency also includes the determination of the centrality-dependent correction factors. The uncertainty on the flow coefficients resulting from uncertainties in the muon ID reconstruction efficiency is evaluated by comparing the results with and without the ID efficiency correction, as the ID efficiency is approximately 99% for all centralities. The muon efficiency systematic uncertainties are correlated between the resulting charm and bottom muon results.
The systematic uncertainty of the ρ fit is due to the uncertainty in the shifts and smearing parameters for single-muon momentum response in simulation, which is evaluated by comparing the nominal results with those 1) without any shifts or smearing, 2) only applying shifts and smearing to the signal muons, and 3) incorporating additional smearing of the background ρ distributions in simulation to match data distributions in the background-enriched region (ρ > 0.2). The changes resulting from these variations are combined in quadrature. The combined uncertainty from the ρ fit is 1%-10% for the charm and bottom muon results, depending on the muon p T and η, but without dependence on centrality. The relative systematic variations are found to be the largest at low p T and large η. The ρ fit systematic uncertainties are correlated between the resulting charm and bottom muon results.

Table 1
Summary of the typical sizes of the absolute systematic uncertainties of all categories for the flow coefficient results. The d 0 related systematic uncertainties are not relevant for the heavy-flavour inclusive flow measurements as the d 0 fit is not utilised for these results. Systematic uncertainties from the event-plane resolution and jet bias are negligible and not included in the final uncertainties, and therefore are not shown in the table. The "<" symbol indicates the provided values are the maximum systematic variation for a given category.
The uncertainty in the d 0 template shift and smearing parameters is tested by comparing results when determining the parameters using 2018 data (as in the nominal results) with results when using the 2015 Pb+Pb data to determine the parameters, which covers the slightly different detector alignment between the two data sets. Sensitivity to the charm and bottom hadron p T spectra reweighting in simulation is covered by a variation in which the p T spectra are reweighted to agree with those from Pythia 8 simulations without any modification due to QGP. The variations in the results due to d 0 template shift and smearing and p T spectra reweighting are considered to be uncorrelated and are combined in quadrature, and the combined systematic uncertainty is 1%-20% for the charm and bottom muon results, depending on muon p T and centrality. The relative systematic variations are found to be the largest at high p T and in more peripheral events. The d 0 fit systematic uncertainties are anti-correlated between the resulting charm and bottom muon results.
The magnitude of the light/onia contribution is held at a fixed fraction relative to the inclusive heavy-flavour muon contribution, based on the Pythia 8 MC sample. The analysis is repeated to study this choice, first ignoring the light/onia contribution and then fixing it to twice the fraction predicted by Pythia 8. As is shown in Ref.
[48], Pythia overestimates prompt quarkonium production at the LHC, and thus these variations of the light/onia contribution are large enough to not underestimate the uncertainty. Each nominal result is assigned a systematic uncertainty equal to the larger of the changes from the two variations. For the nominal results, light/onia muons are assumed to have the same v 2 and v 3 as the inclusive heavy-flavour muons. The analysis is repeated with variations assuming light/onia muons have zero flow coefficients or double the inclusive heavy-flavour muon flow coefficients. The larger of the resulting changes is assigned as the systematic uncertainty. The light/onia systematic uncertainties are anti-correlated between the resulting charm and bottom muon results. The systematic uncertainty associated with the event-plane angle resolution is evaluated by measuring the event-plane resolution in two subregions of the FCal (3.2 < |η| < 4.0 and 4.0 < |η| < 4.8), following a previous ATLAS flow analysis [39]. The systematic uncertainties are evaluated independently for 2 and 3 . The maximum difference between these two variations and the nominal results is considered as a systematic uncertainty. The uncertainty associated with the event-plane angle resolution is found to be negligible compared to other systematic uncertainties, and thus is not included in the results.
The charm and bottom muons are often produced with a recoil jet. The orientation of n could be biased to align with the signal muon if the recoil jet reaches the FCal [7]. The magnitude of this bias in muon v 2 and v 3 is studied with Pythia generator-level muon flow in samples overlaid with Pb+Pb data. The bias is found to be 0.3%-0.4% for v 2 and v 3 , and it is larger in peripheral events than in more-central events. This small bias is negligible and is not included as a systematic uncertainty.    Qualitatively, the charm and bottom v 2 ordering matches theoretical expectations where the heavier bottom quarks have a smaller modification to their initial momentum trajectories due to their larger masses. Light quarks and heavy quarks can lose energy in traversing the QGP via induced gluon radiation [51]; however, heavy quarks with momentum less than or approximately equal to the quark mass (p m) radiate less than light quarks due to a suppression of radiation at small angles relative to the quark direction, referred to as the 'dead-cone' effect [5]. Thus, at high p T >> m, all quark flavors should lose comparable energy in the   expansion, and thus it will be instructive in the future to compare the calculations with a common QGP model to test whether the differences arise from the QGP modelling or the energy-loss implementation.

Conclusion
In summary, a measurement of elliptic and triangular flow coefficients for heavy-flavour decay muons in Pb+Pb collisions at

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgements
We thank CERN for the very successful operation of the LHC, as well as the support staff from our institutions without whom ATLAS could not be operated efficiently. We