Probing Strangeness Canonical Ensemble with $K^{-}$, $\phi(1020)$ and $\Xi^{-}$ Production in Au+Au Collisions at ${\sqrt{s_{\rm NN}} = \rm{3\,GeV}}$

We report the first multi-differential measurements of strange hadrons of $K^{-}$, $\phi$ and $\Xi^{-}$ yields as well as the ratios of $\phi/K^-$ and $\phi/\Xi^-$ in Au+Au collisions at ${\sqrt{s_{\rm NN}} = \rm{3\,GeV}}$ with the STAR experiment fixed target configuration at RHIC. The $\phi$ mesons and $\Xi^{-}$ hyperons are measured through hadronic decay channels, $\phi\rightarrow K^+K^-$ and $\Xi^-\rightarrow \Lambda\pi^-$. Collision centrality and rapidity dependence of the transverse momentum spectra for these strange hadrons are presented. The $4\pi$ yields and ratios are compared to thermal model and hadronic transport model predictions. At this collision energy, thermal model with grand canonical ensemble (GCE) under-predicts the $\phi/K^-$ and $\phi/\Xi^-$ ratios while the result of canonical ensemble (CE) calculations reproduce $\phi/K^-$, with the correlation length $r_c \sim 2.7$\,fm, and $\phi/\Xi^-$, $r_c \sim 4.2$\,fm, for the 0-10\% central collisions. Hadronic transport models including high mass resonance decays could also describe the ratios. While thermal calculations with GCE work well for strangeness production in high energy collisions, the change to CE at $\rm{3\,GeV}$ implies a rather different medium property at high baryon density.

We report the first multi-differential measurements of strange hadrons of K − , φ and Ξ − yields as well as the ratios of φ/K − and φ/Ξ − in Au+Au collisions at √ sNN = 3 GeV with the STAR experiment fixed target configuration at RHIC. The φ mesons and Ξ − hyperons are measured through hadronic decay channels, φ → K + K − and Ξ − → Λπ − . Collision centrality and rapidity dependence of the transverse momentum spectra for these strange hadrons are presented. The 4π yields and ratios are compared to thermal model and hadronic transport model predictions. At this collision energy, thermal model with grand canonical ensemble (GCE) under-predicts the φ/K − and φ/Ξ − ratios while the result of canonical ensemble (CE) calculations reproduce φ/K − , with the correlation length rc ∼ 2.7 fm, and φ/Ξ − , rc ∼ 4.2 fm, for the 0-10% central collisions. Hadronic transport models including high mass resonance decays could also describe the ratios. While thermal calculations with GCE work well for strangeness production in high energy collisions, the change to CE at 3 GeV implies a rather different medium property at high baryon density.

I. INTRODUCTION
Relativistic heavy ion physics is aiming at the detailed investigation of phase structures of strongly interacting matter, governed by quantum chromodynamics (QCD), under extreme conditions of high temperature and density [1][2][3]. Particle production has been studied to investigate properties of the produced QCD matter in heavy-ion collisions. The strange quark mass is comparable to the QCD scale (Λ QCD ∼ 200 MeV), therefore strange quark plays an important role in studying the QCD phase diagram and the Equation-of-State (EoS), particularly in the high density region [4][5][6][7][8][9].
Statistical thermal models have often been used to characterize thermal properties of the produced media [10][11][12][13][14][15][16][17][18][19][20]. In these models, grand canonical ensemble (GCE) and canonical ensemble (CE) statistical descriptions can be applied to conserve electric charge, baryon number, and strangeness number in order to compute the final state particle yields. Both GCE and CE models are able to describe various particle yields including strange particles produced in heavy-ion collisions at RHIC and the LHC at center-of-mass energy ( √ s NN ) greater than 7.7 GeV. It has been argued that at lower energies, strangeness number needs to be conserved locally on an event-by-event basis described by the CE, which leads to a reduction in the yields of hadrons with non-zero strangeness number ("Canonical Suppression") [10,21,22], but not for the φ(1020) meson with zero net strangeness number (S=0). The φ/K − ratio is expected to increase with decreasing collision energy in models using the CE treatment for strangeness, opposite to the trend in the GCE treatment. The canonical suppression power for Ξ − (S=2) is even larger than for K − (S=1). The φ/K − and φ/Ξ − ratios offer a unique test to scrutinize thermodynamic properties of strange quarks in the hot and dense QCD environment. In heavy-ion collisions, the near/sub-threshold production of multi-strange hadrons can be achieved from the multiple collisions of nucleons, produced particles, and short-lived resonances [23]. The particle production in heavy-ion collisions below its free nucleon-nucleon (NN) threshold ( √ s NN ∼2.89 GeV for φ and ∼3.25 GeV for Ξ − ) is expected to be sensitive to the stiffness of the nuclear EoS at high density [24], as it is for single-strange hadrons [6,7]. The near/sub-threshold production further provides the possibility to observe exotic states of QCD matter [25] and signatures of "soft deconfinement" [26]. Previous measurements show that the φ/K − ratio in heavyion collisions stays remarkably flat (∼0.15) at collision ener- * Deceased gies √ s NN > 5 GeV [27][28][29]. At collision energies close to or below the φ and Ξ NN-thresholds, recent measurements of φ/K and φ/Ξ ratios from HADES and FOPI have achieved a significance about 2.2-3.8 sigma in central heavy ion collisions, and the results indicate a relative enhancement compared to those at high energies [30][31][32][33], indicative of the applicability of the CE description for strangeness production at these energies. Measurements from π or proton induced nuclear reactions [34,35] suggest that absorption in cold nuclear matter may play a role in the K − and φ production yields in nuclear collision at low energies. In this Letter, we report high precision measurement of φ/K − and φ/Ξ − ratios in Au+Au collisions at √ s NN = 3 GeV from the STAR experiment.

II. EXPERIMENT AND DATA ANALYSIS
The dataset used in this analysis was collected under the fixed target (FXT) setup [36] in the 2018 RHIC run. A single beam was provided by RHIC with total energy equal to 3.85 GeV/nucleon and incident on a gold target of thickness 0.25 mm, corresponding to a 1% interaction probability. The target is installed inside the vacuum pipe, 2 cm below the center of the beam axis, and located 200 cm to the west of the center of the STAR detector. The main detectors used are the Time Projection Chamber (TPC) [36,37], the Time of Flight (TOF) detector [36,38], and the Beam-Beam Counter (BBC) [39]. The trigger is provided by the signal in the east BBC detector and at least five hits in the TOF detector. To best utilize the detector band-width, the beam-on-target collision rate was tuned to around 1.5 kHz, and the pileup contribution to the triggered event is < 1% [40]. Tracking and particle identification (PID) are done using the TPC and TOF. Both the TPC and TOF detectors have full azimuthal coverage within a pseudorapidity range of 0< η < 1.88 for the TPC and 0< η < 1.5 for the TOF in FXT mode [36][37][38]. Events are selected with the offline reconstructed collision vertex within 1.5 cm of the target center along the beam direction. Approximately 2.6×10 8 minimum bias (MB) triggered events passed the selection criteria and are used in this analysis.
The centrality class is selected using measured charged particle total multiplicity within the TPC acceptance. A Monte Carlo Glauber model, used in conjunction with a negative binomial distribution to model particle production in hadronic collisions, is optimized in order to best match the data and determine the centrality class [40,41]. Due to the trigger inefficiency in the low multiplicity region (corresponding to the most peripheral collisions), we only report the results from the 0-60% centrality class in this paper. In addition, in order to reduce the pile-up contamination, events above the reference multiplicity of 195 are removed from the most central centrality class.
φ mesons are reconstructed via the decay channel φ → K + K − with a branching ratio (BR) of (49.2 ± 0.5)% , while the Ξ − hyperons decay via Ξ − → Λπ − → pπ − π − with a BR of (63.8 ± 0.5)% [42]. Ξ − reconstruction is performed using the KFParticle package based on the Kalman Filter method [44,46]. The charged tracks are reconstructed with the TPC in a 0.5 T uniform magnetic field, and are required to consist of at least 20 TPC hits (out of a maximum of 45) and have a ratio between the number of hit points and the maximum possible number of hit points larger than 0.52 to ensure good tracking and avoid track splitting. The TPC tracking performance with these requirements in the FXT data is similar to that in other data taken in the collider mode. Monte Carlo simulations also reproduces the distributions of various tracking variables. The charged tracks are identified via a combination of the ionization energy loss (dE/dx) measurement with the TPC and the time-of-flight (tof ) measurement with the TOF [47,48]. The resolution-normalized dE/dx or β deviation from the expected values are used for the PID selection. A minimum p T cut of 0.2 GeV/c is required in the analysis. Since the K − /π − ratio is much smaller than the K + /π + ratio at low energies, to reduce the contamination from π − and e − tracks, a strict PID criterion for K − is implemented by requiring the TPC dE/dx and TOF β to be within three standard deviations of the expected values. K + tracks used for the φ analysis are selected with a hybrid algorithm, in which the TPC dE/dx requirement is applied at low momentum p < 0.5 GeV/c while an additional TOF β requirement is imposed at p > 0.5 GeV/c. In the Ξ analysis, proton and π − tracks are identified by requiring the TPC dE/dx to be within three standard deviations of the expect values and the TOF β requirement is only applied when there is a valid measurement. Figure 1 (a) shows the invariant mass distribution of K + K − pairs in the transverse momentum (p T ) region of 0.4-1.6 GeV/c for 0-60% central collisions. The combinatorial background is estimated with the mixed-event (ME) technique in which K + and K − from different events of similar characteristics (centrality, event plane angle) are paired. The mixed-event spectra are normalized to the same-event (SE) distributions in the mass range of 1.04-1.08 GeV/c 2 . After the subtraction of the combinatorial background, the remainder distribution, shown as red solid circles, is fitted with a Breit-Wigner function for the signal plus a linear function which represents the remaining correlated background (< 1%) from a partial reconstruction of strange hadrons. The φ meson raw yields are extracted from the Breit-Wigner function fit within the corresponding ±3Γ mass window (mean value µ ∼ 1.0194 GeV/c 2 , full-width-at-half-maximum Γ ∼ 6.5 MeV/c 2 ). The extracted φ signal shape is consistent with its intrinsic properties convoluted with the detector smearing effect due to finite momentum resolution (< 3% for single track). Note that a Voigt function has been used to extract the signal counts as a cross check, and the extracted yields are consistent with the default value within uncertainties. in the p T region of 0.5-2.0 GeV/c for 0-40% central collisions. The combinatorial background is estimated with the rotating daughter (Rot) method, in which a daughter track of Ξ − is rotated by a random angle between 150 to 210 degrees in the transverse plane. The rotated spectra are normalized to the same-event distributions in the mass ranges of 1.30-1.31 and 1.34-1.35 GeV/c 2 . After the combinatorial background is subtracted, the Λπ − invariant mass distribution is fitted with a Gaussian for the signal plus a linear function for the remaining correlated background (< 1%). The Ξ − raw yields are obtained via histogram bin counting from the invariant mass distributions with all background subtracted within mass windows of ±3σ (µ ∼ 1.3222 GeV/c 2 , Gaussian width σ ∼ 2.0 MeV/c 2 ). The signal-to-background (S/B) ratio for φ and Ξ − within the reconstructed mass windows is about 0.7 and 0.6 respectively. The reconstructed φ and Ξ − acceptances (p T vs. y cm ) in the collision center-of-mass frame are shown in Fig. 1 (c) and (d), respectively. The target is located at y cm = −1.05, using the convention where the beam travels in the positive direction. The red curve represents the TPC and TOF acceptance edge.
Particle raw yields are calculated in each centrality and p T bin within each rapidity slice. The raw yields are corrected for the TPC acceptance and tracking efficiency, the particle identification efficiency, and the TOF matching and PID efficiency. The TPC acceptance and tracking efficiency is obtained using the standard STAR embedding technique [29,49], in which a small number of MC tracks are processed through the GEANT (v3.21) simulation [50], then mixed with the real data and reconstructed using the same algorithm as in the real data. The TPC PID, TOF matching and PID efficiencies are obtained from the data-driven method similar as in Ref. [51]. The final average reconstruction (including acceptance etc.) efficiency is ∼0.30, 0.04, and 0.02 for K − , φ and Ξ − , respectively. MC embedding simulation also reproduces various topological variables used in the Ξ − reconstruction. As a cross-check, we conducted the measurement of Ξ − lifetime from the same data and the result is 164.2 ± 6.6 (stat.) ps, consistent with the PDG value, 163.9 ± 1.5 ps. The corrected p T spectra in symmetric rapidity bins (-0.2,0) vs. (0,0.2) are also consistent.
The systematic uncertainty of the raw yield extraction is estimated by changing the histogram fitting method to bin counting method or by changing the fitting ranges. The maximum difference between these scenarios and the default one is considered as one standard deviation. The contribution varies by p T , rapidity, and centrality and the overall contribution is less than 5% for the invariant yield. The systematic uncertainty in the TPC acceptance and efficiency correction ε TPC is estimated by varying the cuts on track selection criteria [49] and topological variables (for Ξ − only). The contribution to the total yield is 4-5% for K − , 13-16% for φ and 6-10% for Ξ − . This leads to a 10-13% (12-18%) uncertainty in the measured φ/K − (φ/Ξ − ) ratio. The uncertainty of the PID efficiency correction is estimated by varying the PID selection cuts and the contribution is less than 3% to the total yield. For the p T integrated yield, the uncertainty due to the extrapolation to the full p T range is estimated by choosing several fitting functions including Levy, Blast-Wave, m T -exponential, p Texponential [52], and the maximum difference between these scenarios and the default one (m T -exponential) is quoted as one standard deviation. This contribution is 5-7% for K − , 14-17% for φ and 13-15% for Ξ − , respectively. This measurement covers nearly the full rapidity range from y=0 to the target region. The systematic uncertainty due to the rapidity coverage extrapolation is negligible compared to other systematic sources. For each individual φ-meson, K − and Ξ − measurement, some of the uncertainties are correlated or partially correlated (e.g. TPC and PID). To avoid the correlation in the ratio measurement, we vary the above cuts simultaneously for φ, K − and Ξ − , then quote the difference in the final ratios as the systematic uncertainties. Figure 2 and 3 show the acceptance × efficiency corrected K − , φ and Ξ − invariant yields as a function of m T − m 0 (m T = m 2 0 + p 2 T /c 2 , where m 0 is particle rest mass, and c is the speed of light) for various rapidity ranges in 0-10%, 10-40% and 40-60% Centrality Au+Au collisions at √ s NN = 3 GeV. Dashed and solid lines depict fits to the spectra with the m T -exponential function in order to extrapolate to the unmeasured p T ranges (∼20-40% for K − which vary rapidity, ∼33-50% for φ and ∼40-60% for Ξ − ). The fitted inverse slope parameters indicate harder spectra for the φ-mesons compared to the K − and Ξ − within uncertain- ties. The inverse slope parameters gradually decrease from mid-rapidity to forward/backward rapidity and follow the T eff / cosh(y) distribution well. The inverse slope parameter at y = 0, T eff , is extracted to be 177 ± 5(stat) ± 8(sys) MeV for φ meson, 158 ± 3(stat) ± 3(sys) MeV for K − and 156 ± 3(stat) ± 24(sys) MeV for Ξ − in 0-10% central collisions. This agrees with the collision energy dependence trend from other experiments [28,33]. Table I lists the extracted T eff parameter for these particles in different centrality bins from this measurement. The p T integrated rapidity distributions dN/dy are displayed in Fig. 4 for Au+Au collisions at √ s NN = 3 GeV for three different centralities. Solid curves depict Gaussian function fits to the data points with the centroid parameter fixed to zero. They are used to extrapolate to the unmeasured rapidity region (∼5% for K − , ∼9% for φ and ∼6% for Ξ − ) for calculating total multiplicities. The integral yields follow the collision energy trend from other experiments and drop quickly toward the low energies around threshold [28,[30][31][32][33]. The φ/K − and φ/Ξ − ratios are presented in Fig. 5 as a function of collision energy √ s NN , including the midrapidity data in central Au+Au or Pb+Pb data from the AGS, SPS and RHIC BES at higher energies and 4π acceptance data from SIS at lower energies. The black solid circles show our measurements in the 0-10% centrality bin in Au+Au collisions at √ s NN = 3 GeV. The measured φ, K − and Ξ − yields in 4π and the φ/K − , φ/Ξ − ratios in different centrality bins are listed in Tab. II. The φ/K − and φ/Ξ − ratios measured at 3 GeV are slightly higher than, or comparable to, the values at high energies for √ s NN 5 GeV [27][28][29][53][54][55][56][57][58] despite the collision energy being very close to the φ threshold and below the Ξ − threshold in NN collisions. Note that the enhancement of φ/K − and φ/Ξ − were also observed at lower collision energies in √ s NN = 2.4 GeV Au+Au [33] and √ s NN = 2.6 GeV Ar+KCl collisions [30,57], respectively.

III. RESULTS AND DISCUSSIONS
Various curves in Fig. 5 represent the predictions of φ/K − and φ/Ξ − ratios from several model calculations in central A+A collisions. Statistical model calculations, based on the Grand Canonical Ensemble and Canonical Ensemble for strangeness with several different choices of strangeness correlation length (r c ), were calculated using the THERMUS package [59] with energy dependent freeze-out parameters (chemical freeze-out temperature T ch , baryon chemical potentials µ B ) taken from [20], for instance, T ch = 72.9 MeV and µ B = 701.4 MeV for √ s NN = 3 GeV. We noted that the φ/K − and φ/Ξ − ratios from GCE depend on strangeness chemical potential, µ S . From the results of the thermal model fit to the STAR BES-I data [29], there is an empirical relation µ S = µ B /4 in the collision energy region between 7.7 -39 GeV. The same relation was assumed and used in the GCE calculation at lower energies presented in Fig. 5. Our measured φ/K − and φ/Ξ − ratios are larger than this GCE calculation: χ 2 /ndf = 26.0/2 (p-value < 1e −5 ), which indicates the event-by-event strangeness conservation is crucial [17] in such collisions. The exact GCE calculation depends on the precise determination of T ch , µ B , µ S etc, which can be extracted through a global fit to various other particle yields at 3 GeV. In the canonical approach, the correlation length, r c , defines a region of the particle production phase space inside which the production of the strangeness is canonically conserved. Both the φ/K − and φ/Ξ − data from our measurement favor the CE thermodynamics for strangeness with a small strangeness correlation length (r c ∼ 2.7 fm for φ/K − and r c ∼ 4.2 fm for φ/Ξ − ). It is worthwhile to point out that the CE calculations with the same r c parameter cannot describe our φ/K − and φ/Ξ − data simultaneously. The CE calculation with r c ∼ 4.2 fm describes φ/Ξ − well while it deviates by about 3.5σ for φ/K − . r c is an approximation in the CE for reproducing the strange production in heavy-ion collisions. It is unclear if the same value of r c should fit for TABLE I. Inverse slope parameter T eff at y = 0 for the mT spectra of φ, K − , Ξ − in Au+Au collisions at √ sNN = 3 GeV. The first error given corresponds to the statistical one, the second to the systematic error.
The first error given corresponds to the statistical one, the second to the systematic error. Previous measurements from smaller collision systems (Ar+KCl and Al+Al collisions) show comparable or higher φ/K − and/or φ/Ξ − ratios at energies below 3 GeV [30][31][32]57]. The exclusive measurement in p+p collisions at 2.7 GeV shows a much larger φ/K − ratio (1.04 ± 0.23) [64], while the measured ratio at 17.3 GeV (0.11 ± 0.01) [65,66] is comparable to that in central Au+Au/Pb+Pb collisions at similar energies [28,29]. The φ/Ξ − ratio in p+p collisions at 17.3 GeV [66,67], 5.09 ± 0.36, is also significantly larger than that in central Au+Au/Pb+Pb collisions [28,29,55]. In our measurement at 3 GeV, there is no obvious difference in the φ/K − ratio between the 0-10% and 10-40% central bins, while the result in the most peripheral 40-60% central bin shows a hint of a larger value, as shown in Tab. II. Similarly, the φ/Ξ − ratio in mid-central collisions seems to be larger than that in central collisions. Overall, these observations are qualitatively consistent with the expectation that a smaller canonical volume in the smaller system leads to a higher observed φ/K − and/or φ/Ξ − ratio.
Hadronic transport models are widely used in the high baryon density region to study the properties of the produced dense matter [60-63, 68, 69]. In the modified version of the Ultra-relativistic Quantum Molecular Dynamics (UrQMD) model [60], UrQMD 2 , new decay channels from high mass baryon resonances to φ and Ξ − are deployed. The relevant decay branching fraction was determined by fitting the experimental data from p+p collisions [64]. From the comparison shown in Fig. 5, the modified UrQMD 2 calculation for central (b < 5 fm) Au+Au collisions agrees with the data points at low √ s NN , including our new measurement for φ/K − . However calculations from the public UrQMD 1 model [62,63] underestimate our measurements for both φ/K − and φ/Ξ − . The SMASH (Simulating Many Accelerated Strongly-interacting Hadrons) model [61] attempts to incorporate the newest available experimental data to constrain the resonance branching ratios. These data include both elementary hadronic cross sections and dilepton invariant mass spectra. The φ/K − ratio is reasonably reproduced using SMASH in the smaller system and √ s NN below 3 GeV, despite the overestimation of each individual (φ, K − ) transverse mass spectrum measured, e.g. in Au+Au 0-40% system by HADES [33,61]. The predicted φ/K − ratio from the same model is about 2.5σ higher than central Au+Au 0-10% collisions at 3 GeV. This indicates that some important in-medium mechanism for strangeness production and propagation may be missing for the large system in SMASH. Both UrQMD and SMASH calculations reproduced the measured strangeness data highlighting the importance of the contributions of the resonances in the low energies. Furthermore, the φ-meson scattering with the baryonic medium remains an open question from recent measurements of π induced nucleus reactions and the p-φ femtoscopy [35,70]. More detailed investigations are needed in order to understand the dynamics of strange and multi-strange hadrons at low energy nuclear collisions. Our measurement of K − , φ and Ξ production yields in 3 GeV Au+Au collisions demonstrates the necessity of the Canonical Ensemble for strangeness at low energy heavy-ion collisions. In the meantime, hadronic transport model calculations (UrQMD and SMASH) including resonance contributions reproduce the data. These observations suggest a change of the medium properties at 3 GeV compared to those from higher energy collisions. Similar conclusions have been reached from the measurements of collectivity [71] and high moment of protons [40] in 3 GeV Au+Au collisions. The solid black circles show the measurements presented here in 0-10% centrality bin, while empty markers in black are used for data from various other energies and/or collision systems [27][28][29][30][31][32][33]57]. The vertical grey bands on the data points represent the systematic uncertainties. The grey solid line represents a THERMUS calculation based on the Grand Canonical Ensemble (GCE) while the dotted lines depict calculations based on the Canonical Ensemble (CE) with different values of the strangeness correlation radius (rc) [20,59]. The green dashed line, green shaded band and the solid red line show transport model calculations from the public versions UrQMD 1 [62,63], modified UrQMD 2 [60] and SMASH [61], respectively.

IV. SUMMARY
In summary, we report the systematic measurements of K − , φ(1020) and Ξ − production yields and the φ/K − , φ/Ξ − ratios in Au+Au collisions at √ s NN = 3 GeV with the STAR experiment at RHIC. The measured φ/K − ratio is significantly larger than the statistical model prediction based on Grand Canonical Ensemble in the 0-10% central collisions. Both the results of φ/K − (r c ∼ 2.7 fm) and φ/Ξ − (r c ∼ 4.2 fm) ratios favor the Canonical Ensemble model for strangeness production in such collisions. Transport models, including the resonance decays, could reasonably describe our measured φ/K − ratio at 3 GeV and the increasing trend of φ/Ξ − at lower energies. The new results from this paper suggest a significant change in the strangeness production for √ s NN < 5 GeV, providing new insights towards the understanding of the QCD medium properties at high baryon density.
V. ACKNOWLEDGEMENT