Constraints on directionality effect of nuclear recoils in a liquid argon time projection chamber

The direct search for dark matter in the form of weakly interacting massive particles (WIMP) is performed by detecting nuclear recoils produced in a target material from the WIMP elastic scattering. The experimental identification of the direction of the WIMP-induced nuclear re-coils is a crucial asset in this field, as it enables unmis-takable modulation signatures for dark matter. The Recoil Directionality (ReD) experiment was designed to probe for such directional sensitivity in argon dual-phase time projection chambers (TPC), that are widely considered for current and future direct dark matter searches. The TPC of ReD was irradiated with neutrons at the INFN Laboratori Nazion-ali del Sud. Data were taken with nuclear recoils of known directions and kinetic energy of 72 keV, which is within the range of interest for WIMP-induced signals in argon. The direction-dependent liquid argon charge recombination model by Cataudella et al. was adopted and a likelihood statistical analysis was performed, which gave no indications of significant dependence of the detector response to the recoil direction. The aspect ratio R of the initial ionization cloud is R < 1 . 072 with 90 % confidence level.


Introduction
A range of evidences from astronomy and cosmology [1,2,3,4,5] indicates that a substantial fraction of the Universe is made of non-baryonic dark matter, whose nature is still unknown.Weakly Interacting Massive Particles (WIMPs), a common candidate, are actively searched for by many experiments worldwide using different technologies [6,7,8,9].In direct dark matter experiments the expected signal is the a e-mail: ds-ed@lngs.infn.itnuclear recoil (NR) induced by the WIMP elastic scattering.Because of the motion of the Solar system relative to the galactic dark matter halo, an apparent WIMP flux through terrestrial detectors is expected coming from the direction opposite to the Earth's velocity vector, i.e. approximately from the Cygnus constellation.The measurement of the NR angular distribution in a terrestrial detector could hence provide an unmistakeble "smoking gun" for WIMP-induced signals, thus making directional sensitivity a crucial asset for future direct dark matter search experiments [10,11,12].A number of R&D programs is currently in progress for directional direct dark matter search [13,14,15,16,17,18].
One of the most promising approaches for the direct search of WIMPs is based on the argon dual-phase Time Projection Chamber (TPC) [19,20,21,22], whose working principle is briefly described in the following.The TPC contains a volume of liquid argon with a thin layer of gaseous argon, the gas pocket, on the top.The elastic scattering of a WIMP with a Ar nucleus in the TPC would originate a NR of kinetic energy of a few tens of keV, which ionizes the medium along its trajectory.A prompt scintillation light signal (S1 signal) is produced by electron-ion recombination.The residual unrecombined ionization electrons are drifted towards the liquid-gas interface by an appropriate electric field, the drift field E d .They are extracted to the gas phase and accelerated by intense fields, the extraction field E ex and the electrolumiscence field E el , respectively, and emit light by electroluminescence [23], producing the S2 signal.The S1 and S2 signals are separated by the time interval corresponding to the electron drift time from the interaction site to the gas phase.The S2 signal intensity is proportional to the number of extracted electrons.
A dual-phase TPC could potentially offer a directional sensitivity for the events featuring long straight ionization tracks, thanks to the mechanism of columnar recombination [24,25,26].When the track is nearly parallel to E d , drifting electrons pass through the electron-ion column from the track itself and have a higher probability to meet an Ar ion and recombine, compared to a perpendicular track [27,28].Events with tracks parallel to E d are therefore expected to have an enhanced S1 and a reduced S2.The SCENE Collaboration has provided a hint of directional sensitivity in the S1 signal for NRs of about 60 keV [28], and specifically a difference of about 7% on S1 for NRs parallel and perpendicular to the drift field E d =193 V/cm.
The potential directional sensitivity of argon TPCs for future direct dark matter searches motivated the Recoil Directionality (ReD) experiment [29], as a part of the program of the the Global Argon Dark Matter Collaboration (GADMC).To this aim, a miniaturized argon dual-phase TPC was irradiated with neutrons at INFN Laboratori Nazionali del Sud (Italy), to produce NRs at a variety of angles with respect to the TPC drift field.The kinetic energy of NRs was around 70 keV, which falls in the range of interest of WIMP search.This work is organized as follows: Sect. 2 discusses the models to describe the response of an argon dual-phase TPC to NRs of energy relevant for dark matter searches, including the potential directional dependence.The experimental layout of ReD and the description of the individual detectors are given in Sect. 3 and 4, respectively.The data treatment and the subsequent statistical analysis to look for the directional sensitivity are presented in Sect. 5 and 6.The results are discussed in detail in Sect.7, followed by the summary of conclusions in Sect.8.

The response of Ar to nuclear recoils
WIMPs deposit energy in liquid argon (LAr) through elastic scattering on Ar nuclei.The subsequent energy loss of the NR involves nuclear stopping, ionization, charge recombination, and scintillation.Through the series of physical processes, the total energy deposited in the TPC is eventually divided into the detectable photons (S1) and electrons (S2), and the undetectable phonons (heat).
Directional modulation of charge recombination is expected when the spatial charge distribution of ionization is anisotropic.According to Refs.[27,28], this can occur when the ionizing track is longer than the Onsager radius r O , the distance between an ion and a free electron for which the electrostatic potential energy equals the thermal kinetic energy of the electron.As r O = e 2 /(6πε 0 ε r k B T ) is about 80 nm in LAr, argon ions with kinetic energy above ∼ 40 keV, i.e. well within the region of interest for WIMP searches, have a range longer than r O .However, calculations and simulations [30,31] show that the mean thermalization distance of electrons in LAr is about 2.6 µm, which is much longer than the Onsager radius.As recombination mostly takes place when electrons are fully thermalized, the directional sensi-tivity could hence be diluted by electron diffusion during thermalization.
Conventional NR charge recombination models, as the commonly-used Thomas-Imel model [32,33], often assume an isotropic charge distribution.In order to introduce the directionality, the electron distribution after thermalization needs to be included in the model.One approach is to use the Jaffé model [24,25], commonly referred to as the columnar recombination model, which is appropriate for the straight tracks from minimum ionizing particles.Since NR tracks are more localized, a more general and flexible parameterization of the charge distribution q 0 (⃗ r) has been proposed by Cataudella et al. [26], which consists of a three dimensional Gaussian with an elliptical profile where Q 0 is the total charge, σ characterizes the size of the distribution, r0 is the direction of the long axis, and R is the aspect ratio between the long and short axes.The probability of charge surviving recombination is calculated in Ref. [26] as being Li 2 the second order polylogarithm function and which depends on the Langevin recombination coefficient α [34,35] and on the electron mobility µ − .The term f (R, θ ) captures the directionality dependence and it has the functional form being θ the angle between r0 and E d .When R = 1, f (R, θ ) = 1, so directionality vanishes and Eq. 2 reduces to the Thomas-Imel model.Since directionality effects do not occur before recombination, well-established models are used here to describe the S1 and S2 yields, that for NRs also depend on nuclear and electronic quenching.The nuclear and electronic quenching factors, f n and f l , are calculated following the Lindhard [36,37] and Mei [38] models, respectively.The expectation ⟨N 0 ⟩ of the total number of quanta (ionization and excitation) produced by a NR of energy E r in LAr is where W ph = 19.5 eV is the average energy required to produce one scintillation photon in LAr [39].The detectable electron and photon yields eventually generated after recombination are respectively, where N ex /N i is the excitation-to-ionization ratio directly induced by the fast ion and by its secondaries.As a first approximation, N ex /N i is usually treated as an energy independent constant [39,40].However, the distribution of momentum transfer to electrons in the electronic stopping power is energy-dependent, which motivates the introduction of a variable N ex /N i vs. energy.This is corroborated by the SCENE data, which indicate an increase in N ex /N i with respect to the NR energy.The N ex /N i values adopted for this work are taken from Table VIII of Ref. [28], with a linear interpolation between the energy points.The N ex /N i at zero energy is set to the commonly-adopted value of 0.2.
The capability to measure the NR direction can be hidden by random fluctuations in S1 and S2, either intrinsic of the signal generation in LAr or due to detector-related effects.The fluctuation in the total number of quanta N 0 is assumed here to be Gaussian distributed with a Fano factor F = 0.107 [41], namely N 0 ∼ Gaussian(⟨N 0 ⟩, F⟨N 0 ⟩).The partition of N 0 between electrons (N e − ) and scintillation (N ph ) then follows a binomial distribution governed by N ex /N i and by the recombination probability (see Eq. 6).
The TPC signals S1 and S2 are measured in units of photo-electrons (PE) in the photosensor.The stochastic processes of collection of the scintillation light can be described by a binomial distribution, using the gain For S2, the electroluminescence process is described by a Poisson distribution depending on the amplification factor The detector response also includes a position-dependent non-uniformity which could in principle be corrected in analysis.Practically, a small residual error will be present, which can be modeled by an additional Gaussian smearing of standard deviation σ * S1 and σ * S2 for S1 and S2, respectively.Approximating the S1 and S2 distributions with Gaussians, the total contribution from detector response is In conclusion, the argon dual-phase TPC response to a mono-energetic NR follows the probability density function coming from the convolution of the detector and physical terms: Later in Sect.6, a likelihood function is evaluated from the TPC data using this probability density function.An unbinned profile likelihood study is then performed to determine the directionality parameter R.

Experimental setup
The experimental layout is conceived in order to produce and detect Ar nuclear recoils of known energy and direction, by neutron elastic scattering.Neutrons are produced by the primary reaction p( 7 Li, 7 Be)n, by shooting a 7 Li beam on a polyethylene (CH 2 ) target.The neutron energy E n and its direction are kinematically determined by measuring the energy and direction of the accompanying 7 Be nucleus.The neutron can undergo elastic scattering (n,n') with an Ar nucleus inside the TPC, thus producing a NR and a secondary neutron whose energies and momenta are again correlated by two-body kinematics.The scattered neutron is eventually detected by a neutron spectrometer made by an array of liquid scintillator (LSci) detectors; the detection of the neutron by a specific LSci determines the energy and the direction of the Ar recoil.The conceptual layout of ReD is sketched in Fig. 1.The experiment deploys three detector systems: (1) a ∆ E/E telescope made by Si detectors, to identify 7 Be nuclei associated with neutrons; (2) the TPC to detect the Ar NRs; (3) a neutron spectrometer made by seven LSci detectors to detect the neutrons scattered off Ar.The detectors of the neutron spectrometer are placed along the base circumference of a cone with axis corresponding to the target-TPC line (i.e. the direction of the incoming neutron), vertex on the TPC center and opening angle θ lsci .Therefore, all LScis detect neutrons which undergo elastic scattering on Ar at the same angle and hence produce NRs of the same energy E r .While the NRs tagged by the seven individual LScis all have the same energy E r , their momenta ⃗ p r form a different angle θ r with respect to the TPC electric field (z axis in Fig. 1), as required to test the directional effect.As it is important for this work to test the response to NRs also at θ r = 180°, the TPC is placed at a different level with respect to the target, such to Upper panel is the view of the p( 7 Li, 7 Be)n reaction plane, lower panel is the side view.The primary 7 Li beam travels along the x axis and enters the vacuum scattering chamber which hosts the CH 2 target and the ∆ E/E telescope.Neutrons emitted by the p( 7 Li, 7 Be)n reaction undergo elastic scattering inside the TPC and are eventually detected by one of the LScis of the neutron spectrometer, that are deployed within a cone of opening θ lsci = 36.8°withrespect to the target-TPC axis.See text for more details.
provide the incoming neutron with a momentum component along the field direction.
Once the angle θ t pc between the primary 7 Li beam direction and the target-TPC direction and the angle θ lsci are fixed by the setup geometry, ReD is tuned to select monoenergetic Ar recoils of energy E r by the triple coincidence between Si telescope, TPC and neutron spectrometer.The operational parameters chosen for ReD are θ t pc = 22.3°and θ lsci = 36.8°.The target-TPC distance and the TPC-LSci distance are 150 and 100 cm respectively, as a reasonable compromise between angular resolution and solid angle coverage: in both cases the uncertainty on the neutron direction is driven by the dimensions of the TPC and of the LSci, i.e. by the uncertainty on the interaction point within them.Keeping the geometry fixed, the energy E r of the NR can be changed by varying the primary beam energy.The ReD experimental layout was designed to allow for the measurements of NRs in the range of interest for dark matter direct searches, between 20 and 100 keV: this can be achieved by varying the energy of the primary 7 Li beam between 20 and 34 MeV.
3.1 7 Li beam and target The primary 7 Li beam is produced by the 15 MV TANDEM accelerator of the INFN Laboratori Nazionali del Sud [42] at an energy of 28 MeV.The TANDEM offers an excellent resolution in the delivered energy, which is about 1% FWHM in our case.The data reported in this work were collected between January 31 st and February 14 th , 2020.The current of the 7 Li beam ranged between 5 and 15 nA, corresponding to 1 − 3 • 10 10 ( 7 Li/s).The beam is driven to a vacuum scattering chamber, which hosts the CH 2 target and the ∆ E/E telescope.Upstream the target, the 7 Li beam is collimated to obtain a spot of 2 mm diameter at the target position.Neutrons are produced via the p( 7 Li, 7 Be)n reaction.The ∆ E/E telescope detects the 7 Be accompanying the neutrons that travel towards the TPC.As the accelerator does not allow the production of a pulsed beam, the direct detection of 7 Be represents the best solution for event-by-event neutron tagging.The requirement to detect 7 Be drives the choice of inverse kinematics (i.e. 7 Li beam on a hydrogenous target) [43,44], instead of the direct kinematics approach (proton beam on a 7 Li target) employed by other experiments, as SCENE.
The targets of CH 2 have thickness ranging between 150 and 350 µg/cm 2 , which is thin enough to allow for the escape of 7 Be.Due to aging effects, each target was used for about 12 hours of data taking, before being replaced by means of a 12-target holder system placed inside the vacuum scattering chamber.
After the target, the 7 Li beam travels straight forward towards a Ta beam dump placed 3 m downstream (see Fig. 1).Such a long distance is functional to minimize the background on the detectors due to the beam interaction on the beam dump.The beam intensity was precisely measured every few hours of operation by a Faraday Cup deployed about 30 cm downstream the target.However, the Faraday Cup was removed during the data taking, in order to reduce the background radiation close to the TPC.The continuous monitoring of the beam intensity was performed by measuring the rate of the 7 Li elastic scattering on a dedicated Si detector (not shown in Fig. 1) placed at θ = 7°with respect to the beam line, where no 7 Be is allowed by kinematics.

The ∆ E/E telescope
Neutrons directed towards the TPC are produced in association with 7 Be nuclei of energy E Be = 19.0MeV and emitted at angle θ Be = 5.1°. 7Be is detected by a dedicated ∆ E/E telescope placed in the scattering chamber at a distance of 46 cm from the CH 2 target.The telescope is made of two Si detectors manufactured by ORTEC, having thickness of 20 µm and 1000 µm, respectively; the 7 Be loses about 7.6 MeV crossing the thinner stage and it is stopped in the thicker one.The detectors have a 100% efficiency for light charged particles detection and energy resolution of about 1%.The telescope is collimated using an Al shield with a hole of 2 mm diameter.For the fine tuning of the position, the telescope holder is mounted on a two axis remotely-controlled stepper motor which can operate in vacuum.The detectors are readout from a standard spectroscopic chain made by a preamplifier and a charge-sensitive amplifier, with 1 µs shaping time.
The combined measurement of ∆ E and E provides the discrimination in Z, which is necessary to distinguish the interesting Be from the far more abundant elastically-scattered Li.Fig. 2 shows the ∆ E vs. E scatter plot, upon the irradiation of the CH 2 target with the 7 Li beam.The central, and most intense, band is created by Li (Z = 3), mostly by elastic scattering on H and C. The uppermost band is due to Be (Z = 4).As the reaction p( 7 Li, 7 Be)n occurs in inverse kinematics, two different solutions at the same angle θ Be = 5.1°are allowed, with 7 Be having energy of 19.0 MeV ("low energy") and 20.4 MeV ("high energy"), respectively.Neutrons in association with the "low energy" 7 Be are those travelling towards the TPC (θ n = 22.3°), with E n = 7.3 MeV kinetic energy.The "high energy" 7 Be is associated with neutrons of E n = 2.7 MeV emitted at θ n = 44°: these neutrons do not hit directly the TPC, but can contribute to accidental coincidences due to scattering on the floor or on the walls.In Fig. 2 the loci from the two 7 Be solutions are visible and clearly separated; the population between them is due to the inelastic interaction p( 7 Li, 7 Be*)n', which also emits a neutron.Because of the finite extension of the beam spot and of the beam angular divergence, neutrons associated with the 7 Be* detected at θ Be can still travel inside the Fig. 3 ∆ E vs. E distribution obtained from the irradiation of a CH 2 target with a 28-MeV 7 Li beam (color scale).The black dots are the events detected by the Si telescope in coincidence (within 200 ns) with an S1 signal in the TPC having a time profile compatible with a neutroninduced interaction.The dashed red box represents the 7 Be selection cut used in the following analysis and described in Sect.5.2.Inset: distribution of the time difference ∆t between TPC and Si telescope for events within the 200 ns coincidence gate.
TPC and produce an interaction; they also contribute to the diffuse background, e.g.upon scattering on the walls or on the floor of the experimental area.
In order to suppress the dominant contribution from 7 Li elastic scattering, the thresholds for the ∆ E and E detectors shown in Fig. 2 as dashed lines, were used during the data acquisition.Fig. 3 displays the ∆ E vs. E scatter plot, acquired with the thresholds of Fig. 2, without (color) and with (dots) the requirement of coincidence with an event in the TPC compatible with a neutron interaction and within a 200 ns gate.As expected, neutron events in the TPC are mostly associated with a "low-energy" 7 Be nucleus detected by the Si telescope.The dashed red box represents the 7 Be selection cut used in the following analysis and described in Sect.5.2.

The Time Projection Chamber
The heart of the ReD system is the dual-phase Ar TPC, whose detailed description and performance are reported in [29].It is a cubic volume of 5 × 5 × 6 cm 3 , delimited on the top and bottom by two transparent windows that are operated as anode and cathode, respectively.The lower part of the TPC contains LAr: the liquid fills the entire volume between the cathode and the extraction grid, plus 3 mm above the grid, while the gas pocket occupies the remaining 7-mm thick region up to the anode.
The TPC electric fields which are set for this work are: drift field (E d ) of 152 V/cm; extraction field (E ex ) of 3.9 kV/cm; and electroluminescence field (E el ) of 5.9 kV/cm.The maximum drift time is about 66 µs: this is the time required for an electron produced at the cathode to travel until the liquid surface.The extraction field is strong enough to give a 100% extraction efficiency of the electrons from the liquid to the gas phase [45].
The scintillation and electroluminescence signals are detected by two 5 × 5 cm 2 tiles of Silicon PhotoMultipliers (SiPMs), each containing 24 cryogenic SiPMs [46] arranged in a 4 × 6 array.The tiles are placed behind the top and bottom windows of the TPC.As the position of the S2 event in the gas phase can be used to estimate the x − y coordinate of the original interaction point in the TPC, the SiPMs of top tile are readout in 22 channels for improved resolution: 20 SiPMs are readout individually, while 4 lateral SiPMs are summed in pairs and grouped into two readout channels.The SiPMs of the bottom tiles are summed in groups of twelve, hence giving two readout channels.The SiPMs are operated at +7 V of overvoltage with respect to the breakdown voltage.Due to the presence of resistors in the bias chain, the effective overvoltage of the SiPMs gets smaller than the nominal +7 V when the bias current of the devices is high.This typically happens when the SiPMs are exposed to a significant amount of light, e.g.due to the high interaction rate under beam irradiation, and causes a change in the SiPM response (see Sect. 5.1).
More details about the cryogenic setup, the TPC, the photosensors and the readout system can be found in [29].

The neutron spectrometer
The neutron spectrometer used in ReD is made of seven 3inch liquid scintillator (LSci) cells, individually read-out by photomultipliers (PMTs).The assembly includes the liquid scintillator cell, a ETL-9821B PMT and the front-end electronics with the amplifier.The cells are filled with the EJ-309 liquid scintillator by Eljen Technologies, which features a very powerful neutron-γ pulse shape discrimination (PSD) based on the time pattern of the scintillation pulse.The neutron detection efficiency of the detectors was measured individually by using a 252 Cf source [47,48] and found to be about 28% for the 7-MeV neutrons of interest for this work.The calibration of the energy scale was performed with γ-ray sources ( 241 Am, 137 Cs and 22 Na).Dedicated measurements taken with the annihilation γ-rays from the 22 Na source confirmed the time resolution to be better than 1 ns (rms).

Data acquisition and control infrastructure
The output signals from all of the detectors are sent to CAEN V1730 Flash ADC Waveform Digitizers and digitized with 14-bit resolution at a sampling rate of 500 MHz.In total a signal of 100 µs (50k samples) is acquired at each trigger: this is sufficiently long to contain the S1 and S2 signals of the TPC, given the maximum drift time of 66 µs.Two 16channel CAEN V1730 boards were used for the measurement, synchronized with a daisy chain.
The data acquisition (DAQ) software was built upon a package developed for the PADME experiment [49] and based on the CAEN Digitizer Libraries.The trigger logic is implemented by means of an external NIM logic module as the coincidence of the Si telescope with any other detector of the setup (i.e.either the TPC or one of the LScis of the neutron spectrometer).This logic allows to collect a large sample of neutron events in the TPC and to maximize the trigger efficiency for the rare triple-coincidence signal events.The Si telescope trigger is built as the coincidence of the ∆ E and E detectors, with the thresholds displayed in Fig. 2. The TPC trigger consists in the logical AND between the two readout channels of the bottom tile within a coincidence gate of 200 ns, in order to suppress the dark rate [29].The individual thresholds are set to approximately 2 PE.The TPC is expected to trigger with 100% efficiency on S1 signals from the E r = 72 keV NR events (S1 ∼ 190 PE) which are of interest for this work, although trigger inefficiencies can possibly come from pile-up.The neutron spectrometer trigger is produced by the logical OR of the five readout channels of the seven scintillators.The energy threshold of each cell is set to approximately 20 keV ee (electron equivalent), which corresponds to about 200 keV for a proton recoil [47].This is sufficient to have a nearly-100% trigger efficiency for the neutron events of interest, as their elastic scattering on the scintillator produces protons of average energy ∼ 3.6 MeV, giving a 1.1 MeV ee signal.
All detectors and sensors of the setup can be operated and readout remotely by means of a slow control system made of a suite of LabVIEW-based [50] applications.All parameters under control (e.g.temperatures, bias voltages, leakage currents) are monitored continuously, and readings are stored in a database every 10 s.

Event reconstruction and calibrations
The raw data from the TPC are the digitized waveforms of each of the SiPM channels, from which the event type, time, and 3D position were reconstructed following the procedure described in [29].A dedicated pulse-finder algorithm searched for possible S1 and S2 signals.Each pulse was classified as either S1 or S2 by using the pulse shape parameter f p , defined as the ratio of the charge in the first 700 ns over the total charge: pulses with f p < 0.2 are classified as S2.The pulse-finder algorithm is fully efficient for S1 signals above a few keV.The time delay between the S1 and S2 pulses, i.e. the electron drift time t drift , was used to estimate the z coordinate of the interaction below the liquid-gas interface.Events with a single S1 pulse and a S2 pulse with t drift between 6 µs and the maximum drift time were kept for the subsequent analysis.The cut t drift > 6 µs removes the events produced just below the extraction grid of the TPC, in which the S1 and S2 pulses are piled-up.The approximate x − y position of the event was evaluated as the charge-weighted center of the S2 signal in the top SiPM array.The parameter f p defined above was also used to perform the NR/ER discrimination: S1 pulses with f p > 0.4 are selected as from NR.This simple cut was shown to allow for a NR/ER separation better than 2 σ for S1 above 50 PE [29].
The Single Electron Response (SER) and the cross-talk and afterpulsing effects were studied by irradiating the SiPMs with a 403-nm laser source and by modeling the photon counting statistics according to the geometric chain process model by Vinogradov [51,52].The calibration was performed channel by channel, as described in [29].The final PE gain is corrected to remove cross-talk and afterpulsing.Dedicated laser calibrations were taken every 12 hours throughout the beam time to monitor the stability of the SiPMs.
As mentioned in Sect.4.2, the voltage drop in the bias resistor chain causes a reduction of the bias voltage of the SiPMs, which is proportional to the bias current and must be properly accounted for in the data analysis.The bias current registered during the laser calibrations by the slow control system was < 0.5 µA.During the beam irradiation, because of the much higher interaction rate and the much higher amount of light hitting the SiPMs, the bias current ranged up to 90 (150) µA for the bottom (top) SiPMs, depending on the intensity of the primary 7 Li beam, which was not constant in time.To derive the corrections to the SER for each individual SiPM, three dedicated laser runs in which the TPC was simultaneously irradiated with high-activity radioactive sources were performed.The typical correction is of the order of 0.5% • I, where I is the bias current in µA.For this reason, the SER and the cross-talk and afterpulsing corrections were time-dependent and calculated using the closest reading of the bias current registered by the slow control.Besides the SER, the photon detection efficiency also changes with bias voltage: a set of runs with 241 Am 60 keV γ and the 7 Li beam irradiation was performed to calibrate the additional bias current dependency in PE yield.
Additional calibrations with 241 Am were taken daily during the campaign, to evaluate the dependence of the TPC response on the interaction position, and to determine the correction factors for S1 and S2.The events featuring one sin- gle S1 and one single S2 and having S1 compatible with the full energy deposition of the 60 keV γ-ray from 241 Am were grouped in a 22 × 11 mesh, according to the interaction position in the TPC.The mesh has 22 entries in x − y, based on the top SiPM channel detecting the largest fraction of the S2 signal, and 11 bins in z, equally spaced between t drift = 6 µs and 72 µs.Firstly, S2 was corrected to account for the presence of impurities in LAr, which can cause the absorption of electrons during their drift path.The electron life time was typically > 1 ms, i.e. much longer than the 66-µs maximum drift time, and it was estimated with an exponential fit of the S2 vs. t drift profile, restricted to the events in the central eight x − y bins.The z dependency of S1 and S2 was further corrected by using a set of 5th-order polynomials S1 i (t drift ) and S2 i (t drift ): they are calculated by interpolating over the z-points within each bin i in x−y.Three examples are shown in Fig. 4: the correction vs. z is within 10-15%, for both S1 and S2.No significant variation in the position correction was found throughout the sequence calibration runs.Position dependencies mostly result from non-uniformity in the light collection efficiency within the TPC: as a consequence, the same corrections for S1 and S2 derived from 241 Am (ER) events were also applied to NR events.The calibrations with 241 Am were also used to evaluate the light yield of the TPC at 60 keV and at E d =150 V: it is (8.53±0.19)PE/keV, which is very well consistent with the expectation of 8.6 PE/keV based on the parametrization obtained in the pre-irradiation campaign [29].
A simpler processing was performed for the digitized waveforms from the liquid scintillators and from the Si detectors of the telescope.The signal in the LSci detectors was processed by calculating the total charge, integrated within a gate of 600 ns.The ratio between the charge in the first 80 ns and the total was used as the discrimination parameter, resulting in a neutron-γ discrimination better than 3σ above 200 keV ee [48].The signals from the E and ∆ E detectors of the telescope were evaluated by taking the maximum of the digitized shaped waveforms from the charge-sensitive amplifier.
The time signal of all three kinds of detectors in the setup is critical for the coincidence event selection.The time stamp of a TPC event was defined as the zero-crossing time of the pulse obtained by passing the S1 pulse through a digital constant fraction discriminator (CFD).The ∆ E/E telescope generates two time stamps, one for the ∆ E detector and one for the E detector, which were both evaluated with CFDs.The reference time for the ∆ E/E telescope used for the coincidence was taken as the average of the two time stamps.Finally, the time stamp for the neutron spectrometer was defined as the zero-crossing CFD time of the digitized waveforms.

Selection of signal events
The events of interest are triple coincidences between a 7 Be nucleus detected in the ∆ E/E telescope, and the two subsequent neutron scatterings in the TPC and in the neutron spectrometer.
A clean sample of signal events with the proper topology was selected through a sequence of cuts.Firstly, unambiguous TPC events were selected according to same criteria of Sect.5.1: events with only one S1 and only one S2, separated by a t drift within the range [6,66] µs.An additional S2 "echo" signal, namely a secondary event due to photoionization of the cathode from the main S2 electroluminescence, is allowed in the time window [67.5, 72] µs after the primary S2.
Afterwards, events in the TPC were selected by requesting that S1 is in time coincindence within a gate of 200 ns with the ∆ E/E telescope and with one single LSci detector of the neutron spectrometer.In addition, neutron-induced (n,n') events in the neutron spectrometer were efficiently selected by PSD against the dominant γ-ray background.The ER/NR discrimination based on the S1 signal of the TPC was not applied.This was meant to avoid an undesirable S1dependent selection efficiency, given the fact that the discrimination based on f p gets progressively worse for S1 signals below 100 PE.
The 7 Be ion which accompanies the neutron traveling towards the TPC was selected by a combined cut on ∆ E and E, which is shown in Fig. 3 (red dashed contour).The selection is not sensitive enough to resolve between the 7 Be emitted at the ground state in the p( 7 Li, 7 Be)n reaction and 7 Be* in the first excited state coming from the p( 7 Li, 7 Be*)n' reaction.Therefore, the neutron energy distribution consisted of two different mono-energetic components.
The data sample was further selected by using the timeof-flight (ToF) of the TPC with respect to the ∆ E/E telescope, namely by keeping the events in which the delay between the telescope and the TPC (see inset in Fig. 3) is consistent with the flight time of the neutrons.The coincidence window in ToF was set to be S1-dependent, in order to ensure a S1-independent selection efficiency.The boundaries of the coincidence window were defined as the 1% and 99% quantiles in each S1 slice of 10 PE, after the subtraction of the constant background due to random neutrons and γ-rays.The random background contributes to about 1 % of the events in the coincidence windows.
The coincidence windows for the delay ∆t(LSci−SiTel) between the LSci and the telescope in triple-coincidence events were set with very stringent cuts, so to guarantee the selection of pure single-scattering neutron interactions.The timing of the individual LScis was calibrated by using as a reference the γ-rays produced in the TPC by inelastic interactions (n,n'γ) and then detected in the LScis: all γ peaks were aligned to ∆t(LSci − SiTel) = 0, as displayed in Fig. 5, where the effect of used cuts applied sequentially is shown.The single-scattered neutron events of interest form the peak around 20 ns.The low-statistics peak at about 25 ns comes from the lower-energy neutrons produced in the p( 7 Li, 7 Be*)n' interactions, while the tails at longer times are mostly due to multi-scattered neutron background.Monte Carlo simulations indicate that the hump around 60 ns is originated by the neutrons associated with the "high energy" 7 Be, which reach the TPC after scattering on the floor or other passive structures.The peaks around −35 ns and −20 ns are γ-rays emitted by p( 7 Li, 7 Li*)p inelastic scattering.Gaussian fits to the peak around 20 ns determined the position and width of the window, individually for each scintillator.As mentioned in Sect.4.3, the LSci channels which selected NR events at θ r = 20 • and 40 • were each made from the analogue sum of the signals of two different detectors.Since the cable lengths for the two detectors at 20 • were not properly matched, this introduced a split in the timing: the ∆t(LSci − SiTel) distribution for the channel at 20 • was hence fitted with a double Gaussian.The coincidence windows were defined according to the position µ and width σ of the peaks from the Gaussian fits, as summarized in Table 1 and they are used to select the triple coincidence events.The coincidence windows were further extended by 5 ns in order to include the slower neutrons from p( 7 Li, 7 Be*)n'.Side-bands were also defined to estimate the random coincidence rate in each channel, see Table 1.The triple coincidence events eventually considered for the statistical analysis of Sect.6 are those which pass the sequence of cuts displayed in Fig. 5 and the additional selection in the ∆t(LSci − SiTel) ToF from Table 1.

Statistical analysis
The S2 vs. S1 distribution of the NR events in the TPC which pass the selection procedure of Sect.5.2 is displayed in Figure 6: the pink dots represent the events selected requiring the triple coincidence (TPC, Si telescope and neutron spectrometer); the colour-coded distribution includes the events in double coincidence (TPC and telescope).The triple coincidence sample contains about 650 NR events with S1 above 20 PE, which were collected during 10.7 live days of beam run.The double coincidence events constitute a large sample of about 70000 TPC NR events in all directions: they were hence used as a calibration data set to constrain the nuisance model parameters in the global fit below.
The data samples were statistically analysed in order to evaluate the best estimate of the directionality parameter δ R = R − 1, which measures how much the shape of the initial ionization charge cloud differs from a sphere.As the number of events is relatively modest, an unbinned profile likelihood was applied.Fig. 6 S2 vs. S1 distribution of NR events in the TPC.The color-coded histogram includes the event in double coincidence (TPC and telescope), namely at all angles θ r with respect to the electric field.The pink circles are the events in triple coincidence (TPC, telescope and spectrometer).All corrections and cuts are applied.The white contour is the fit range in the (S1,S2) plane used for the statistical analysis described in Sect.6.

Likelihood function and fit parameters
The global likelihood L is written as a product of three likelihood terms: where the product over i refers to the five samples taken at the five angles θ 1, each containing the observed array of events X i = (S1, S2); δ R is the parameter of interest (POI); ν is the array of nuisance parameters; X cali is the array of calibration data set.The POI is constrained in this work to δ R ≥ 0, as negative values of δ R are not physically allowed by the recombination model [26].The three likelihood terms of Eq. 11 are described in detail below.
L i is the extended likelihood of each sample of NR events at the recoil angle θ (i) r : where n i and ni are the size of X i and its mean, respectively, and P i is the joint probability density function (PDF) of the events (S1,S2).The PDF is made as the combination of three The first component is the energy spectrum for the signal F sig (E r ), which depends on the recoil energy E r , convolved with the response function P of the TPC to mono-energetic NR events, as defined in Eq. 10.The parameters λ 1i are the fractions of random coincidences within each data sample: they were estimated from the data, using the counting rate in the side-band in ToF and are listed in Table 2. Similarly, λ 2 is the scaling factor for multi-scattering background, namely the fraction of those events with respect to all NR events in the coincidence window.The other two components are the energy distributions of the backgrounds due to random coincidences, F bkg1 (E r ), and to multiple neutron scattering, F bkg2 (E r ).They are also convolved with the response function P of the TPC.As the angular distribution for background events is approximately random, the θ r dependence of f (θ r , R) is averaged out by using the equivalent angle θr calculated analytically for an isotropic distribution and the functional dependence on the angle is approximated as The factor λ 2 and the three energy spectra (F sig , F bkg1 , and F bkg2 ) were evaluated by means of a dedicated Monte Carlo simulation using the Geant4-based framework g4ds [53,54,55,56].The events from the simulations underwent the same sequence of selection cuts used for the real data.The energy distributions derived by the Monte Carlo are displayed in Fig. 7.The three energy distributions were then analytically parametrized in order to optimize the calculation of the CPU-intensive PDF P i .F sig consists of two Gaussian peaks corresponding to the NR induced by neutrons from p( 7 Li, 7 Be)n and p( 7 Li, 7 Be*)n'.F bkg1 and F bkg2 were approximated by a double-exponential and a single exponential, respectively, whose parameters were calculated by fits to the Monte Carlo distributions.
The factor L cali of the global likelihood of Eq. 11 is the constraint term on the nuisance parameters and it depends on the events X cali in the calibration set (i.e.colour-coded histogram in Fig. 6).While the energy spectrum of the calibration events is a broad and featureless distribution, the joint  7 Signal and background spectra from g4ds Monte Carlo simulation.Solid and dotted black histograms show the distribution of the NR energy E r for signal events associated with 7 Be and 7 Be * neutrons, respectively.Solid red line shows the random coincidence background spectrum F bkg1 and the dashed red line the spectrum F bkg2 from multiscattered triple coincidence events.The error bars represent the Monte Carlo statistical uncertainties.distribution of the NR band in the (S1,S2) plane can set a strong constraint on the nuisance parameters.Since the fraction of signal events in the calibration sample is negligible, the energy distribution is well approximated by the random background F bkg1 .The calibration term is hence written as: In order to avoid any analysis bias, δ R should be decoupled from the nuisance parameters as much as possible.The explicit occurrence of the POI δ R in Eq. 14 is due to the fact that the parameter ξ m in the modified Thomas-Imel model in Eq. 2 is dependent on δ R because of the track length.To remove such undesirable degeneracy, the angular dependence term and the Thomas-Imel parameter of Eq. 2 were re-defined as and respectively.In this way the angle-averaged position of the NR band in calibration data does not depend on δ R and the POI δ R is left as a pure representation of directionality.Furthermore, the degenerate nuisance parameters were re-cast into a unique nuisance parameter which represents the recombination probability of one electron-ion pair.
The last factor of the global likelihood, L contraint (ν), is the pull term for the nuisance parameters which were known by prior independent measurements.Those parameters are constrained by Gaussian terms Table 3 List of the parameters used in the model.δ R is the parameter of interest, while all others are nuisance parameters, constrained by the calibration data and/or by a Gaussian pull term.The error bars are the standard deviation which is taken in the Gaussian pull terms.The parameters reported without uncertainties are fixed.The gains g 1 and g 2 come from the previous TPC performance study [29].The S1 resolution of the TPC of Eq. 10 is parametrized as σ based on the previously-measured values ν 0 i of the parameters ν i and on their corresponding uncertainties.
The parameters and their reference values are summarized in Table 3.The recombination probability A depends on σ , the size of the ionization cluster of Eq. 1, which is dominated by the electron diffusion during thermalization.Due to their high mobility and long thermalization time, electrons diffuse for a few µm in LAr [30,57].It is found that A = 0.04/e − , which corresponds to σ = 1.8 µm, was an appropriate initialization parameter for the likelihood fit.The ratio N ex /N i was treated as a function of recoil energy, as discussed in Sect. 2. The TPC gains g 1 and g 2 were estimated according to the TPC characterization in [29], and were treated as nuisance parameters in order to accommodate for possible variations in the TPC performance.The parameters W ph , k e , λ 1 and λ 2 were fixed in order to limit the degeneracies in the fit: their effect on the POI is minor and is accounted below as a systematic uncertainty.

Results
Experimental data of Fig. 6 (calibration and five triple-coincidence samples) were fitted against the model of Eq. 11.In order The fit result is shown in Fig. 8 and reported in Table 4.The positions of the signal peak in both S1 and S2 (middle and bottom rows of Fig. 8) are mutually consistent among the five samples at different θ r .The best-fit for the POI is δ R = 0.037 ± 0.027, which is less than 2σ away from a null result; the uncertainty on δ R is largely driven by statistics.The upper limit of δ R is calculated by a toy Monte Carlo approach, in order to guarantee the correct coverage: it results to be δ R < 0.072 at 90% CL.The best-fits of the nuisance parameters are in good agreement with the central values of their estimates used for the constraints.In particular, the smallness of the best-fit for the parameters σ * S1 /S1 and σ * S2 /S2, which are the extra (non-statistical) contributions to the experimental resolution in S1 and S2, demonstrates that the spatial inhomogeneities of the detector response were properly corrected.Furthermore, the proper convergence and the absence of a significant bias for all fit parameters, notably including the POI δ R, were checked by running a dedicated set of toy Monte Carlo simulations.
The uncertainties on δ R related to the nuisance parameters are automatically accounted in the fit.All other systematic uncertainties on δ R, e.g.those related to the values of W ph , k e , λ 1 and λ 2 , to the spectral shapes F sig , F bkg1 and F bkg2 , and to the approximation of θr from isotropic distribution, are globally evaluated to be an order of magnitude smaller than the statistical term and are hence neglected in this work.

Discussion
The results of this work suggest that the charge recombination in NRs in the energy range of interest for WIMP dark matter searches has a limited directional dependence.A possible explanation is that the directional effect is washed out in the isotropic thermalization process of the electrons: the range of 70 keV argon ions in LAr, 0.18 µm [58], is much shorter than the electron thermalization radius ∼2.5 µm [30,31].If all electrons were confined within the Onsager radius, the recombination probability A would be 8/e − , namely, two orders of magnitude higher than measured in this work.This indicates that the extension of the thermalized electron cloud is much bigger than the Onsager radius, thus weakening any initial directional effect.Other non-local processes at the length scale of a few µm can also contribute to the size of the electron cloud, including the emission of Auger electrons and fluorescence X-rays from excited Ar atoms [59,60].
The strongest constraint on δ R from the fit comes from the position of the S2 peak, since g 2 ≫ g 1 .In fact, the SCENE hint for directional sensitivity was primarily given by the 7% variation in S1 for NRs parallel and perpendicular to E d : no variation of S2 vs. direction was observed.While the SCENE data were never analyzed according to the directional model of [26], an asymmetry δ R ≈ 2 would be re-quired to generate a 7%-effect on S1.However, given the anti-correlations of Eqs. 6 and 7, such a large δ R would produce a much more significant variation in S2 (∼ 80% between parallel and perpendicular directions), which is not observed in SCENE.The lack of a variation in the S2 signal, which is further confirmed in this work, rules out the directional modulation in charge recombination as the explanation of the effect and sets an upper limit on δ R. Furthermore, the ReD data, with an improved signal yield and resolution in S1, do not confirm the variation in S1 at different directions which was reported by SCENE.As for S2, no statistically-significant variation was found for S1.
The LAr signal model adopted in this work has two major upgrades comparing to the models commonly used in the literature.The first modification is about charge recombination, by the introduction of the directional term of Eq. 2. The second modification is the use of an energy-dependent ratio N ex /N i , which allows for a better fit of the NR band shape and improves significantly the performance of the likelihood fit.If N ex /N i is kept constant to the value 1, which is commonlyadopted for NRs [39,40], the fit still returns a value of δ R compatible with zero, but the model fails to reproduce the shape of the S2 vs. S1 band and the S2 distribution for NRs; furthermore, the best-fit for g 2 in this case is 29.9±0.1 PE/e − , which is in tension with the prior measurement of Table 3.While the SCENE data also support the energy dependence of N ex /N i , the physical motivation of it requires further study.One possibility is that this is the apparent effect of energy dependences in the nuclear quenching, electron quenching and recombination processes, which are unaccounted by the Lindhard, Mei and Thomas-Imel models used in this work.Specifically: the Thomas-Fermi screening function used in the Lindhard model is known to have a bias in the O(10)keV range [37,61]; the Mei model simplifies the average electronic stopping power by taking the value at the initial electron kinetic energy; the charge recombination model does not consider the dependence on the charge cloud size on the recoil energy.All these energy-dependent factors are not accounted in the models and they could eventually show up as an effective energy dependence in N ex /N i .It has anyhow a small effect on the systematic uncertainty on δ R, due to the weak correlation reported in Table 4.

Conclusions
The Recoil Directionality (ReD) experiment was designed within the GADMC to explore the possible directional sensitivity of an Ar dual-phase TPC to nuclear recoils in the energy range of interest for WIMP dark matter searches.The ReD TPC was irradiated with neutrons of known energy and direction at the INFN Laboratori Nazionali del Sud, in order to produce Ar recoils of about 70 keV kinetic energy.Nuclear recoils traveling in five different directions with respect to the drift field E d of the TPC were selected using a neutron spectrometer made by liquid scintillation detectors.A statistical analysis based on the Cataudella et al. model [26] was performed to assess the TPC response for those samples of NR events.
The data from this work do not show any statisticallysignificant dependence of either S1 or S2 on the direction for NRs of ∼ 70 keV.The upper limit for the parameter of interest R, which measures the aspect ratio between the long and short axes of the initial electron cloud, is R < 1.072 at 90% CL.The absence of significant deviations from the spherical symmetry of the electron cloud indicates that the electron thermalization process likely plays a significant role in weakening any initial direction-induced anisotropy of the charge cloud.

Appendix A: A data-driven analysis approach
The directionality analysis presented in this work depends on the specific model by Cataudella et al. which was adopted to describe the phenomenon.One possibility to release such a model depencence is to employ a data-driven approach based on Machine Learning (ML) techniques.ML techniques are in fact very effective in revealing possible correlations between quantities in the study of phenomena for which large data samples are available, even if a model for their description is lacking.
Supervised learning algorithms were used to try to highlight the signature for possible directionality effects in the electron-ion recombination in the ReD data [62,63].Due to the limited size of the triple coincidence event samples, an indirect approach was adopted, which makes use of all TPC calibration events.The data set of the double coincidence events provides a two-order-of-magnitude larger amount of data for the training of the model, and this is a desirable condition when working with ML algorithms.
In an ideal LAr TPC, the S2 signal is expected to be related to S1 through some functional form S2= f (S1).The basic assumption of this strategy is that the function f does not depend on the direction of the Ar recoil, namely that the angle θ r between the recoil and the drift field does not affect the balance between S1 and S2.Deviations from this trend would highlight a possible effect of the recoil directionality.
The model was derived by using the calibration data set, which is made of NR events characterized by a wide distribution of angles θ r .The data set contained about 72000 events and it was randomly split (70:30) in a training and testing set, on which the model was trained and tested, respectively.During the training phase, the algorithm built the function f used to predict S2, event-by-event, based on the patterns which are learnt from the training set.Each pattern consists of a vector of features: S1 signal [PE], x − y position [cm], and t drift [µs] as the z coordinate of the event, within the appropriate ranges, and the measured S2 value as a target.The derived model aims to predict the value of the ionization signal S2, for each of the events, from the knowledge of S1 and of the reconstructed interaction point within the TPC.
The Extreme Gradient Boosting algorithm (xgboost) was used to derive the model [64].To evaluate the accuracy of the model, the metric of the relative prediction error was adopted.This is defined, for the i-th pattern (i.e. the i-th event in the TPC), as The trend of ε i pred was investigated for each event, and also against each feature describing the patterns, to verify that there were no regions in the feature space in which the model has a worse response that could introduce any bias in the predictions.At the end of the training phase, the model was able to provide a satisfactory prediction of the experimental S2 of the events in the testing set: the relative errors ε i pred resulted to be approximately Gaussian-distributed with mean 0.0043 (6) and standard deviation 0.09.Subsequently, the model was used to make predictions on the triple coincidence data set.For these data, the known θ r values are used to check for possible directional-dependent deviations of the predicted S2 values compared to those measured experimentally.ε i pred was initially calculated for each event in the triple coincidence data set and the corresponding values were subdivided into four subsets, according to the angle θ r determined by the coincident neutron detection.The mean value of ε pred in each data set and the corresponding uncertainty are displayed in Fig. 9 as a function of the recoil direction θ r .The point at θ r = 0 • is lower than the others, as expected in the case of directionality effects, since traces parallel to E d would result in enhanced S1 signals and reduced S2.Nevertheless, experimental data are compatible with the null hypothesis of no directionality effect: the pvalue calculated from the χ 2 test is 23%.Therefore, the datadriven analysis carried out using ML techniques on the data collected in the ReD TPC is compatible with the absence of any directional effect1 , in agreement with the analysis based on the model by Cataudella et al.

Fig. 1
Fig. 1 Schematic layout of the ReD experimental setup (not in scale).Upper panel is the view of the p( 7 Li, 7 Be)n reaction plane, lower panel is the side view.The primary7 Li beam travels along the x axis and enters the vacuum scattering chamber which hosts the CH 2 target and the ∆ E/E telescope.Neutrons emitted by the p( 7 Li, 7 Be)n reaction undergo elastic scattering inside the TPC and are eventually detected by one of the LScis of the neutron spectrometer, that are deployed within a cone of opening θ lsci = 36.8°withrespect to the target-TPC axis.See text for more details.

Fig. 2 ∆
Fig. 2 ∆ E vs. E scatter plot obtained from the irradiation of a CH 2 target with a 28-MeV 7 Li beam.The bands identify nuclei of different Z (α, Li and Be), as discussed in the text.Neutrons traveling towards the TPC are produced in association with the 7 Be nuclei of the locus labeled as " 7 Be low".The dashed lines show the thresholds used in the normal operating conditions for the ∆ E and E detectors, and which are meant to suppress the dominant contribution from Li.

Fig. 4
Fig. 4 Examples of correction factors for the z position dependence for events located below two central SiPMs (blue and green) and below one corner SiPM (red).Upper (lower) panel: correction factor for S1 (S2).

Fig. 5
Fig. 5 Timing spectra of ∆t(LSci − SiTel).All channels of the neutron spectrometer are aligned with the γ peak at t = 0. Cuts are applied accumulatively starting from the raw distribution (dark blue dashed histogram).
Fig.7Signal and background spectra from g4ds Monte Carlo simulation.Solid and dotted black histograms show the distribution of the NR energy E r for signal events associated with7 Be and7 Be * neutrons, respectively.Solid red line shows the random coincidence background spectrum F bkg1 and the dashed red line the spectrum F bkg2 from multiscattered triple coincidence events.The error bars represent the Monte Carlo statistical uncertainties.

Fig. 8
Fig. 8 Experimental data superimposed with the best-fit model.The fit is performed for S1 ∈ [120, 400]PE and within the white contour of Fig. 6.Upper row: S2 vs. S1 distribution for the calibration data set (color-coded histogram) and for the triple-coincidence data sets at different angles θ r (pink circles).Middle row: projection on S1 for the triple-coincidence samples.The dashed vertical lines at 120PE mark the left edge of the fit range.Bottom row: projection on S2 for the triple-coincidence sample in the range of S1 ∈ [120, 400]PE.The blue, red and pink curves are the total spectrum, the random coincidence background F bkg1 , and multi-scattering background F bkg2 , respectively.

Fig. 9
Fig.9Mean relative prediction error for each ε i pred distribution obtained by splitting data into the four data sets.The red dashed line marks the level ε pred = 0.0043 from the testing set (i.e.null directional effect).

Table 2
Fraction of random coincidence events, λ 1i , in the range S1 ∈ [120, 400]PE and S2 ∈ [800, 2800]PE in the five samples of triplecoincidence events at different θ r .Uncertainty is about 2% for all samples.

Table 4
Best-fit of the parameters and correlation coefficients between the nuisance parameters and the POI δ R. [120, 400]PE, as represented by the white contour in Fig.6.The S1 range corresponds to NR energies between approximately 35 and 150 keV, and hence comfortably includes the expected NR signal at ∼ 72 keV.The low-S1 edge S1 > 120 PE was set in order to avoid any inefficiencies in the event reconstruction and selection.The center of the NR band was empirically parametrized with the function S2 /[PE] = 455 ln(S1/[PE]) − 535 and the cut was set as ±500 PE in S2.The fit region globally contains 529 triple coincidence and 42340 calibration events.