Proton, deuteron and triton flow measurements in Au+Au collisions at $\sqrt{s_{NN}} = 2.4$ GeV

High precision measurements of flow coefficients $v_{n}$ ($n = 1 - 4$) for protons, deuterons and tritons relative to the first-order spectator plane have been performed in Au+Au collisions at $\sqrt{s_{NN}} = 2.4$ GeV with the High-Acceptance Di-Electron Spectrometer (HADES) at the SIS18/GSI. Flow coefficients are studied as a function of transverse momentum $p_{t}$ and rapidity $y_{cm}$ over a large region of phase space and for several classes of collision centrality. A clear mass hierarchy is found for the slope of $v_{1}$, $d v_{1}/d y^{\prime}|_{y^{\prime} = 0}$ where $y^{\prime}$ is the scaled rapidity, and for $v_{2}$ at mid-rapidity. Scaling with the number of nucleons is observed for the $p_{t}$ dependence of $v_{2}$ and $v_{4}$ at mid-rapidity, which is indicative for nuclear coalescence as the main process responsible for light nuclei formation. $v_{2}$ is found to scale with the initial eccentricity $\langle \epsilon_{2} \rangle$, while $v_{4}$ scales with $\langle \epsilon_{2} \rangle^{2}$ and $\langle \epsilon_{4} \rangle$. The multi-differential high-precision data on $v_{1}$, $v_{2}$, $v_{3}$, and $v_{4}$ provides important constraints on the equation-of-state of compressed baryonic matter.


Introduction
Heavy-ion collisions are a tool to investigate the properties of strongly-interacting matter under extreme conditions, such as high temperatures typical for the early phase of the universe and high baryon number densities occurring in compact stellar objects [1].Especially the latter conditions cannot easily be addressed by ab-initio calculations based on Quantum Chromo-Dynamics (QCD), the theory of strong interaction, but are to be addressed with effective theories.Therefore, measurements are indispensable to determine the properties of dense matter, which can be -in local equilibrium -encoded in the equation-of-state (EOS) [2][3][4].
For non-polarised and equal-size projectile and target nuclei, the azimuthal uniformity of the momentum tensor in the final state is broken by the spatial asymmetry of the initial state at finite impact parameter which is transferred into momentum space via pressure gradients generated by multiple interactions of the matter constituents.The resulting structure of the azimuthal distributions of produced particles is conveniently parametrised by a Fourier decomposition [5] with the coefficients v n (p t , y): ( Here, ϕ is the azimuthal angle of the particle and Ψ RP stands for the azimuthal angle of the reaction plane, defined by the beam direction ⃗ z and the direction of the impact parameter ⃗ b of the colliding nuclei.p z = p cos θ is the momentum component along the beam direction ⃗ z with the laboratory momentum p and the polar angle θ, p t = p sin θ the one perpendicular to it, and y = tanh −1 (p z /E) denotes the rapidity of a given particle with energy E in the laboratory frame.
The rapidity in the centre-of-mass system is denoted by y cm = y− y proj /2, with the projectile rapidity y proj .The coefficients v 1 and v 2 quantify the so-called directed and elliptic flow, respectively.More generally, the shape of the anisotropy is quantified by odd coefficients v 1 , v 3 , v 5 . . .v 2n+1 and even coefficients v 2 , v 4 . . .v 2n .Due to momentum conservation, and assuming that initialstate fluctuations with large eccentricities are absent, the values of even coefficients are expected to be symmetric around mid-rapidity, while the ones of the odd harmonics should change their sign when going from forward to backward centre-of-mass rapidities.Flow data in the few GeV energy range at BE-VALAC and SIS18 have been reported for pions, charged kaons, hyperons, neutrons, as well as protons and many light nuclei.For reviews see [6][7][8][9] and references therein.Most of these data are for a limited phase space or integrated over transverse momenta.High-statistics, multidifferential data on v 1 and v 2 for identified particles measured over a large region of phase-space is a valuable extension of the existing world data.In addition, the study of higher-order flow coefficients can provide information about the various contributions to the bulk properties of dense nuclear matter.At RHIC and LHC energies, these were employed to determine the ratios of the shear and bulk viscosity to the entropy density Fig. 2 The phase-space coverage of identified protons (left panel), deuterons (middle panel) and tritons (right panel) as accepted in the HADES experiment for Au+Au collisions at √ s NN = 2.4 GeV as a function of the centre-of-mass rapidity ycm and transverse momentum p t (the curves correspond to the given polar angles θ in the laboratory system).
η/s, and ζ/s, respectively, of high-temperature matter [10,11].Attempts to extract η/s for dense baryonic matter have been made by comparing various model approaches to the available data [12][13][14][15][16][17], however with large uncertainties.Data on higher-order flow coefficients are essential in order to disentangle the effects of shear and bulk viscosity [17] and can also provide important information on the EOS.A comparison of the proton v 3 , measured by HADES, with UrQMD transport model calculations indicates an enhanced sensitivity to the EOS of the hadronic medium [18,19].
Other transport model calculations suggest that a nonvanishing value of v 4 , measured at center-of-mass energies of a few GeV, can constrain the nuclear mean field at high net-baryon densities [20].The E877 collaboration has reported a non-vanishing v 4 at 10.1 AGeV [21], but measurements of higher coefficients are generally scarce at low (i.e.AGS and SIS18) energies.A multi-differential measurement of several Fourier coefficients allows for a three-dimensional characterization of heavy-ion collisions in different representations [22][23][24].
The scaling properties of the flow coefficients v n of different order n with the number of nucleons A of the respective nucleus can provide information on the production mechanisms of light nuclei, e.g.via nucleon coalescence [25,26].The relation of the v n to the shape of the initial eccentricity of the collision system can shed light on the reaction dynamics and the transport properties of the produced medium.
In this article we present results on flow of protons, deuterons and tritons in Au+Au collisions at a centreof-mass energy of √ s NN = 2.4 GeV, corresponding to a kinetic beam energy of 1.23 AGeV.We extend our pre-vious study [22] to multi-differential data on the flow coefficients v 1 − v 4 as a function of transverse momentum and rapidity over a large region of phase-space and for several classes of reaction centrality.The paper is organised as follows.Section 2 introduces the experimental set-up and the particle reconstruction methods, while section 3 discusses the procedures used to determine the flow coefficients.Section 4 presents the results on directed, elliptic and higher harmonics for semi-central collisions.In section 5 our results of v 1 and v 2 are compared with existing world data, while in section 6 the scaling properties of the data are discussed.Section 7 presents comparisons to several model calculations.We summarize in section 8.

Experimental set-up
HADES is a charged-particle detector consisting of a six-coil toroidal superconducting magnet centered around the beam axis with six identical detection sections located between the coils and covering polar angles between 18 • and 85 • (see Fig. 1).Each sector is equipped with a Ring-Imaging Cherenkov (RICH) detector followed by four layers of Multi-Wire Drift Chambers (MDCs), two in front of and two behind the magnetic field, as well as a scintillator Time-Of-Flight detector (TOF) (polar angle coverage: 44  -45 • ).Hadron identification is based on the time-of-flight measured with TOF and RPC, and on the energy-loss information from TOF as well as from the MDC tracking chambers.Electron candidates are in addition selected via their signals in the RICH detector.Combin-ing this information with the track momentum, as determined from the deflection of the tracks in the magnetic field, allows for an identification of charged particles (e.g.pions, kaons or protons).
The spectrometer set-up is complemented by the Forward Wall (FW) detector.It consists of 288 scintillator elements of 2.54 cm thickness and three different front area sizes (innermost region: 4×4 cm 2 , intermediate region: 8×8 cm 2 and outermost region: 16×16 cm 2 ) which are read out with photomultipliers.The FW is placed downstream at a 6.8 m distance from the target and covers the polar angles 0.34 • < θ < 7.4 • .It thus allows for a measurement of the emission angles and the charge states of projectile spectators and is used to determine the event-plane angle.A detailed description of the HADES experiment can be found in [27].

Data set and event selection
Several triggers are implemented to start the data acquisition.The minimum-bias trigger is defined by a signal in a 60 µm thick mono-crystalline CVD 1 diamond detector (START) in the beam line [28].In addition, online Physics Triggers (PT) are used, which are based on hardware thresholds on the TOF signals, proportional to the event multiplicity, corresponding to at least 5 (PT2) or 20 (PT3) hits in the TOF.
By comparing the measured hit multiplicity distribution with Glauber and transport model simulations, it has been estimated that the PT3 trigger is selecting (43 ± 4) % (PT2 trigger: (72 ± 4) %) of the total inelastic cross section [29].The selection of centrality classes is based on the summed hit multiplicity of the TOF and RPC detectors.Four classes are defined, which together cover the 40 % most central collisions in steps of 10 % of the total Au+Au cross section of 6.83 ± 0.43 b.Events are selected offline by requiring that their reconstructed global event vertex lies inside the target region, i.e. between z = −65 mm and 0 mm along the beam axis.Additionally, only events with at least four hits in the FW with a charge Z ≥ 1 are used for the reconstruction of the event-plane.It was verified that this criterion does not introduce any significant bias to the centrality selection.The mean number of accepted proton, deuteron and triton candidates are summarized in Tab. 1.

Track selection
The flow coefficients are determined for the charged particles detected by the detectors MDC, TOF and  RPC.Their trajectories are reconstructed using the MDC information.The resulting tracks are selected according to the quality parameter provided by the employed Runge-Kutta track fitting algorithm χ 2 RK and a maximal Distance of Closest Approach (DCA) of the extrapolated track to the reconstructed primary vertex position.In order to assure a good matching of the tracks to the hits measured in the particle identification detectors TOF and RPC, an additional selection criterion is applied.It involves an upper limit on the quality parameter Q MM = dx/σ x , which is defined as the deviation of the intersection point of a given reconstructed track from the position of the associated hit in the RPC and TOF detectors, dx, normalized to the corresponding measurement uncertainty, σ x .The nominal selection values on χ 2 RK , DCA and Q MM are summarized in Tab. 2.

Particle identification
The Particle IDentification (PID) is based on a combined measurement of time-of-flight and energy loss.The time-of-flight, as determined by the TOF and RPC detectors, allows for a separation of particles in different momentum dependent regions of velocity β.To select protons, deuterons and tritons windows with widths of n σ β (p) with n = 2.5 are placed around the corresponding expected β values (see also Tab. 2).The respective resolutions σ β (p) depend on the particle momenta p and are parametrised accordingly.
In addition, the energy loss (dE/dx) measurements in the MDCs are employed for PID.This is particularly important to suppress the 4 He contamination in the deuteron sample, as the two nuclei cannot be separated by time-of-flight alone due to the same Z/A ratio.The variable Z MDC is constructed from the energy loss measured in all four MDC layers, dE/dx exp , and the theoretically expected value, dE/dx th , The selection windows applied to this variable are momentum independent, but different for protons, deuterons and tritons, see Tab. 2. The purity of the particle identification procedure is determined by analysing simulated data (see Sect. 3.2) and by fitting the mass distributions, calculated from the measured values of β and momentum for different rapidity and transverse momentum intervals, with a function that describes the signal as well as the background component.Phase-space intervals, in which the purity of the particle identification is found to be lower than 80 %, are excluded from further analysis.Also, intervals at the edges of the detector acceptance, i.e. on the borders of the polar angle range 16 • < θ < 85 • and the gaps in azimuth between the detector sectors, are excluded.This translates into a rapidity dependent lower transverse momentum cut-off.At very high momenta, phase-space regions are rejected if the accuracy of the momentum measurement is not sufficient.The phase-space coverage for the identified particles is shown in Fig. 2.

Determination of flow coefficients
The flow coefficients v n of order n are defined in their relation to the reaction plane angle Ψ RP as [30,31] Here, ⟨. . .⟩ denotes the average over all selected particles and events in a given sample.As the reaction plane is not accessible experimentally, it is replaced by the event-plane angle (Ψ EP ) constructed from measured event anisotropies as described in the following.

Reaction plane reconstruction
For the determination of the event-plane angle Ψ EP , hits of projectile spectators recorded by the FW are used as illustrated in Fig. 3. Projectile spectators are thus measured in the polar angle interval 0.34

HADES
Fig. 4 The distribution of the first-order event-plane angles Only those hits are used for which the energy deposit in the scintillator cells and their flight time corresponds to the values expected for spectators with charges of Z ≥ 1.
From the azimuthal angles ϕ FW of the FW cells hit by spectators, a vector ⃗ Q n = (Q n,x , Q n,y ) is calculated event-by-event: Here, N FW is the number of detected FW cell hits.The weights w i are here chosen to be w i = |Z i |, where Z i is the charge of a given hit as determined via the signal amplitude seen by the FW cell.Because of nonuniformities in the FW acceptance, caused by few dead cells and by deviations of the beam position relative to the nominal centre of the experimental set-up, the distribution of FW hits averaged over many events is not centred around the origin.To correct for this, the individual positions of the FW-hit X FW,i and Y FW,i are recentred by the corresponding first moments (⟨X FW ⟩, ⟨Y FW ⟩) and scaled by the second moments (σ X FW , σ Y FW ), which are calculated for each day of data-taking and centrality class separately.To remove the residual non-uniformities in the event-plane angular distribution an additional flattening procedure was applied [32].
The corresponding event-plane angle of order n is then defined as: Figure 4 shows distributions of the first-order eventplane angles before and after applying the above described correction procedure.As a result of the corrections, Ψ EP,1 is distributed uniformly in all centrality classes.The comparison to a flat distribution results in χ 2 /N DF values in the range 0.83 − 1.09 for the centralities 0 − 40 %.Generally, Ψ EP,n can be determined for each order n.As the reaction plane orientation is mainly connected to the deflection of the projectile spectators, n = 1 provides the highest resolution and therefore the firstorder event-plane angle Ψ EP,1 is used in the following for the extraction of the observable flow coefficients of all orders: As the computed event-plane angles will fluctuate around the true reaction plane angles, the observed flow coefficients will come out smaller than the true ones.This can be corrected using the resolution correction of the event-plane (Ψ EP,1 ):

HADES
Fig. 5 The resolution ℜn of the first-order spectator event plane for flow coefficients of different orders n as a function of the event centrality [22].The circles correspond to centrality intervals of 5 % width and the squares to 10 % width (curves are meant to guide the eye).
For the first-order event-plane, assuming that non-flow contributions can be neglected, the resolution can be expressed as [30,31,33,34] where I ν are the modified Bessel functions of the order ν and χ is the resolution parameter.For the determination of ℜ n , the two-sub-event method is employed.For this purpose the FW hits in a given event are randomly divided into two sub-events A and B of equal multiplicity.From the correlation of the two the resolution for the sub-events is then calculated as By replacing χ in Eq. ( 8) with χ sub and inverting the equation, the value for χ sub can be calculated.The value of the resolution parameter for the full FW is then χ = √ 2 χ sub , which yields the full resolution ℜ n after inserting it into Eq.(8).
The resulting values for the resolution of different order n are exhibited in Fig. 5.In the case n = 1, it is found to be around 80 % and higher in the centrality range 10 − 40 %, while it drops towards a value of ∼ 50 % for very central collisions.
Alternatively, two methods proposed in [30] are used.The resolution parameter χ is obtained by fitting Eq. ( 12) of [30] to the measured distribution of the differences between the two sub-event-plane angles ∆Ψ = |Ψ EP,A − Ψ RP,B | or by using the approximate relation where the effect of these differences on the systematic uncertainties is found to be negligible.

Correction for reconstruction inefficiencies
In the high-multiplicity environment of Au+Au collisions the reconstruction of tracks is affected by ambiguities in the assignment of firing MDC drift cells to a given track.This results in reconstruction inefficiencies which depend on the local track multiplicities N tracks .Anisotropies in the event shape, as caused by flow effects, will in turn generate local modulations of the track densities and thus of the reconstruction inefficiencies, which consequently distort the determination of the flow coefficients.Therefore, any efficiency correction must also account for the track orientation relative to the event-plane.
With the help of simulated data, generated using Geant3.21[35] in combination with a detailed description of the detector geometry and response, the multiplicity dependence of the reconstruction efficiency ϵ was studied.It was found that it can be described by the following function: From simulations, a maximal efficiency of ϵ max = 0.98 is determined.In the phenomenological data-driven approach used here the parameter c ϵ is adjusted such that v 1 = 0 for y cm = 0, as required by the symmetry of the reaction system.In a next step, the average local track multiplicity ⟨N loc.tracks ⟩ is calculated from data in intervals of the track polar angle θ, of the difference between its azimuth angle ϕ and the one of the eventplane ϕ EP = ϕ−Ψ EP,1 and of the event centrality.Using these three-dimensional matrices as input to Eq. ( 11), relative efficiency tables ϵ(θ, ϕ EP , cent.) are determined.Examples are shown in Fig. 6.These are then used to weight all tracks used to calculate the flow coefficients according to

Systematic uncertainties
The systematic uncertainties of the measured flow harmonics v n can be separated into global ones, i.e. those which affect all data points in the same way, and those which depend on phase-space position.The latter include effects of the reconstruction and selection of the tracks, particle identification, correction procedures for reconstruction inefficiencies and time-dependent changes of the acceptance.They are determined as a function of y cm and p t , separately for each particle species, the order of the flow harmonics v n , and the centrality class.
The systematic uncertainties due to the track reconstruction are estimated by varying the track selection cuts.Table 2 lists, in addition to the nominal selection criteria, also two values for each cut used for the determination of systematic effects.Impurities in the selected particle samples, i.e. a background of misidentified particles, will also modify the corresponding flow result.Their contribution to the systematic uncertainty is evaluated by varying the PID selection criteria (see also Tab. 2).
The parameter c ϵ used in the correction for multiplicity dependent inefficiencies in Eq. ( 11) is modified relative to its nominal value to evaluate the influence of the correction procedure on the resulting v n .This variation covers all values of c ϵ which are still compatible within errors with v 1 = 0 at mid-rapidity.
In larger periods of the data-taking time, one sector of the MDC was not fully operational.As this introduces an azimuthal asymmetry into the acceptance and therefore increases the sensitivity to an imperfect re-centring of the event-plane, it can be the cause of an additional systematic uncertainty.This is estimated by comparing the results obtained for a fully symmetric detector (i.e. six operational sectors) with those for only five sectors.In addition, configurations were analyzed with only four or even three active sectors, corresponding to the upper or lower part of the detector.
The total systematic uncertainly is derived by independently analyzing all different variations and then evaluating the overall distributions of the resulting flow coefficients.It is found that for the even coefficients all the effects described above contribute roughly at the same level to the point-by-point systematic uncertainties, whereas azimuthal anisotropies, like efficiency losses in whole sectors, dominate the systematic uncertainties of the odd flow coefficients.A summary of the different systematic uncertainties is given in Tab. 3.
In order to verify that the higher flow harmonics are not artificially generated by acceptance holes, a toy MC study was performed.This simulation mimics corresponding effects by passing tracks through an accep-   |ycm| < 0.05) and the rapidity dependence in a selected region of p t (1.0 < p t < 1.5 GeV/c) the minimal and maximal values are given.In addition to the total systematic uncertainties, also the individual contributions are shown.These are due to the procedures for particle identification ("PID"), the quality selection criteria applied to the tracks ("Track Quality"), the correction for inefficiencies due to high track densities ("Occupancy") and the effects of an azimuthally non-uniform detector acceptance ("Acceptance").
tance filter.This filter includes the gaps between the sectors for the support structures and in addition one entirely missing sector.Furthermore, also the effect of a non-uniform event-plane distribution was included.No significant differences between the input values v n and the ones extracted after filtering are observed, see left panel of Fig. 7.
Another systematic check is performed by analysing data that was recorded with a reversed magnetic field setting.In this configuration, the bending directions of positively and negatively charged particles are swapped such that they are measured by different areas in the outer two MDC layers, as well as TOF and RPC.No significant differences between the two settings are found, see right panel of Fig. 7.In addition, the analyses are also performed for each day of data-taking separately, in order to investigate whether any systematic trends appear in the course of the whole data-taking time.Also in this case, no deviations beyond the systematic uncertainty are observed.
Residual systematic effects can also be assessed by investigating whether the Fourier decomposition of the azimuthal particle distributions contains sine terms in addition to the cosine terms in Eq. ( 1).These are found to be of smaller or similar magnitude than the system-atic uncertainties estimated via the methods discussed above.Therefore, no additional systematic uncertainty is assigned due to these differences.
The main contribution to the global systematic uncertainty arises from the event-plane resolution.This is mainly caused by so-called "non-flow" correlations which can distort the event-plane determination.The magnitude of these systematic effects is evaluated using the three-sub-event method, i.e. by determining the event-plane resolution for combinations of different subevents separated in rapidity.It is found to be below 5 % for the centralities 10 − 40 % [36].While v 1 of protons is consistent with zero at midrapidity, it rises towards forward and decreases towards backward rapidities (see left top panel of Fig. 8).This  rapidity dependence is stronger at higher than at lower transverse momenta.The p t dependence of the proton v 1 is shown in the left top panel of Fig. 9 for four exemplary rapidity intervals.Its absolute value exhibits an almost linear rapid rise in the region p t < 0.6 GeV/c and then increases only moderately or even saturates for p t > 1 GeV/c.A comparison of the absolute v 1 values measured in forward and backward rapidity intervals, chosen symmetrically around mid-rapidity, results in an agreement well within systematic errors.
A p t and y cm dependence similar in shape is observed for the v 1 of deuterons and tritons (middle and right panels in the top row of Figs. 8 and 9).However, there are quantitative differences, namely that the saturation behaviour sets in at higher p t values (above p t ≈ 1.2 GeV/c for deuterons and p t ≈ 1.4 GeV/c for tritons) and reaches higher absolute values of v 1 (e.g.| ≈ 0.68 for the |y cm | interval 0.55−0.65).This also implies that the dependence of v 1 on rapidity gets more pronounced with increasing particle mass.

Elliptic flow: v 2
Figures 8 and 9 show in the second row a compilation of v 2 values for protons, deuterons and tritons as a function of p t and y cm .Their rapidity dependence is opposite to v 1 , i.e. the absolute value of v 2 is largest at mid-rapidity and decreases towards forward and backward rapidities for all investigated particles.The v 2 values around mid-rapidity decrease continuously with p t , and an indication for a saturation behaviour is seen at relatively high p t for protons only.The drop with p t is most pronounced for protons and gets weaker with increasing particle mass.Also the rapidity distribution of v 2 in a fixed p t interval strongly depends on the particle type.While for protons it reaches zero at rapidities of |y cm | ≈ 0.70, the distributions for deuterons is significantly narrower, such that it crosses zero already at |y cm | ≈ 0.50 and v 2 changes sign for larger centre-ofmass rapidities.For tritons this change of sign already happens around |y cm | ≈ 0.35.
The shape of the p t dependence for deuterons and tritons in the rapidity region, where a positive v 2 is observed, is clearly different to the one in the regions with negative v 2 .In the region y cm < −0.5 (y cm < −0.35), v 2 rises with p t for deuterons (tritons) towards a maximum, whose position seems to move towards higher p t with decreasing y cm , and then starts to drop again.

Higher flow harmonics: v 3 and v 4
In addition to the directed and elliptic flow coefficients also higher moments of the azimuthal distributions of particle emission relative to the reaction plane have been extracted.In ref. [22] data on flow coefficients up to the sixth order were presented for a limited region of phase-space.A systematic multi-differential analysis of higher orders over a larger p t − y cm range with satisfactory accuracy turned out to be possible only for the coefficients v 3 and v 4 .
Figures 8 and 9 exhibit in the third row the results on v 3 for protons, deuterons and tritons.The rapidity and p t dependences are similar for the three analysed particle species (see also Fig. 2 in ref. [22]).Generally, the rapidity dependence is comparable in shape to the one of v 1 , however, the v 3 values have opposite signs.Taking a closer look, one finds that the y cm distributions start with a steeper slope at midrapidity than v 1 and exhibit a turn around away from mid-rapidity.The positions of the corresponding maxima depend slightly on the particle mass and are found at |y cm | ≈ 0.5 (protons), ≈ 0.4 (deuterons) and ≈ 0.3 (tritons).Also, in distinction to v 1 , no clear evidence for a saturation is seen at high p t for v 3 (see Fig. 9).
Figures 8 and 9 present in the bottom row the results on v 4 for protons, deuterons and tritons.The rapidity distributions are similar in shape to the ones measured for v 2 for the corresponding particle, but have opposite signs.Also, they are narrower for v 4 than for v 2 and cross the v 1 = 0 line at smaller values of |y cm |.For the different particle species this is found to be at |y cm | ≈ 0.35 (protons), ≈ 0.3 (deuterons) and ≈ 0.25 (tritons).The increase of the absolute v 4 values with p t is also significantly less pronounced as in the case of v 2 .Therefore, in contrast to the case of v 2 , no saturation or even a maximum is observed at higher values of p t .

Centrality dependences
The directed flow at mid-rapidity can be quantified by its slope dv 1 /dy ′ | y ′ =0 which is defined relative to the scaled rapidity y ′ = y cm /y mid , with y mid = 0.74 as mid-rapidity in the laboratory system.It is determined as the linear term, dv 1 /dy ′ | y ′ =0 = a 1 , of a cubic ansatz v 1 (y ′ ) = a 1 y ′ + a 3 y ′ 3 which has been fitted to the measured data points.Similarly, the slope of v 3 dv 3 /dy ′ | y ′ =0 is extracted.The upper panels of Fig. 10 displays the corresponding values as determined for two different p t intervals (0.6 < p t < 0.9 GeV/c and 1.5 < p t < 1.8 GeV/c) and for the four centrality classes investigated in this analysis.The slope of v 1  exhibits no significant centrality dependence for all particles and p t intervals, except for the very central class where dv 1 /dy ′ is smaller than for the other centralities.This is distinctly different to the centrality dependence of the slope of v 3 , where the absolute value |dv 3 /dy ′ | is continuously increasing with centrality.Also, the values are almost identical for the different particles at all centralities, while for dv 1 /dy ′ a significant mass hierarchy is observed.
Similar to the triangular flow, also v 2 and v 4 at midrapidity depend on the reaction centrality.While the absolute value |v 2 | increases roughly linearly with centrality (see lower left panel of Fig. 10), v 4 exhibits a stronger centrality dependence (see lower right panel of Fig. 10).The mass hierarchy is visible for v 2 and v 4 in the lower p t interval for all centrality classes.In the higher p t region, however, only v 2 of tritons is different decreasing passage time of the colliding nuclei, which reduces the pressure transfer onto the light nuclei at higher energies.
Also v 2 at mid-rapidity exhibits a very distinct energy dependence.In the region 0.1 ≲ E beam /A ≲ 5 GeV v 2 is negative, i.e. the particle emission is out-of-plane, as the passage time of the spectator matter is long enough to cause the squeeze-out effect [9,57], i.e. the fireball pressure pushes particles preferentially into the direction which is not shadowed by spectators.At higher energies the passage times are too short and particle emission is in-plane as the pressure gradients are steepest in this direction.The integrated v 2 obtained for Au+Au collisions at √ s NN = 2.4 GeV in this analysis is in the region where out-of-plane emission is still very strong.It is also well in accordance with other measurements by EOS [42] and FOPI [38] (see right panel of Fig. 11).The p t dependence of v 2 at mid-rapidity measured by HADES is compared with results of other experiments in the same energy region (KaoS [56] and FOPI [38]) in Fig. 12. Within uncertainties and considering the slight differences of beam energies, good agreement with the other data sets is found.The new HADES data extend the phase-space coverage significantly in comparison to previous measurements with clearly improved accuracy.

Coalescence scaling
A comparison of the p t dependences of v 2 measured at mid-rapidity for protons, deuterons and tritons (see Fig. 13) demonstrates a clear mass ordering v prot. 2 2 for p t < 1.5 GeV/c.Within a naive nucleon-coalescence scenario one would expect that the observed flow coefficients scale with the nuclear mass number A according to the relation v n,A (A p t ) = A v n (p t ), where p t is the momentum of a single nucleon and p t,A = Ap t the one of the composite nuclei.The correspondingly scaled p t dependences of the proton v 2 are shown in Fig. 13 as solid curves for A = 2 and 3.The agreement with the measured v 2 values for deuterons and tritons is already quite good.However, this approximation only holds for small flow values and, as v 2 measured at high p t is quite sizeable, a correction term has to be taken into account [25,26]: as a function of p t at mid-rapidity (|ycm| < 0.05).The solid curves represent the proton distribution after scaling according to v n,A (A p t ) = A vn(p t ).The coloured bands depict the results as calculated for the higher-order nucleon-coalescence scenario given in Eq. ( 13).Systematic uncertainties are displayed as boxes.
In fact, the correspondingly scaled proton v 2 values agree well with the ones measured for deuterons and tritons up to the highest p t , as depicted by the coloured bands in Fig. 13.It should be noted though, that this kind of scaling is only observed in the region around mid-rapidity, where the elliptic flow is the predominant component of the azimuthal distributions.Similar to the case of v 2 (see Fig. 13) also for v 4 a nuclear mass scaling behaviour is observed at midrapidity.This is demonstrated in Fig. 14, which shows a comparison of v 4 at mid-rapidity as a function of p t for protons, deuterons and tritons.A mass hierarchy can be observed,

4
, at least in the lower p t region.In order to test whether this ordering of v 4 is also compatible with a nucleon coalescence scenario, we use an extension of Eq. ( 13) which takes combinations of different orders into account2 : If the higher-order correction is omitted, this results in the simple approximation of v 4,A (A p t ) = A 2 v 4 (p t ), which therefore should only be valid for small flow values.Figure 14 includes a comparison of these different approximations to the data.While the relation given in Eq. ( 13) does not provide a good match with the data, its extended version given in Eq. ( 15) results in a very good description of the deuteron and triton data.Also, the simple relation v 4,A (A p t ) = A 2 v 4 (p t ) is quite close to the data points, indicating that the higher-order corrections are small.While the above discussed scaling properties can be indicative for nucleon coalescence as the main process responsible for light nuclei formation, a more refined discussion would involve the comparison to various models.Examples for implementations of the coalescence approach within transport models to describe HADES and STAR flow data, respectively, can be found Table 4 Parameters describing the initial nucleon distribution for the different centrality classes as calculated within the Glauber-MC approach [29].Listed are the corresponding average impact parameter ⟨b⟩ and the average participant eccentricities ⟨ϵ 2 ⟩ and ⟨ϵ 4 ⟩.
in [19,58].These studies should be extended in the future in a more systematical way using several transport models in order to arrive on firmer conclusions on this topic.

Initial eccentricity
In order to investigate to what extent the spatial distribution of the nucleons in the initial state of the collision system determines the observed flow pattern, we use the eccentricity ϵ n of order n of the participant nucleon distribution in the transverse plane as calculated within the Glauber-MC approach [29,59]: with r = x 2 + y 2 , ϕ = tan −1 (y/x) and x, y as the nucleon coordinates in the plane perpendicular to the beam axis, where x is oriented in the direction of the impact parameter.The values calculated for the different centrality classes are given in Tab. 4.
Figure 15 shows the elliptic flow measured at midrapidity for all three investigated particle species after dividing it by the event-by-event averaged second-order participant eccentricity, v 2 /⟨ϵ 2 ⟩.Remarkably, this scaling results in almost identical values for all centrality classes at high transverse momenta, indicating that the centrality dependence of the elliptic flow of particles emitted at early times is to a large extent already determined by the initial nucleon distribution.However, as the elliptic flow at these beam energies is due to the so-called squeeze-out effect, caused by the passing spectators, it is not immediately clear how the flow pattern can be directly related to the initial participant distribution.A possible explanation might be that the distribution of the spectators forms a negative image of the one of the participants and thus could imprint its shape onto the emission pattern of the light nuclei.The scaling of v 2 works less well at lower p t , which suggests that particles emitted at later times are less affected by the initial state geometry.Also, we observe a scaling of v 4 with ϵ 2 2 , as depicted in the left panel of Fig. 16 which presents v 4 /⟨ϵ 2 ⟩ 2 for different centralities in two transverse momentum intervals.This points to a fixed relation between v 2 and v 4 , such that the latter is a second order correction ∝ ϵ 2 2 to the overall emission pattern defined at mid-rapidity by v 2 .This is contrary to the case at very high collision energies, where higher-order flow coefficients are related to initial state fluctuations and thus, to a large extent, should be independent of one another.In this scenario one would also expect v 4 to scale rather with ϵ 4 .While this might be observable also here at lower p t , v 4 /⟨ϵ 4 ⟩ is not independent of centrality in the high p t region, i.e. for particles emitted at early times, as demonstrated in the right panel of Fig. 16.

Model comparisons
In the following, several transport model calculations are compared with the measured flow data.These models provide the possibility to test the effect of the EOS of dense nuclear matter on the flow coefficients by implementing different density dependent potentials.Usually, these are parameterised such that the dependence on the baryon density ρ results in either a weak ("soft EOS") or a strong ("hard EOS") repulsion of compressed nuclear matter.The comparison with data then allows for a discrimination between these two scenarios.While previous investigations were only based on measurements of the directed and elliptic flow [60,61], the  information from higher-order flow coefficients will provide additional discriminating power.Ultimately, the multi-differential high-statistics data presented here should enable a direct extraction of the EOS parameters via a Bayesian fit of the models to the data.However, as a prerequisite it is important to establish that the various model approaches do not differ significantly in their predictions in order to allow for a consistent determination of the EOS.For a detailed review of the different approaches used for transport simulations see [62][63][64].As examples the predictions by two QMD models, JAM 1.9 [65] and UrQMD 3.4 [18], and one BUU model, GiBUU 2019 [66] are considered here.The JAM code is used with three different EOS implementations: hard momentum independent NS3, hard momentum dependent MD1 and soft momentum dependent MD4.The UrQMD code is employed with the "hard EOS", and GiBUU with the "soft EOS" (Skyrme 12).
Comparisons of the model predictions to the proton flow of different order measured at mid-rapidity are presented, as a function of centrality, in Fig. 17 (low p t interval 0.6 < p t < 0.9 GeV/c) and Fig. 18 (high p t interval 1.5 < p t < 1.8 GeV/c).As most models do not include a dedicated mechanism for the generation of light clusters, which would be needed for a realistic prediction of deuteron and triton flow, we restrict the comparison here to protons only.
Generally, all models roughly capture the overall magnitude and trend of the measured data.In the lower p t region the differences between the models are relatively small.JAM with MD4 provides the overall best reproduction of the data points, with the exception of  v 2 where MD1 is closer to the data.UrQMD is close to the data for v 2 , but deviates for v 1 , v 3 and v 4 at several centralities, while GiBUU reproduces generally better the exception of v 4 .In the higher p t interval the deviations between the models are a bit larger.
Here JAM with MD1 yields the best match to the data, while MD4 and NS3 do not provide a consistent description of the measurements.Also, for UrQMD and GiBUU, systematic deviations are observed for some orders of the flow coefficients.Nevertheless, the models presented here should form a useful basis for further, more detailed data comparisons and and consistent determination of the EOS.It should be noted that not all model calculations include the effects of momentum and isospin dependent potentials, which, however, will be essential for this purpose.Furthermore, a common treatment of cluster formation should be implemented which will allow for an usage of the data on deuteron and triton flow as an additional constraint.

Conclusions
In summary, we have presented a detailed multi-differential measurement of collective flow coefficients of protons, deuterons and tritons in Au+Au collisions at √ s NN = 2.4 GeV using the high-statistics data set collected with the HADES experiment.The directed (v 1 ), elliptic (v 2 ) and higher order (v 3 and v 4 ) flow coefficients were determined with respect to the first-order event-plane measured at projectile rapidities.The centreof-mass energy of √ s NN = 2.4 GeV is close to the region where dv 1 /dy ′ | y ′ =0 is maximal and v 2 is minimal.All flow coefficients were extracted as a function of transverse momentum p t and rapidity y cm over a large region of phase-space and in four centrality classes.The p t and y cm dependences of v 1 are very similar in shape for protons, deuterons and tritons.A clear mass hierarchy is observed for v 1 values measured away from mid-rapidity at higher p t , as well as for dv 1 /dy ′ | y ′ =0 , which both increase with the mass of the particle.The elliptic flow coefficient v 2 has a Gaussian shaped rapidity distribution, whose width narrows with increasing particle mass, such that the rapidity value for which v 2 changes sign moves closer towards mid-rapidity for increasing mass number.Both, the proton directed flow dv 1 /dy ′ | y ′ =0 and the elliptic flow v 2 at mid-rapidity are in line with the established energy dependence.The p t and y cm dependences of v 3 (v 4 ) are relatively similar to the one of v 1 (v 2 ), but have the opposite sign.At mid-rapidity a nucleon number scaling is observed in the p t dependence of v 2 (v 4 ), when dividing the values of v 2 (v 4 ) by A (A 2 ) and p t by A. This might be indicative for nuclear coalescence as the main process responsible for light nuclei formation.Such a straightforward scaling is not seen at more forward and backward rapidities.The elliptic flow measured at mid-rapidity at higher p t is found to be independent of centrality for all three investigated particle species after dividing it by the event-by-event averaged second order participant eccentricity v 2 /⟨ϵ 2 ⟩.A similar scaling is observed for v 4 after division by ϵ 2 2 .The new multi-differential high-precision data on v 1 , v 2 , v 3 , and v 4 provides important constraints on the equation-of-state of compressed baryonic matter as used in models of relativistic nuclear collisions [4].In particular, the higher moments provide more discriminating power than the directed and elliptic flow alone.The general features of the data on proton flow at midrapidity are qualitatively captured by several transport models.A consistent and exact description of all flow coefficients over the whole phase-space and at all investigated centralities is not yet possible.With further progress in the theoretical developments it should be feasible to use the data shown here, together with other measurements, to directly extract a precise parametrization of the equation-of-state of compressed nuclear matter.

1
Cross section of the HADES set-up during the measurement of Au+Au collisions at √ s NN = 2.4 GeV.Shown is the arrangement of the different detectors and a magnet coil on one side of the beam pipe.

4 Flow coefficients 4 . 1 1 Figures 8
Figures 8 and 9  present in the uppermost row an overview on the directed flow coefficient v 1 measured for protons, deuterons and tritons in various p t and y cm intervals in semi-central (20 − 30 %) Au+Au collisions.While v 1 of protons is consistent with zero at midrapidity, it rises towards forward and decreases towards backward rapidities (see left top panel of Fig.8).This

Fig. 8
Fig. 8 The flow coefficients v 1 , v 2 , v 3 , and v 4 (from top to bottom panels) of protons, deuterons and tritons (from left to right panels) in semi-central (20 − 30 %) Au+Au collisions at √ s NN = 2.4 GeV as a function of the centre-of-mass rapidity ycm in transverse momentum intervals of 50 MeV/c width.Systematic uncertainties are displayed as boxes.Lines are to guide the eye.

Fig. 9 1 .
Fig.9The flow coefficients v 1 , v 2 , v 3 , and v 4 (from top to bottom panels) of protons, deuterons and tritons (from left to right panels) in semi-central(20 − 30 %) Au+Au collisions at √ s NN = 2.4 GeV as a function of p t in several rapidity intervals chosen symmetrically around mid-rapidity.The values measured in the forward hemisphere (open symbols) have been multiplied by −1.Systematic uncertainties are displayed as boxes.Lines are to guide the eye.

Fig. 15
Fig. 15 Scaled elliptic flow of protons, deuterons and tritons in two transverse momentum intervals at mid-rapidity in Au+Au collisions at √ s NN = 2.4 GeV for four centralityclasses.The values are divided by the second order eccentricity, v 2 /⟨ϵ 2 ⟩, as calculated within the Glauber-MC approach for the corresponding centrality range (see Tab. 4).Systematic uncertainties are displayed as boxes.

Fig. 16
Fig.16 Same as in Fig.15, but for the scaled quadrangular flow.The values are divided by the square of second order eccentricity, v 4 /⟨ϵ 2 ⟩ 2 (left panel), and by the fourth order eccentricity, v 4 /⟨ϵ 4 ⟩ (right panel).

Fig. 17
Fig. 17Directed dv 1 /dy ′ | y ′ =0 (left top panel), elliptic v 2 (right top panel), triangular dv 3 /dy ′ | y ′ =0 (left bottom panel) and quadrangular v 4 (right bottom panel) flow of protons in the transverse momentum interval 0.6 < p t < 0.9 GeV/c at mid-rapidity in Au+Au collisions at √ s NN = 2.4 GeV for four centrality classes.The data are compared to several model predictions (see Fig. 17Directed dv 1 /dy ′ | y ′ =0 (left top panel), elliptic v 2 (right top panel), triangular dv 3 /dy ′ | y ′ =0 (left bottom panel) and quadrangular v 4 (right bottom panel) flow of protons in the transverse momentum interval 0.6 < p t < 0.9 GeV/c at mid-rapidity in Au+Au collisions at √ s NN = 2.4 GeV for four centrality classes.The data are compared to several model predictions (see text for details).The width of the bands reflect the statistical uncertainties of the model calculations.

Fig. 18
Fig.18 Same as in Fig.17, but for the transverse momentum interval 1.5 < p t < 1.8 GeV/c.
Fig.18 Same as in Fig.17, but for the transverse momentum interval 1.5 < p t < 1.8 GeV/c.

Table 2
List of applied track selection and PID criteria.In addition to the nominally applied values, also two variations are given for each selection criterion, which are used for the determination of the systematic uncertainties.

Table 3
Summary of the contributions to the absolute systematic uncertainties for different particle species and flow coefficients of different orders.For the p t dependence in different regions of rapidity (v 1 and v 3 : −0.25 < ycm < −0.15, v 2 and v 4 : ) Quadrangular flow v 4 of protons, deuterons and tritons in semi-central (20 − 30 %) Au+Au collisions at √ s NN = 2.4 GeV as a function of p t at mid-rapidity (|ycm| < 0.05).The dashed curves represent the proton distribution after scaling according the higher order nucleon-coalescence scenario given in Eq. (13).The coloured bands depict the results as calculated with Eq. (15) which includes the additional contribution of v 2 assuming the relationv 2 = − √ 2 v 4 .The solid curves show the result for the approximation v n,A (A p t ) = A 2 vn(p t ).Systematic uncertainties are displayed as boxes.