Quenching factor measurements of sodium nuclear recoils in NaI:Tl determined by spectrum fitting
L. J. Bignell, I. Mahmood, F. Nuti, G. J. Lane, A. Akber, E. Barberio, T. Baroncelli, B. Coombes, W. Dix, J. T. H. Dowie, T. Eriksen, M. S. M. Gerathy, T. J. Gray, B. P. McCormick, A. J. Mitchell, M. S. Rahman, F. Scutti, N. J. Spinks, A. E. Stuchbery, H. Timmers, P. Urquijo, L. Wang, Y. Y. Zhong, M. Zurowski
PPrepared for submission to JINST
Quenching factor measurements of sodium nuclearrecoils in NaI:Tl determined by spectrum fitting
L. J. Bignell, a, I. Mahmood, b F. Nuti, b G. J. Lane, a A. Akber, a E. Barberio, b T. Baroncelli, b B. Coombes, a W. Dix, b J. T. H. Dowie, a T. Eriksen, a M. S. M. Gerathy, a T. J. Gray, a B. P. McCormick, a A. J. Mitchell, a M. S. Rahman, c F. Scutti, b N. J. Spinks, a A. E. Stuchbery, a H. Timmers, c P. Urquijo, b L. Wang, a Y. Y. Zhong, a and M. Zurowski b a Department of Nuclear Physics, The Australian National University, Canberra, ACT 2601, Australia b School of Physics, University of Melbourne, Melbourne, Victoria 3010, Australia c School of Science, University of New South Wales, Canberra, ACT 2610, Australia
E-mail: [email protected]
Abstract: We have performed measurements of sodium nuclear recoils in NaI:Tl crystals, followingscattering by neutrons produced in a Li(p,n) Be reaction. Understanding the light output fromsuch recoils, which is reduced relative to electrons of equivalent energy by the quenching factor, iscritical to interpret dark matter experiments that search for nuclear scattering interactions. We havedeveloped a spectrum-fitting methodology to extract the quenching factor from our measurements,and report quenching factors for nuclear recoil energies between 36 and 401 keV. Our results agreewith other recent quenching factor measurements that use quasi-monoenergetic neutron sources.The new method will be applied in the future to the NaI:Tl crystals used in the SABRE experiment.Keywords: Scintillators, scintillation and light emission processes (solid, gas and liquid scin-tillators), Dark Matter detectors (WIMPs, axions, etc.), Detector modelling and simulations I(interaction of radiation with matter, interaction of photons with matter, interaction of hadrons withmatter, etc), Instrumentation for heavy-ion accelerators Corresponding author. a r X i v : . [ phy s i c s . i n s - d e t ] F e b ontents Weakly Interacting Massive Particles (WIMPs) are a well-motivated dark matter candidate thatare expected to interact with ordinary matter via nuclear scattering. There is a large and ongoingsearch for such interactions from galactic WIMPs in terrestrial detectors, with most experimentsreturning null results [1]. A notable exception is the DAMA experiment, which for many years hasobserved a modulation in the count rate of low energy events in their NaI:Tl detector consistent– 1 –ith that expected from WIMP dark matter [2]. This observation has motivated numerous otherNaI:Tl-based searches that aim to independently measure the DAMA modulation [3–5], thoughnone have unambiguously confirmed or refuted this result to date.Understanding the scintillation response of NaI:Tl to nuclear recoils is essential to interpretany potential dark matter signal. The light yield from nuclear recoils is reduced relative to electronrecoils of the same energy due to energy losses to lattice excitations and additional non-radiativecarrier de-excitation pathways in regions of high excitation density [6]. Furthermore, in practice,detectors are calibrated using gamma-ray sources that produce electron recoils, so to set the correctenergy scale for nuclear recoil interactions the quenching factor Q = L NR L ER (1.1)must be known. Here L NR and L ER are the light yields due to a nuclear recoil and electron recoilof equal energy, respectively.The DAMA experiment has measured values of 0.3 and 0.09 for Na and I recoils, respectively,using a Cf source and a simulation of the expected distribution of recoils [7]. More recentexperiments using quasi-monoenergetic neutrons from deuteron-deuteron neutron generators and Li(p,n) Be reactions have indicated an energy dependence where the Na quenching factor decreaseswith energy, though there is some inconsistency between these results [8–11]. This study providesan additional measurement of the quenching factor of Na recoils in a NaI:Tl crystal, and exploits afull spectrum fitting methodology as a means to improve the modelling of the recoil spectrum.
Figure 1 illustrates the experimental setup. The choice of geometry was guided by a toy model thatused cross-section data and kinematic calculations to predict coincidence rates, neutron times offlight, and nuclear recoil energies (including their spread) based on a simplified angular acceptance.
An encapsulated 40 mm diameter, 40 mm height, cylindrical NaI:Tl crystal was employed tomeasure the nuclear recoils. The crystal enclosure includes a transparent window, allowing readoutwith photomultipliers. Two Hamamatsu H11934-200 1 x 1 inch ultra-bialkali (43% peak quantumefficiency) photomultipliers were coupled to the window using optical grease before wrapping withlayers of teflon tape and black adhesive tape. The light yield from this setup was measured, usingcalibration sources and the peak charge from single photoelectrons, to be 2 PE/keV, which waslower than expected. We attribute this low sensitivity to a combination of poor optical coupling,owing to the mismatch between the crystal window and the PMT size, and a low intrinsic light yieldfrom the crystal. A µ g/cm LiF foil on a 20 mg/cm Ta backing was placed at the target location. This wasirradiated by a pulsed beam of protons. A range of beam energies were chosen, sufficient togenerate quasi-monoenergetic neutrons via Li(p,n) Be. The target backing was sufficient to stop– 2 – igure 1 : A schematic drawing of the nuclear recoil experiment (not to scale). Neutrons wereproduced in the tantalum-backed LiF target using Li(p,n) Be. These were collimated with acylinder of HDPE (grey) onto the NaI:Tl (red) and the time-of-flight to events in the liquid scintillatordetectors (blue), along with particle identification via pulse shape discrimination, were used to tagnuclear recoils. The polar angle with respect to the beam axis and the distance from the face of eachliquid scintillator to the centre of the NaI:Tl is shown. The dashed lines indicate through-holes inthe HDPE collimator.the primary proton beam, but had a negligible effect upon the neutrons. The pulsed beam widthhad a FWHM less than 2 nanoseconds, repeated every 747 ns. The beam current varied throughoutthe measurements, and also with beam energy, in a range between approximately 0.2 and 2 nA.Neutron production was confirmed at the time of measurement using the time difference spectrumbetween a BaF scintillation detector placed close to the target and the beam pulse timing signal.This time spectrum exhibited a prompt peak from gammas emanating from the target and a delayedpeak due to neutrons (section 3.6). Two measurement campaigns were carried out, with protonbeam energies of 2.44 MeV, 3 MeV, and 5.2 MeV in the first and 3 MeV and 6 MeV in the second.Neutrons exited the vacuum chamber through a thinned aluminium vacuum flange. A cylin-drical HDPE collimator was placed after the target so that the 25 mm aperture was collinear withthe beam axis. The NaI:Tl detector was aligned with the aperture, 30 cm from the LiF target.This collimator included 25 mm diameter through-holes to allow the transport of neutrons to thebackward-angle scattered neutron detectors. – 3 – .3 Liquid Scintillator Detectors Seven liquid scintillator detectors were placed at the angles and distances from the centre of theNaI:Tl crystal outlined in figure 1 to tag scattered neutrons. The detectors were aligned using atheodolite, so the error in the mean angle is expected to be negligibly small.The scintillation detectors were made with a diameter of 38 mm and a length of 15 cm. Thelength was chosen to give a balance between a sufficiently high macroscopic scattering cross-section(> 90% for the neutrons traversing the cylindrical axis) and a small neutron time of flight over thelength of the sensitive volume. The neutron time of flight over the liquid scintillator volume variedfrom 14 ns to 5 ns at the lowest and highest beam energies, respectively. The inner surfaces of thedetectors were fabricated from 5 mm thick teflon, aside from a borosilicate glass window sealedwith a Viton o-ring. The scintillation light was read out via a 1.5 inch diameter Hamamatsu H10828photomultiplier with a super bialkali photocathode (35% peak quantum efficiency), coupled to thedetector window using optical grease. A thin aluminium housing covered the outer surfaces of theliquid scintillator detectors as light shielding, supplemented with duct tape and aluminium foil.We used EJ-309, a low-hazard liquid scintillator with excellent particle identification capabili-ties, as the detection material. The scintillator was sparged with nitrogen before filling the detectorsin a nitrogen glove box. A small gas volume was included with the liquid scintillator to allow forthermal expansion.
An XIA Pixie16 digitizer [12] was used to capture digitized waveform signals from all 9 detectors(7 liquid scintillators, 2 NaI:Tl channels) at 500 MSPS with 12 bit resolution in both runs. A beamsync signal was generated as a logic signal in phase with the beam bunches. This was captured asa constant-fraction discriminator (CFD) time relative to the Pixie16 clock, not as a waveform. Inthe second run, the CFD timing signal from the BaF neutron monitor was also captured and usedto correct for the drift in the phase relationship between the beam timing signal and the arrival ofthe beam burst at the target (section 3.6). The data acquisition system did not have a global trigger,instead each channel was digitized whenever the signal exceeded a threshold. The thresholds forthe detector channels were set low enough to measure single photoelectrons.To prevent excessive triggering, a signal over threshold in the liquid scintillator or NaI:Tlchannels required a coincident signal in one of the NaI:Tl or liquid scintillator channels, respectively.An ideal trigger for the beam timing signal would require a coincident liquid scintillator and NaI:Tlchannel. However, due to a limitation of our digitizer control implementation, we instead required atwo-fold multiplicity of coincident signals amongst all of the liquid scintillator and NaI:Tl channels.This meant that events where both NaI:Tl channels triggered dominated the data stream, howeverthese events could be readily rejected in later analysis. When it was acquired, the rate in the BaF detector was low enough that self-triggering was sufficient. The digitized waveforms were processed to estimate the pulse arrival time, charge, and a particleidentification metric for each detector channel.– 4 – .1 Trigger Time
A significant amount of RF noise was observed in the NaI:Tl waveforms. To mitigate this noise andreduce the analysis threshold, the waveforms were digitally filtered using a Wiener filter. A Wienerfilter requires an estimate of the power spectral density of the signal ( S ), and the noise ( N ), whichtogether give an estimate of the signal to noise ratio in Fourier space. The filter is then defined as H = S/N
S/N so that frequencies where the signal-to-noise ratio is high are passed, and frequencies where it islow are attenuated.Since noise was most problematic for low amplitude pulses, the filter was optimised for singlephotoelectron events. A population of single photoelectron waveforms was selected using a timeover threshold cut. These were then aligned in time by up-sampling the waveforms and aligningto the maximum. The averaged power spectral density in the trigger region about the singlephotoelectron signal was used as an estimate of S , and the averaged power spectral density of thepre-trigger region was used as an estimate of N .For the adopted analysis threshold of ∼ The low charge NaI:Tl events tend to be composed of single photoelectrons that are well-separatedin time owing to the 250 ns NaI:Tl decay time, which is slow relative to the few ns signal width froma single photoelectron. Thus, a simple fixed integration approach performed poorly for low-chargeNaI:Tl events and was unable to resolve the single photoelectron peak due to integrating the noisebetween well-separated single photoelectron events.To avoid this, only regions around where the waveform exceeded a threshold were integrated.The threshold was chosen to be approximately 0.5 photoelectrons, to optimise the energy resolution.This algorithm is outlined schematically in figure 2a. This approach performed much better than afixed integration window at low charge, allowing single photoelectron events to be identified (figure2b), while maintaining identical performance at high charge. This algorithm was adopted for allNaI:Tl energy estimates in this analysis. – 5 – a) (b)
Figure 2 : (a) The pulse charge estimate integrated the charge only in regions where the waveformexceeded a ∼ The NaI:Tl detection efficiency drops below 1 at low energies due to there being few scintillationphotons per event. Thus, the numbers of events measured near the threshold is an underestimationof the true spectrum, which can lead to incorrect predictions of the quenching factor [13]. We haveestimated the detection efficiency using coincidences between the two NaI:Tl photomultipliers; theefficiency of an individual channel is given by (cid:15) i = N C N j , (3.1)where N C is the number of observed coincidences between channel i and j and N j is the number ofcounts in channel j evaluated in an identical energy range. The coincident efficiency is the productof the two efficiencies and the efficiency of the union of the two channels is given by: (cid:15) U = (cid:15) + (cid:15) − (cid:15) c , (3.2)where (cid:15) c is the coincident efficiency. The union of all events passing the cuts was used for theanalysis of our measurements, and the resulting spectra were corrected by an energy-dependent fitto this efficiency given by f ( E ) = 1 − e − λE . (3.3)This model function is motivated by the fact that the detection efficiency ought to be 1 minus thePoisson probability of detecting zero photoelectrons. The fit is given in figure 3. The model fit isnot perfect, and there are differences in the higher energy range that may to stem from the difficultyof measuring unit efficiency with limited statistics. We imposed an analysis threshold at 3 keV inpart to mitigate the degraded fit at low energies. The effect of the fit to the efficiency does not have– 6 – Energy (keV) D e t e c t i o n E ff i c i e n c y Figure 3 : The measured energy-dependent detection efficiency, as well as its fit (see section 3.3).The fit values were used above the vertical dashed line.a significant systematic effect on our quenching factor results: a reanalysis of the quenching factorswithout any efficiency correction resulted in a ≤ . change to the quenching factor. A well-known attribute of liquid scintillators is that the time distribution of the emitted scintillationlight depends upon the stopping power of the particle that deposited energy in the material. Thisproperty can form the basis of particle identification.There are a number of techniques possible for particle identification. We used the chargecomparison method, with the metric Q f /Q t , where Q f is the charge in the first part of the pulse andQ t is the total pulse charge. Fixed integration windows about the pulse arrival time were used toestimate the charge. The free parameters for the algorithm are the integration lengths, which weredetermined using a parameter grid search to be optimal at 10 ns for Q f and 100 ns for Q t . Figure4 illustrates the excellent particle identification performance of the scintillator. We applied a quality cut to all measurements to remove spurious noise events in the NaI:Tlphotomultiplier channels. This charge asymmetry cut required that the difference in the proportionof charge detected between the two photomultipliers not be too great, such that | Q − Q | Q + Q < . where Q and Q are the charges measured in the two NaI:Tl photomultipliers.– 7 – Q t (AU) Q f / Q t Figure 4 : Particle identification in the liquid scintillation detector at 67.5 degrees, with a protonbeam energy of 5.2 MeV. See text for details of the definition of Q f and Q t . The shaded red regionindicates the selected nuclear recoil events. The upper band is from electron recoils, predominantlyfrom γ rays.We also applied a separate cut requiring that the two NaI:Tl signals arrive within 250 ns ofeach other. An analysis of the first measurement campaign revealed that the phase relationship between thebeam sync signal and the arrival of the proton pulse at the target changed over time. In the secondcampaign this relationship was tracked by recording the timing signal from the BaF signal, whichyielded time spectra as shown in figure 5a. The two spectra, which represent 20 second slices ofthe run, show an obvious shift in the arrival time of the pulse relative to the beam sync signal.To track this shift, we have performed fits to every analogous timing spectrum representingsequential 20 second slices of beam irradiation. The fit model used two Gaussian functions, andcuts were made to ensure the fit converged and achieved an acceptable goodness-of-fit. For timeperiods where the fits were not of sufficient quality, the mean of adjacent fits was used instead. Ingeneral, poor fits were only returned when the beam was not irradiating the target. The resultingdistribution of the position of the prompt gamma signal during the 3 MeV beam measurementis shown in figure 5b. There appears to be continuous small-scale variability as well as suddenstep-like changes. The origins of these timing variations are not known, however they have beencorrected for in our analysis. – 8 – s (a) (b) Figure 5 : The drift of the beam pulse arrival time relative to the sync signal. (a) The timespectrum of prompt gamma rays and neutrons measured in the BaF monitor, relative to the syncsignal (so delayed events occur at earlier times). The larger and smaller peaks in each spectrumare due to prompt gamma rays and neutrons, respectively. The blue and red histograms representarbitrary 20 second acquisitions within the 3 MeV beam measurement, and illustrate the level of driftencountered. (b) The fitted location of the beam arrival time across the entire 3 MeV measurementperiod. The combined energy estimate for all events was based on the summed charge contributions ofboth NaI:Tl photomultipliers in both nuclear and electron recoil measurements. The electron recoilenergy scale was set using gamma ray sources:
Ba,
Cs,
Eu, and
Am. The spectral peakswere fit to a Gaussian function on a linearly interpolated background.We corrected the measured charge in the NaI:Tl detector for the known non-linear electronrecoil response [14, 15] to obtain a linearised electron-equivalent energy scale. A small non-linearity remained evident at high energies, which was accounted for by using a piecewise linear fitbased on a linear calibration of the low-energy peaks below 80 keV and a separate calibration basedon the high energy peaks above this energy.A shift in the energy scale was evident during the irradiation, based on the movement of thelocation of the gamma ray associated with the first excited state of
I at 57.6 keV. These couldhave been gain variations in the photomultipliers or due to differences in the spatial location ofsource calibration and beam events. The shifts were accounted for by tracking the position ofthis peak at the different beam energies and applying a uniform energy scale correction. Only thelow-angle liquid scintillator detectors with iodine recoil energies less than 5 keV were used in thegain correction, so that the contribution from the iodine recoil energies did not strongly affect thepeak location. (Note that the low quenching factor reduces the correction to approximately 250 eV.)The uncertainty of the corrections was limited by the number of inelastic scattering events collectedat these angles, which translates to an additional 1% uncertainty in the quenching factor results toaccount for the variability in the energy scale. – 9 – a) T1/2=1457±484 s (b)
Figure 6 : Background measurements of the NaI:Tl crystal post-irradiation. (a) The energy spectrumin the first 600 and last 600 seconds of the 2600 second acquisition. The peak at 1461 keV is due tonaturally occurring K background. (b) The time dependence of the event rate between 450 and2000 keV, with a fitted exponential decay.The energy resolution of the linearised NaI:Tl energies was fitted with a function proportionalto the square root of the energy, and this function was used to smear the simulated spectra (section5) with the experimental resolution.Following the irradiation, we captured background measurements of the NaI:Tl to investigatethe activation of the crystal. The background rate varied throughout the 2600 second measurement,especially in the 450 keV to 2000 keV range (figure 6a). A fit of the decay to an exponential functionwith constant background in this region gives a half-life of 1457 ±
484 seconds, which is consistentwith the
I half-life of 1499.4 ± I is expected to be produced abundantly vianeutron capture on I. Nuclear recoils were selected using timing cuts and particle identification cuts in the liquid scintil-lator detectors.Figure 7 shows a time correlation plot (relative to the beam pulse time) for the 3 MeV protonbeam dataset at 90 ◦ recoil angle. The population lying along the diagonal line represents fast coin-cidences from events such as Compton-scattered gamma rays and cosmic ray showers. The excessalong this line (green circle, figure 7) is due to prompt gamma ray interactions in both detectors. Thenuclear recoil events are delayed in both the NaI:Tl detector and the liquid scintillator detectors, andthe cut region used to select these events is shown on the figure. Given the measurement geometryused to generate the plot, the neutron arrival time should be delayed by ∼ ns and ∼ ns relativeto the prompt gamma burst for the NaI:Tl and liquid scintillator detectors, respectively. Thesepredicted times are consistent with the observed delays in figure 7). Analogous cuts were made forall detectors and beam energies. – 10 – igure 7 : The time correlation of events in the NaI:Tl and nuclear recoil-like events in the 90 ◦ liquid scintillator detector, for the 3 MeV proton beam. Larger values are more delayed with respectto the beam arrival for both detectors. The green circle highlights the early beam-correlated eventscaused by gamma ray scattering, and the red triangular section illustrates the cuts used to select thelater nuclear recoil events. We used the Geant4 toolkit, version 10.5 [16], to simulate the neutron interactions within theexperimental geometry. In order to properly simulate the low energy neutrons generated for thisexperiment, we used the "QGSP_BIC_HP" physics list. This list uses nuclear data to transportneutrons with energies less than ∼
20 MeV with high precision. Range cuts for neutrons and protonswere reduced to zero in order to invoke the Geant4 processes needed to accurately simulate lowenergy behaviour. Primary neutrons were emitted in a 1 mm radius beam spot from the centre ofthe LiF target.The energies and momenta of the neutrons were calculated using SimLiT [17]. SimLiTis a specialized Monte Carlo computational tool designed to calculate the neutron energy andangle distributions produced from Li(p,n) Be and Li(p,n) Be ∗ reactions for given proton beamand target parameters. SimLiT natively supports proton energies above the Li(p,n) Be reactionthreshold up to 3.5 MeV. This is lower than some of the energies used for our experiment. Therefore,we have extended the SimLiT energy range using neutron production cross sections from Liskien et al. [18] and LiF stopping powers from SRIM [19].Due to neutron emission from the Be excited state, the proton beam produces two quasi-– 11 –onoenergetic, forward-directed neutron beams. Neutrons from the Li(p,n) Be ∗ reaction makeup 0.8%, 20%, and 10% of the total neutrons produced for the 2.44, 3 and 5.2 MeV proton beamsrespectively. Hence the Li(p,n) Be ∗ reaction contributes a non-negligible portion of the totalneutrons produced, for all but the lowest energy. We simulated the unquenched recoil spectra by tagging neutrons that deposit energy in the NaI:Tland exit at a polar angle subtended by the front face of one of the liquid scintillator detectors usedin our geometry. Figure 8 shows the fraction of the simulated neutrons that scattered in the NaIdetector towards the scintillator placed at 67.5 degrees from the beam line for the initial 3 MeVproton energy; this represents 0.007% of the simulated events. Most of these events, approximately73%, involve neutrons that travel straight into the NaI detector. However, the remaining 27% areneutrons that interact with other materials (predominantly the collimator) prior to reaching the NaIdetector. Thus, they can have significantly different energies, as shown by the spread of eventsabove the blue dashed line in figure 8. Of these, neutrons that are not emitted in the forwarddirection towards the collimator can reach the detector by scattering off beamline and vacuumfittings surrounding the target. Of the neutrons that scatter into the 67.5 degree scintillator, 7% arefrom the Li(p,n) Be ∗ reaction.Figure 9 shows the energy spectrum associated with the neutrons that were scattered into theliquid scintillator placed at 67.5 degrees from the beam line. The non-resonant backgrounds infigures 9a and 9b are due to multiple neutron scattering, which makes up 23.8% of events. Smallerpeaks seen to the left of their larger counterparts (for example the small elastic Na peak centredat ∼
40 keV before the larger one centred at ∼
60 keV) in figure 9a are due to neutrons from the Li(p,n) Be ∗ reaction. The elastic Na and I peaks in the low energy region <10 keV, are mainlydue to events where neutrons lose energy in the collimator or aluminium casing before scattering inthe detector. This was determined by comparing simulations where the collimator and aluminiumcasing were removed from the geometry. Previous quasi-monoenergetic neutron studies of the quenching factor in NaI:Tl have fit the nuclearrecoil peaks with Gaussian functions to extract the quenching factor. In our analysis, we havedecided to fit the spectrum of simulated results to the experimentally measured recoil spectra. Thisapproach allowed us to leverage the information contained in the inelastic recoil peaks, and alsohelped to incorporate complicating effects such as multiple scattering within the NaI:Tl, scatteringfrom the collimator, and lower energy neutron emissions from the excited Be state. A comparativeanalysis involving fits with a simple Gaussian fit to the visible Na recoil peaks and linear backgroundgave quenching factor results that differed from the spectrum fit analysis by between 7% and 43%.The fitted peak location also depended upon the choice of the fit range. Furthermore, not all recoilspectra could easily be fitted with a Gaussian. This is especially true for the higher beam energydata, where the elastic scattering peak was broadened by lower energy components associated withemissions from the excited Be state and inelastic scattering events where the gamma ray escapedthe NaI:Tl. – 12 – igure 8 : Fraction of total emitted neutrons that scatter in NaI detector and proceed to the scintillatorplaced at 67.5 degrees from the beamline. 73% of neutrons are below the angle of acceptance at2.2 degrees; shown in blue. Events in the left group (labelled Be ∗ ) are neutrons produced fromthe Li(p,n) Be ∗ reaction, while in the right group (labelled Be) neutrons are produced from the Li(p,n) Be reaction
We applied multiplicative parameters to the simulated Na and I recoil energies to model thequenching factors. These parameters were adjusted so as to minimise the negative log-likelihoodof the resulting simulated distribution, given the measured data.An important feature in the recoil spectra is the peak at ∼
58 keV from the inelastic excitationof
I. The energy deposit from such events is the sum of the 57.6 keV gamma ray and therecoiling nuclear energy deposit. Since the I recoil energies were low due to the experimentalkinematics, and the I quenching factor is low, the location of this peak did not move during ourmeasurements. However, this feature was consistently under-predicted by the Geant4 model at allneutron energies, which led to poor fits. To mitigate this effect, a third parameter was introducedin the likelihood model which weighted these inelastic scattering events. Figure 10 shows the fittedinelastic scattering weights, which were broadly consistent with being near 2 for the 2.44 MeV and3 MeV beam energies, and higher for the higher beam energies.The fitted detector resolution was applied to the simulated spectra to allow comparison withthe experimental results. Our early analysis of the resolution smearing showed that smearing thequenched nuclear recoil energies tended to under-predict the widths of the associated features in the– 13 – a) (b)
Figure 9 : Simulation of the unquenched energy deposited in the NaI detector due to neutrons thatscatter to the angle subtended by the liquid scintillator placed at 67.5 degrees from the beamline,for the 3 MeV proton beam energy. (a) and (b) show the Na and I spectra, respectively. The redand blue lines show events where the neutron scatters elastically or inelastically, respectively. Wehave applied an energy threshold of 0.1 keV to remove thermal neutron events. The error bars showthe standard Poisson error. Each of the spectra exhibits peaks that correspond to single scatteringevents. For the inelastic spectra, the peaks at lower energy correspond to events where the gammaray escapes the detector.
Figure 10 : The event weighting factors for the inelastic
I excitation. Only the fits used in thefinal analysis are shown. The different series denote the different proton beam energies.– 14 – .0000.0050.0100.0150.0200.0250.030 I Q u e n c h i n g F a c t o r Na Quenching Factor I n e l a s t i c E v e n t W e i g h t i n g F a c t o r Figure 11 : Every Markov Chain Monte Carlo sample for the fit to the NaI:Tl recoil spectrumobserved in coincidence with the 135 ◦ liquid scintillator with a 3 MeV proton beam energy. 10000Monte Carlo trials were used in the sampling. The distribution of samples obtained by this fit wereused as an estimator of the posterior distribution of the fit parameters.spectrum. Instead we have applied the energy resolution to the unquenched nuclear recoil energies;this gives a better fit to our measurements.Once a maximum likelihood parameter estimate was made, we sampled via Markov ChainMonte Carlo in order to assess the posterior distribution of quenching factors and inelastic peakbiasing allowed by the data. Figure 11 demonstrates this sampling for a single beam energy andliquid scintillator detector. Typical fits showed only weak correlations between the fit parameters,with the strongest correlation between the Na quenching factor and the peak biasing shown in thefigure, with a correlation coefficient of 0.8. This correlation arose from a trade-off between thelocation of the Na elastic recoil peak and the inelastic I peak intensity, when the peaks were closein energy.Recoil spectra were excluded from our analysis where the fitted quenching factors predictedrecoil peaks below the 3 keV analysis threshold, since there is less informative data available toproperly constrain the fit. Fits where the Na recoil peak coincided with the 57.6 keV inelastic– 15 –xcitation peak were also excluded, since this degeneracy also led to poor fitting constraints. Thehighest energy nuclear recoil distributions above the maximum fitted energy of 401 keV were alsoexcluded, since these gave poor fits to the simulated spectra – the Na elastic recoil peak tended tobe broader than expected from simulation, and this flat spectrum did not adequately constrain thefit parameters. Figure 12 shows some of the fitted recoil spectra across the range of Na recoil energies studied. Thefits describe the data reasonably well around the recoil peak and the inelastic excitation lines. Thelow-energy rise is due to the tail of the I elastic scattering peak.We chose the fit range to be as large as possible to make full use of the measured data, includingthe inelastic excitation peaks, which carry information about the iodine recoils that would otherwisefall below threshold. An exception is the highest energy fitted Na recoils, associated with the 112.5 ◦ detector at 5.2 MeV beam energy, and the 90 ◦ detector at 6 MeV beam energy. In these casesthe elastic scattering distributions were quite broad and these fits tended to find quenching factorswhich optimised the Na elastic scattering to other features of the spectrum. In these cases the fitwas performed using a restricted energy range to avoid this problem (figures 12(g) and 12(h)).The Na recoil quenching factors are given in figure 13 and table 1. The Na recoil energies aretaken to be the mean of the elastic scattering peak predicted by SimLiT.Our fits also included the I quenching factor. These were constrained by the location of the mainiodine inelastic excitation peak. However, these quenching factors tended to have large uncertaintiesand are not reported. The covariance between the I and Na quenching factors was not large for anyof the fits, so that the poor constraints on the I quenching factors do not strongly influence the Naresults.The systematic errors that could affect our results have been minimised as much as possible.The event selection for tagging nuclear recoil events was deliberately chosen to be relativelyconservative to avoid contamination from electron recoil events, and we assume that the smallremaining proportion of such events does not affect the quenching factor in a significant way. Animperfect knowledge of the resolution and energy calibration at low energies can also give riseto systematic errors. The resolution model was extrapolated below the energy threshold, whichreproduced the observed rise at low energy associated with the tail of the iodine elastic scatteringpeak. The effect of these sources of systematic error is difficult to assess quantitatively, however weassume such errors are small owing to our relatively high analysis threshold. Our results are self-consistent across two measurement campaigns and different neutron beamenergies. They are also consistent with previous measurements insofar as they show a decreasingquenching factor with recoil energy. There is some tension amongst the previous studies, whichmay arise from a crystal-dependence to the quenching factor or systematic differences between themethodologies employed in those studies. Our measurement results show a mild preference for thelower quenching factors amongst those previously measured, though we were limited in our abilityto examine these differences by the relatively high energy threshold.– 16 –a)
20 40 60 80 100
Energy (keV) (b)
20 40 60 80 100
Energy (keV) (c)
20 40 60 80 100
Energy (keV) (d)
20 40 60 80 100
Energy (keV) (e) (f)(g)
40 60 80 100 120 140 160
Energy (keV) (h)
Figure 12 : A selection of spectrum fits, where the fit parameters are taken to be the mean of theposterior distribution given by the Markov Chain Monte Carlo sampling. The red and blue stackedhistograms denote the Na and I recoil contributions, respectively, while the green points are themeasured spectra. All spectra are normalised to unit area. The Na elastic scattering peak can beseen to move to higher energy and become broader as the recoil energy is increased. Note thedifferent energy scale in the highest recoil energy plots (g) and (h).– 17 – igure 13 : The sodium recoil quenching factors in NaI:Tl inferred from our measurements, togetherwith previous results [9–11]The main outcome of this work is to demonstrate a spectrum fitting methodology for themeasurement of the Na quenching factors in NaI:Tl with Li(p,n) Be neutrons. Our fits exhibitexcellent agreement with the experimental results over a wide energy range. We intend to performsimilar measurements using this methodology on high quality NaI:Tl off-cuts from the SABRE ex-periment’s crystals [3] in an optimised enclosure and measurement geometry. These measurementswill allow a lower energy quenching factor determination and help to probe any crystal dependence.
We would like to thank the technical staff at the ANU Heavy Ion Accelerator Facility, who fabricateddetector housings, built the beamline hardware, and carefully aligned the detectors in the correctgeometry. We would also like to thank A. Duffy from Swinburne University for the use of hisNaI:Tl crystal. This work was supported by the Australian Research Council Discovery Programthrough project number DP170101675.
References [1] P. A. Zyla et al.,
Review of particle physics , Progress of Theoretical and Experimental Physics (2020) 083C01. – 18 –cattered neutronangle Proton BeamEnergy (MeV) Recoil Energy (keV) Quenching Factor67.5 ◦ ± ± ◦ ± ± ◦ ± ± ◦ † ± ± ◦ ±
10 0.235 ± ◦ ± ± ◦ ±
12 0.194 ± ◦ ± ± ◦ ± ± ◦ † ± ± ◦ ± ± ◦ † ± ± ◦ ± ± ◦ † ± ± ◦ ±
24 0.309 ± ◦ ±
15 0.305 ± Table 1 : The measured Na quenching factors at the elastic scattering recoil energies predictedby SimLiT. The measurements with the 3 MeV proton beam were repeated across two differentmeasurement campaigns. The uncertainties are statistical. For the repeated measurements, † indicates results from the second measurement campaign. [2] R. Bernabei, P. Belli, A. Bussolotti, F. Cappella, V. Caracciolo, R. Cerulli et al., First modelindependent results from DAMA/LIBRA-phase2 , Nuclear Physics and Atomic Energy (2018)307–325.[3] M. Antonello, E. Barberio, T. Baroncelli, J. Benziger, L. J. Bignell, I. Bolognino et al., The SABREproject and the SABRE Proof-of-Principle , The European Physical Journal C (2019) .[4] G. Adhikari, P. Adhikari, E. B. de Souza, N. Carlin, S. Choi, M. Djamal et al., Search for a darkmatter-induced annual modulation signal in NaI(Tl) with the COSINE-100 experiment , PhysicalReview Letters (2019) .[5] J. Amaré, S. Cebrián, I. Coarasa, C. Cuesta, E. García, M. Martínez et al.,
First results on dark matterannual modulation from the ANAIS-112 experiment , Physical Review Letters (2019) .[6] A. Hitachi,
Liquid rare gases and inorganic scintillators for WIMP searches , in
The Identification ofDark Matter , pp. 344–349, 2005, DOI[ ].[7] R. Bernabei, P. Belli, V. Landoni, F. Montecchia, N. W. Di, A. Incicchitti et al.,
New limits on WIMPsearch with large-mass low-radioactivity NaI(Tl) set-up at Gran Sasso , Physics Letters B (1996)757 .[8] J. I. Collar,
Quenching and channeling of nuclear recoils in NaI(Tl): Implications for dark-mattersearches , Phys. Rev. C (2013) 035806. – 19 –
9] J. Xu, E. Shields, F. Calaprice, S. Westerdale, F. Froborg, B. Suerfu et al.,
Scintillation efficiencymeasurement of Na recoils in NaI(Tl) below the DAMA/LIBRA energy threshold , Phys. Rev. C (2015) 015807.[10] T. Stiegler, C. Sofka, R. C. Webb and J. T. White, A study of the NaI(Tl) detector response to lowenergy nuclear recoils and a measurement of the quenching factor in NaI(Tl) , 2017.[11] H. Joo, H. Park, J. Kim, J. Lee, S. Kim, Y. Kim et al.,
Quenching factor measurement for NaI(Tl)scintillation crystal , Astroparticle Physics (2019) 50 .[12] X. LLC.,
Pixie-16 , 2020 (accessed Sept 15, 2020).[13] J. I. Collar,
Light WIMP searches: The effect of the uncertainty in recoil energy scale and quenchingfactor , arXiv [astro-ph.IM] (2010) 1010.5187.[14] B. D. Rooney and J. D. Valentine, Scintillator light yield nonproportionality: calculating photonresponse using measured electron response , IEEE Transactions on Nuclear Science (1997) 509.[15] I. V. Khodyuk, P. A. Rodnyi and P. Dorenbos, Nonproportional scintillation response of NaI:Tl to lowenergy x-ray photons and electrons , Journal of Applied Physics (2010) 113513[ https://doi.org/10.1063/1.3431009 ].[16] S. Agostinelli, J. Allison, K. Amako, J. Apostolakis, H. Araujo, P. Arce et al.,
Geant4—a simulationtoolkit , Nuclear Instruments and Methods in Physics Research Section A: Accelerators,Spectrometers, Detectors and Associated Equipment (2003) 250 .[17] M. Friedman, D. Cohen, M. Paul, D. Berkovits, Y. Eisen, G. Feinberg et al.,
Simulation of the neutronspectrum from the Li(p,n) reaction with a liquid-lithium target at Soreq Applied ResearchAccelerator Facility , Nuclear Instruments and Methods in Physics Research Section A: Accelerators,Spectrometers, Detectors and Associated Equipment (2013) 117 .[18] H. Liskien and A. Paulsen,
Neutron production cross sections and energies for the reactions Li(p,n) Be and Li(p,n) Be ∗ , Atomic Data and Nuclear Data Tables (1975) 57 .[19] J. F. Ziegler, M. Ziegler and J. Biersack, SRIM – the stopping and range of ions in matter (2010) , Nuclear Instruments and Methods in Physics Research Section B: Beam Interactions with Materialsand Atoms (2010) 1818 .(2010) 1818 .