Variance reduced monte-carlo simulation for neutrons in XENON

(1)

MSc Physics and Astronomy

GRAPPA

Master Thesis

Variance reduced monte-carlo simulation for neutrons in

XENON

by

Gijs Leguijt

11000279

June 2020

60EC

September 2019 – June 2020

Supervisor/Examiner:

Examiner:

Prof. dr. A.P. Colijn

Dr. M. Vreeswijk

(2)

Neutron vrt-simulation 2

Summary

Multiple independent observations indicate there is more mass than we can see. These ob-servations can be explained by imposing a new particle: dark matter. As the additional mass has only been seen through its gravitational effect on ordinary matter, it must have low cross-sections for reacting with Standard Model particles.

To make the first non-gravitational discovery of a dark matter particle big detectors, like XENONnT, are build which consist of large target volumes that are constantly monitored. As the probability that a dark matter particle interacts within the detector, they are placed in low-background places, e.g. underground or within mountains, to make the detectors as quiet as possible. In addition, they are constructed using radio-quiet materials and they include systems to actively remove background signals from the detector.

Still, the majority of the events in the detector are background events, which makes it hard to spot the rare dark matter signals. To increase the probability that a dark matter signal is found, background simulations need to be improved. This would allow the identification of smaller excesses and thus might enable a dark matter discovery.

A method by Kesber et al. enables a considerable acceleration of the background simula-tion for photons, meaning better models can be achieved. This thesis will extend this method to include neutrons. As neutron-interactions are more complicated than photon-interactions this work will use the Geant4 simulation package. It is important to include neutrons in the improved background model as they closely resemble potential dark matter signals.

It is shown that the method of Kesber et al. can be implemented in Geant4 to yield the correct energy spectrum. For 1 MeV neutrons this results in a 12 times more effective simulation. By improving the energy of the neutrons and restricting the allowed energy deposits this factor can be increased further.

(3)

3.3 Geant4 . . . 21 4 Photon validation 22 4.1 Photon interactions . . . 22 4.1.1 Compton-effect . . . 22 4.1.2 Doppler broadening . . . 23 4.1.3 Photo-electric effect . . . 24 4.1.4 Pair production . . . 24 4.2 Validation photon-spectrum . . . 25 4.3 Energy cut . . . 27 4.3.1 Maximum angle . . . 27 4.3.2 Look-up table . . . 28 5 Neutrons 32 5.1 Origin of neutrons . . . 32 5.1.1 Cosmogenic neutrons . . . 32 5.1.2 Radiogenic neutrons . . . 33

(4)

5.2 Neutron interactions . . . 34

5.2.1 Elastic scatter . . . 34

5.2.2 Neutron capture . . . 36

5.2.3 Inelastic scatter and neutron induced fission . . . 36

5.2.4 Isotope dependent cross-section . . . 37

5.3 Validation (neutrons) . . . 39

6 Acceleration 41 6.1 Relative error . . . 41

6.2 Acceleration result . . . 43

(5)

1 Introduction

In 2018 XENON1T, at that moment the largest liquid xenon detector or earth [1], stopped taking data. Data acquisition was not stopped due to a lack of interest, instead, the detector had to make room for its successor: XENONnT, an even larger liquid xenon detector, to continue the search [2]. It will be the fourth in a line of detectors searching for dark matter [1], [3], [4]. At the same time, similar detectors (e.g. Lux Zeplin [5]) are being constructed trying to make the first direct detection of a dark matter particle. Evidently, the dark matter search is an active field of research.

Even though there is no dark matter in the Standard Model, there is both astrophysical and cosmological evidence for its existence as discussed in section 2 [6]. From these observa-tions the total dark matter fraction of the universe can be calculated and it turns out that for every kilogram of (ordinary) Standard Model particles there are four kilograms of dark matter particles [6]. Despite this large amount of dark matter, all proof for dark matter is based on gravitational effects of dark matter on ordinary matter. The first non-gravitationally based measurement of a dark matter particle is still to be made.

The lack of a direct detection of a single dark matter particle is due to the low cross-section between dark and ordinary matter. Regardless of the detection mechanism, subject of section 2, the probability for a dark matter particle to react within the detector is really small. On the other hand, the probability for ordinary particles to leave a signal in the detector is high. Detectors therefore have to put great effort into reducing the background signal to a minimum to be able to see the low dark matter signal above the noise [2].

To this end, the XENON-detector is build within a mountain which shields it from the outside world. Also a screening campaign was issued to ensure the detector was build with radio-clean material [2]. Still, the background signal exceeds the dark matter signal by orders of magnitude. To reduce the background further, a good understanding of the detector’s response is needed. Simulations enable a better identification of background signals as well as an estimation of the number of background events in the detector. This means excesses of events above the background can be recognised better. With and increased understanding of the background signal, lower excesses can be identified [7].

As searching for an excess above the estimated background signal has not yet resulted in a dark matter detection, exceedingly large detectors are build increasing the probability for a dark matter particle to react in the detector. Another, more practical, way to increase the probability of noticing a dark matter particle in the detector is by increasing the number of particles in the background models. As these models are based on large simulations [2], [8], [9], it seems infeasible to make the number of particles even larger. However, if instead the simulations can be made more effective, the same simulation time can lead to a more accurate model.

(6)

Increasing the efficiency of the standard background simulations is the goal of this thesis. It is a continuation of the work by Kesber et al. , who tested the method on photon simulations [10]. The working principle of the method is the realisation that most simulated particles are not important for the dark matter search. For example, the probability that a dark matter particle scatters within the detector is small, meaning that the probability for the same particle to scatter twice within the detector is negligible. As these unwanted events form a large fraction of the simulated events, the efficiency of the simulations can be increased by forcing the particles to be interesting events (e.g. scatter only once within the detector). To keep a correct prediction of the number of events, particles are given a weight which keeps track of how likely the specific event is. A more thorough explanation of this variance reduced Monte-Carlo technique (vrt MC) is given in section 3.

As opposed to [10], which uses a python-based code, this work will use the Geant4 sim-ulation software [11]. Geant4 is a dedicated software package commonly used in the field of particle physics [12], [13]. It embodies an integrated neutron-physics package, enabling the simulation of neutrons which behave more complicatedly than photons [14]. In addition it allows for the simulation of complex geometries which can not be achieved with the current python simulation code. However, before discussing the results of the neutron simulations in section 5, the variance reduced code code will be validated with a photon simulation in section 4. Finally, section 6 will show the effect of reducing the variance of the neutron simulations.

(7)

2 The XENON-experiment

This section will give a short overview of the dark matter problem in section 2.1. After giving evidence for the existence of dark matter in section 2.1.1, section 2.1.2 will discuss the most common dark matter model. Section 2.2 will explain the working principle of the detector, which mainly focuses on the model of section 2.1.2.

2.1 Dark matter

2.1.1 Evidence for dark matter

Radial velocity distribution in galaxies and galaxy clusters

As mentioned in the introduction, all evidence for dark matter is based on gravitational effects. One phenomenon that can be explained by dark matter is the radial dependence of velocities in galaxies and galaxy-clusters [15]. In 1957 measurements of stars in the galaxy M31 showed higher rotational velocities than could be explained by the visible matter [16]. Far away from the center of the galaxy, where increasing the distance no longer increases the encapsulated mass one would expect the velocity to drop off with increasing radial distance. Instead measurements show a flattening of the velocity distribution as can be seen in figure 1.

Figure 1: Rotation curve of stars in M33 taken from [17]. Both the expected and measured curve are shown. The expected curve starts to fall off once the majority of the galaxy mass is encircled. Measurements however indicate a flat distribution which indicates the encircled mass still increases. An extended dark matter halo could describe the shape of the measured curve as it would mean the encircled mass increases even at high distances from the galaxy center. [15]

(8)

in an extended halo around the center of the galaxy [15]. This would ensure a growth of encircled mass even at large distances. As this additional mass is not visible it would be dark matter.

This argument is not limited to lone galaxies. It turns out a similar effect is measured in rotational curves in galaxy clusters [18]. The velocity of galaxies with respect to the center of the cluster indicate additional mass in between the different galaxies, extending the previous argument.

Cosmic Microwave Background

A different proof for dark matter comes from precise measurements of the Cosmic Microwave Background (CMB) [19]. In the early universe photons and matter were strongly coupled. However, due to the weak coupling between dark matter and ordinary matter, dark matter could decouple from ordinary matter at an early point in time [20]. As a result, dark matter was able to clump together due to its gravity leading to potential wells. As photons were coupled to ordinary matter they fell into the gravitational well, increasing the photon pressure. As the pressure increased it started to push matter out of the gravitational wells, creating waves. These acoustic oscillations left an imprint on the photons when they decoupled [21].

The resulting CMB-spectrum shown in figure 2 is hard to explain without potential wells from dark matter.

Figure 2: Power spectrum of the Cosmic Microwave Background taken from [22]. The oscillation is the result off a balance between photon pressure and gravity. It is hard to explain the shape of the spectrum without imposing gravity wells originating from dark matter.

Structure formation

As well as leaving an imprint on the CMB, dark matter potential wells play an important role in structure formation [15]. Quantum density fluctuations in the early universe form the seeds of galaxy clusters. Over-fluctuations will attract surrounding matter which in turn leads to higher differences and pile up of matter. However, to grow from a quantum fluctuation to a galaxy cluster requires time. At the same time the universe expands, lowering the matter

(9)

density, slowing down this growth. Simulations show that gravity from ordinary matter is not enough to lead to a structure comparable to the current universe [6].

As dark matter decouples earlier from photons, it can clump together before ordinary matter can. The resulting gravitational wells can result in a head-start for structure formation through ordinary matter. Simulations that include dark matter can successfully produce universes comparable to ours [23].

Bullet Cluster

Although dark matter would solve all previous phenomena, there has been an effort to explain these observations without the need for dark matter [24]–[27]. Instead of adding additional matter modifications can be made to nature’s laws to yield similar effects. However, explaining the bullet cluster by modifying the law of gravity proves difficult [28]. The bullet cluster is formed by the collision of two clusters. As the spacing between stars is orders of magnitude larger than their radii there are no collisions between stars and they continue without slowing down. The interstellar gas, which carries the majority of the mass of a cluster, is however much denser meaning it is slowed down heavily when the clusters collide.

As a result, one would expect two galaxies moving apart after the collision while the interstellar gas trails behind. This would mean that most mass would become separated from the stars. However, measuring the gravity in the galaxies using gravitational lensing, most mass is centered with the stars [28]. This means there has to be another mass component next to the gas that does not feel any effect from a galaxy collision. A dark matter halo would make a perfect candidate due to the small cross-section of dark matter with other particles.

Figure 3 shows a photo of the bullet cluster. The gas, identified by infra-red emission is shown in red while the majority of the mass is shown in purple. Explaining the separation of gravity from the visible mass is hard to achieve by modifying gravity. Similar collisions of galaxies form strong evidence for dark matter [29].

(10)

Figure 3: Photo of the bullet cluster by [30]. The bullet cluster is a collision of two galaxies. As opposed to stars, which are not affected by the collision, the interstellar gas (shown in red) is slowed down. Although the interstellar gas forms most visible mass of the cluster, gravitational lensing measurements show that most mass is centered around the stars, shown in purple. Dark matter, which would not be affected by the collision of galaxies and continue moving with the stars, is a perfect candidate to explain this feature. However, trying to explain the separation between the gas and gravity is hard to accomplish by modifying gravity.

2.1.2 Weakly Interacting Massive Particle

As no dark matter has been directly observed yet, there is a variety of theories describing the dark matter particles [31]. The main focus of the XENON-experiment are Weakly Interacting Massive Particles (wimps) [32]. This is a class of potential dark matter particles that exhibit the following properties:

1. they are massive compared to SM-particles,

2. they have a long half-life as they originated from the Big Bang and are present still, 3. they interact weakly with SM-particles, where weakly does not necessarily imply through

the weak-force.

Candidate wimps follow from different theoretical frameworks (e.q. super-symmetry and universal extra dimension models [33], [34]). These models predict the existence of one or multiple additional particles where the lightest new particle is stable enough to form a wimp-candidate.

Although wimps do not necessarily interact through the weak force, cosmological models

imply annihilation cross-sections of hσvi ≈ 3 · 10−26 _cm3_s-1 _{[35]. A wimp reacting through}

the weak force would have the correct cross-section for a mass of the order of 100 GeV/c2

which is about the mass expected to explain the phenomena mentioned in section 2.1.1 [36]. This so-called wimp-miracle makes the wimp an attractive dark matter candidate.

(11)

2.2 The XENON-detector

In 2006 the XENON collaboration commissioned their first detector: XENON10 [3]. It was followed up by XENON100, XENON1T and XENONnT which is currently in an advanced state of construction [2]. As XENONnT will use 8.3 tonnes of xenon, it will be orders of magnitude larger than XENON10 which used only 14 kg of xenon. Although there are plans for another upgrade, named DARWIN, these are only in an early stage of planning which is why this work will focus on the XENONnT-detector.

Like its predecessors, XENONnT is positioned in the underground facilities of the Labo-ratori Nazionali del Gran Sasso (LNGS) [37]. The underground laboratory is build within a mountain which serves as a natural shield against cosmic rays. However, as the rock-material itself emits particles the detector is submerged in a water tank which reduces the background signal further [1].

The detector itself is shown in figure 4. It contains both of liquid and gaseous xenon which will be explained in section 2.2.1. One of the reasons to choose xenon as a target material is that it has a relatively high mass which raises the probability of wimp-collisions in the detector. At the same time this high mass makes it optimal for shielding purposes [35]. Therefore, 4.3 tonne of xenon is used as a background shield (partly passive and partly active) [2]. This leaves 4 tonne of xenon in the inner part of the detector to search for dark matter signals which is called the Fiducial Volume (FV). This region is in no way separated from the rest of the xenon and is solely a region defined within the xenon. Although its exact shape depends on multiple factors (e.g. field strength and background activity) it is represented by a cylinder of radius 57 cm and height of 134 cm in this work.

(12)

Figure 4: Schematic view of the XENONnT detector taken from [2]. The 8.3 tonnes of xenon are contained within an outer vessel. However, 4.3 tonnes are used as a (partly) active shield against background signals. At the top and bottom a PMT-array is used to detect the light from signals in the detector. The sides are covered with field shaping wires which create an electric field to drift electrons upwards as explained in section 2.2.1.

It is important to note that the xenon in the detector is natural xenon. This means it consists of 9 different isotopes of which the natural abundances are shown in table 1. Different isotopes have different cross-sections for neutrons with will be shown more thoroughly in section 5.2.4. It is therefore important to take all the different isotopes into consideration while simulating particles in the detector.

(13)

Isotope Abundance Xe-124 0.10 % Xe-126 0.09 % Xe-128 1.91 % Xe-129 26.4 % Xe-130 4.1 % Xe-131 21.2 % Xe-132 26.9 % Xe-134 10.4 % Xe-136 8.9 %

Table 1: Relative abundance of the various isotopes in natural xenon taken from Geant4.

2.2.1 Detection principle

As mentioned in the previous section both liquid and gaseous xenon are present within the detector. This separation of phases serves to split each event into two signals. When a particle hits a xenon-atom it will ionise the atom; 12.13 eV is needed to partly ionise the

xenon [35]. The ionised xenon will combine to form an excited Xe2-molecule called a dimer.

When the dimer de-excites it emits a photon. As the dimer has different energy levels than the individual xenon-atoms, the photon can move through the detector-volume without being absorbed [38]. Upon reaching the top or bottom it is detected by an array of photo-multiplier tubes which is the first detected signal of the event (S1).

A second signal results from the electrons that were freed upon ionisation of the xenon-atom. These electrons are drifted upwards towards the liquid-gas interface by an applied electric field. Upon reaching the gaseous xenon, the electrons have gained enough energy from the electric field to ionise the gas which leads to a second signal (S2). The temporal difference between the S1 and S2-signal can be used to estimate the detector depth at which the event happened as electrons originating lower in the detector need to drift longer before they reach the gas-layer [35]. Figure 5 shows how one event leads to two signals.

(14)

Figure 5: Schematic view of the signals originating from an event in the detector taken from [1]. As a particle hits a xenon-atom it will ionise the xenon-atom. This will lead to emitting both photons and electrons. The photons form the first detected signal (S1). The electrons are drifted towards the liquid-gas interface by an electric field. When they reach the gas-layer, the gas gets ionised leading to a second signal: S2. Both the S1 and S2-signals are detected by photo-multiplier tubes situated at the top and bottom of the detector.

(15)

3 Monte-Carlo simulations

As mentioned before, understanding the background in the detector is paramount to the search for dark matter. Overestimating the number of background events can cause a signal to disappear while a too optimistic model of your background can lead to a fake discovery. To ensure a good understanding of the detector, extensive simulations are run [2]. These simulations encapsulate detailed models of the geometry and propagate background particles through the detector until they have deposited all their energy or leave the detector again.

However, as will be shown in section 3.1, a large part of the available parameter space is cut out of the acquired data to enhance the signal-to-noise ratio [39]. This means simulated particles that deposited energy in this region of parameter space will have no impact on the final understanding of the signal. As a result, most simulated events are redundant.

If instead simulated events are forced to lie within the region of interest, the number of interesting events will increase which enables higher statistics and thus a better understanding of the background. To make sure this enforcement does not bias the simulation, particles will have to be assigned a weight that keeps track of how likely the simulated effect is. The weighting procedure is explained in section 3.2.

3.1 Signatures of dark matter events

3.1.1 Fiducial volume

Because the xenon in the experiment is used both for detection and shielding purposes, the detector volume is divided into two sections. An outer part which is used both for active and passive shielding and an inner part. As mentioned in section 2.2 this inner section is called the Fiducial Volume.

To simulate the same target mass as in XENONnT, a cylindrical FV was used with a radius of 57 cm and a height of 134 cm centered around the origin. In this way, particles originating at the edge of the detector (a cylinder with radius 65 cm and height 150 cm) will travel at least 8 cm before reaching the FV. As a 1 MeV neutron has an attenuation length of 10.8 cm, over 50 % of all 1 MeV neutrons are cut away by selecting only events in the FV.

3.1.2 Energy range of recoil spectrum

As opposed to photons, neutrons and other particles, a dark matter particle can not deposit large amounts of energy in the detector. If it were able to, it would have visible cosmological effects [40]. Therefore, all energy deposits above a certain threshold can be removed from the detector data without losing dark matter signals. For this work, particles were not allowed to deposit more than 200 keV.

(16)

3.1.3 Number of scatters

Just like high energy deposits would have a visible effect, a large cross-section for a reaction between a dark matter particle and an ordinary particle is ruled out by experiment [41]. Predecessors of XENONnT have set increasingly stringent constraints on the dark matter-matter cross-section showing that it is unlikely for a dark matter-matter particle to scatter inside the detector [42].

A particle that scatters more than once in the detector is bound to be background and will be removed from the data. For this reason, simulated particles were forced to react only once in the detector. However, to enable a better identification of double scatters the simulation code incorporates a handle to tune the number of times a particle will scatter. In this way, a closer look can be taken at e.g. neutron double scatters such that neutrons that react twice but at small spatial separation can be more accurately distinguished from single scatters that deposit more energy.

3.2 Accelerating Monte-Carlo

As emphasized in section 3.1, the majority of events in the detector is removed. This means also the majority of simulated particles is removed and thus the efficiency of the simulations is rather low, simulating 10000 1 MeV neutrons would only result in ≈ 32 useful events (which will be shown in section 6).

If the particles are forced to deposit a low amount of energy in the FV and leave the detector afterwards, the same 10000 1 MeV neutrons would result in 10000 useful events which gives improved statistics.

Still, it is important to keep the link between simulated particles and useful events. One does still expect 0.32 % of the particles to be interesting, and increasing that factor will result in an overestimation of your background. By adding weights that keep track of how likely an event is, the 10000 1 MeV neutrons will result in 10000 useful events with a total weight of 32.

To give the correct weight to an event, five contributions have to be taken into account, as shown in figure 6.

(17)

Figure 6: Schematic side view of the simulated detector, including a single scatter event. The different steps in the simulation are indicated by colour changes and the numbers indicate the corresponding weights. 1: reaching the FV. 2: scattering in the FV. 3: Selecting a scatter process. 4: energy deposit below a threshold. 5: leaving the detector. Also the three distances important for calculating the different weights are included.

3.2.1 Intersection of particle with FV

To keep track of the various distances needed to describe the weights, the following is a short description of how the intersections with the FV are calculated.

The problem is split in two; intersections with the top/bottom are derived separately from

intersections with the FV-edge. As the FV has a height zFV the top/bottom are situated at

z = ±1

2zFV and the edge is a distance rFV away from the center in the (x, y)-plane.

The particle position can be parametrized as

~r = ~r0+ s · ~p, (1)

where ~r is a three-component vector, ~r0a randomly picked position at the edge, top or bottom

of the detector, s > 0 interaction length which resembles time and ~p the normalized initial

direction.

Using this parametrization, calculating the intersections with the fiducial volume is straight-forward. For the edges one has to solve

(x0+ s · px)2+ (y0+ s · py)2 = r2FV (2)

and check whether the resulting z-coordinate lies within

−1

2zFV < z <

1

(18)

Similarly for intersections with the top (bottom) s is found by solving

z0+ s · pz = +(−)

1

2zFV (4)

and checking whether the corresponding (x, y)-coordinate satisfies

x2+ y2 < r_FV2 . (5)

If the initial position and direction of the particle are aimed towards the FV there will be two intersections. The (interaction) distance to the closest intersection will be called spre-FV and the distance between the two interaction points is sFV. As the particle scatters

in the fiducial volume, it changes direction and the distance towards the outside world is

recalculated. This distance will be called sout, see figure 6.

3.2.2 Reaching the FV

If λ(E) is the energy dependent attenuation length of a particle in xenon, the chance for the particle to react within a distance s is equal to

P (E) = 1 − e−s/λ(E). (6)

Thus, to force a particle to reach the FV it has to be given a weight

w1(E) = e−spre-FV/λ(E). (7)

3.2.3 Generating interaction point

Upon reaching the fiducial volume, the chance for a particle to scatter before leaving the FV again is simply

w2(E) = 1 − e−sFV/λ(E). (8)

To retain the correct distribution of travelled distances, one can not select a random point along the path. Instead, a random number r is uniformly drawn from the interval

h

0, 1 − e−sFV/λ(E)

i . Inverting the equation

r = 1 − e−sint/λ(E) ₍₉₎

gives

sint= −λ(E) log(1 − r), (10)

where sint is the interaction length to the scatter. The position can be found using equation

(19)

3.2.4 Generating an interaction

The weighting due to forcing the interaction is two-fold. First, the correct process must be

chosen and secondly the energy can not exceed a maximum deposit Emax.

If all processes are allowed, a process will be chosen based on the relative cross-sections and the corresponding weight is equal to 1 as no process was forced. However, sometimes there are processes that should not be simulated. For example, the photo-electric effect will result in

the photon depositing all his energy, so if the photon’s energy is higher than Emaxthis process

should be switched off. For neutrons, all but the elastic scatters are easily distinguishable from dark matter signal and thus are removed through data-cuts. Therefore there is no need to simulate anything but the elastic scatters for neutrons.

The corresponding weight is then

w3(E) =

P

iσi(E)

σtot(E)

, (11)

where σi(E) are the energy-dependent cross-sections of the allowed processes and σtot(E) the

summed cross-section of all processes.

After selecting the process, the energy deposition must not exceed Emax, resulting in a

weight w4(Ecut) = 1 σtot Z E0=Ecut E0₌₀ dσ(E0) dE0 dE 0_, ₍₁₂₎

where is the total integrated cross-section for the selected process and σ(E0) the energy

dependent cross-section for that particular process.

3.2.5 Leaving the detector

Similar to equation 7 the chance to leave the detector is

w5(E) = e−sout/λ(E), (13)

with sout the distance towards the outside of the detector.

The total weight is now simply the multiplication of the individual ones:

wtot =

i=5

Y

i=0

wi (14)

The benefit of this method is easily seen in figure 7 and 8 which show the location of single-energy deposits in the detector both by a conventional and a vrt simulation respectively. Both

simulations consist of 107 _{neutrons with an energy of 1 MeV. The left panels show the deposits}

before any data cuts are applied, the right panels show the same deposits post-data cut. The majority of simulated events in the conventional simulation lie outside the FV. The events that do show up within the FV are too few in number to create a smooth distribution resulting

(20)

in the freckled pattern on the right-hand side of figure 7. Do notice that the color-scale of the left panel is multiplied by a factor of 10. The distribution of events in figure 8 is totally different. No energy-deposits are removed by the data-cuts as all particles were forced to obey these cuts. This results in higher statistics and a smoother distribution. Even though the right-hand side panel in figure 8 includes more data points than its counterpart in figure 7, the total weight is equal. This is where the benefit of the vrt-simulation lies.

0 1000 2000 3000 4000

r

2

_[cm

2

_]

60 40 20 0 20 40 60

z [

cm

]

Pre-cut

0 50 100 150 200 250 300 350 400 0 1000 2000 3000 4000

r

2

_[cm

2

_]

60 40 20 0 20 40 60

z [

cm

]

Post-cut

0 5 10 15 20 25 30 35 40

Figure 7: Spatial distribution of lone energy deposits from a 107 _{1 MeV neutron non-vrt simulation.}

The colour scale indicates the number of deposits in the detector. The x-axis lists r2 _{to account for}

the increased circumference at larger radii. The left-hand side plot shows all deposits before any data-cut is applied while the right-hand side panel shows the data post-data-cut. The boundary of the fiducial volume is shown with a dashed black line.

0 1000 2000 3000 4000

r

2

_[cm

2

_]

60 40 20 0 20 40 60

z [

cm

]

Pre-cut

0 5 10 15 20 25 30 35 40 0 1000 2000 3000 4000

r

2

_[cm

2

_]

60 40 20 0 20 40 60

z [

cm

]

Post-cut

0 5 10 15 20 25 30 35 40

Figure 8: Same plot as figure 7 for a vrt-simulation. Although the left panel is pre-data cut, it is exactly the same plot as the post-data cut plot on the right. This is because all particles were forced to obey the data-cut resulting in a much smoother distribution than in figure 7.

(21)

3.3 Geant4

As Geant4 is also used in experiments in for example the LHC [13], the particles in this work, with energies ranging up to 1 MeV are considered low energetic particles [14]. Therefore, dedicated low energy-models were used. For photons LowEPComptonModel was used which is a dedicated model for photons below 20 MeV. It incorporates Doppler-broadening and takes into account that the electrons are bound to the nucleus. This gives a shift in energy compared to scattering of free electrons. The neutron were simulated using the NeutronHP-package (high precision) [14], [43]. Like LowEPComptonModel, this NeutronHP-package is valid up to energies of 20 MeV. The four types of reaction that are part of this package, elastic scattering, neutron capture, inelastic scattering and neutron induced fission will be covered in section 5.

(22)

4 Photon validation

As simulating photon-interactions is less complicated than neutron interactions it is not needed to use a dedicated simulation package like Geant4. A python-based model, created and validated by Kesber et al. successfully describes propagation of photons through a liquid xenon detector and the resulting energy deposits, using the same methods as described in section 3 [10]. The comparison between this work and the data of Kesber et al. had a twofold purpose. First, comparing the results of Kesber et al. with a standard Geant4 monte-carlo simulation served as another validation of the variance reduced method. Secondly, it served as an independent test of the Geant4 vrt code.

To compare the different simulations, section 4.1 will first discuss the different photon-interactions that have to be modelled. Section 4.2 will validate the implementation of these interactions. Consequently, section 4.3 discusses and validates the energy-cut.

4.1 Photon interactions

4.1.1 Compton-effect

If a photon scatters off a charged particle (in our case an electron) it can transfer some of its energy to the particle, this process is called the Compton-effect [44]. By considering a situation as shown in figure 9 it is straightforward to show that the scattered photon has an energy

E0 = E

1 + E

me(1 − cos θ)

, (15)

where E is the energy of the initial photon, me the electron mass and θ the angle between

the initial and final photon’s direction [44]. This argument only uses conservation of energy and momentum combined with the standard relation

(23)

Figure 9: Schematic visualisation of Compton-scattering. The incoming photon with energy E trans-fers some of its energy to the electron. Therefore the outgoing photon has an energy E0 _{< E which}

depends on the angle θ under which it scatters.

As can be seen from equation 15 E0 has a minimum for backward scattering (θ = π)

E_max0 = E

1 +2E_m

e

. (16)

This corresponds to a maximum energy deposit of Edep,max= E − Emax0 =

2E2

me+ 2E

(17)

4.1.2 Doppler broadening

Although figure 9 shows a photon that scatters of a static particle, the photons that scatter in the detector are not standing still, instead they are circling the nucleus. A hand-wavy argu-ment, that ignores most of the quantum nature of the electron and uses only non-relativistic dynamics, gives a rough estimation of their velocity.

Using the uncertainty relation one finds that

p · r ∼ ~ → mev ∼ ~

p.

Then setting the kinetic energy of the electron equal to its coulomb energy yields 1

2mev

2_∼ Ze2

4π0r

,

where Z is the atomic number of xenon. Combining these two equations finally gives

v ∼ 2Z e

2

(24)

with α = ₁₃₇1 the fine-structure constant. The right hand side of equation 18 is close to the

speed of light and even though it is a hand-wavy argument it shows the electrons are moving at considerable speed.

As a result, the photon will be either blue- or red-shifted depending on whether the electron is moving towards or away from the photon. This means a mono-energetic photon beam will have a broader spectrum than expected from equation 15.

4.1.3 Photo-electric effect

Instead of scattering of an electron, the photon can also be absorbed by the xenon-atom. All of the photon’s energy will be transferred to an electron which gains enough energy to be emitted as shown in figure 10. This process is called the photo-electric effect [45].

Figure 10: Schematic version of the photo-electric effect. The photon deposits all of its energy which is used to emit an electron from the electron cloud.

4.1.4 Pair production

For high photon-energies (MeV and higher) it is also possible for photons to decay into a particle and its anti-particle, e.g. an electron and a positron. However, as the process

γ → e++ e−

violates momentum-conservation, the process has to involve an additional particle. In xenon pair production will look like

γ + Xe → e++ e−+ Xe (19)

(25)

Figure 11: Pair production in which a photon decays into an electron and a positron. As this process would violate momentum-conservation another particle has to participate in the reaction, in this case a xenon-atom.

When the photon-energy increases, pair production becomes the dominant process. How-ever, at the energies simulated in this work (up to 1 MeV), pair production is a negligible effect [46]. It is therefore completely omitted from the simulation.

4.2 Validation photon-spectrum

Comparison between the vrt python code and the vrt Geant4 code is straightforward as they use the same principles. To compare these spectra with the standard monte carlo Geant4 code (from now on referenced to as the non-vrt code) is slightly more complicated. First of all the non-vrt code has a different way of tracking energy deposits. As opposed to the vrt-code, that only tracks the photon and considers all energy transferred from the photon as an energy deposit, it tracks the secondary particles. As a result, energy deposits are split into small energy transfers and tracked all the way down to zero. So a Compton-scatter of 600 keV would be split up in e.g. 40 deposits of (on average) 15 keV.

As the time between two scatters of the primary photon is about an order of magnitude larger than the time-scale in which the secondary particles lose their energy the hits can be clustered based on time. All energy deposits between two primary photon interactions are added to find the total energy deposit as is schematically shown in figure 12.

(26)

Figure 12: Clustering of the individual energy deposits from the non-vrt Geant4 code. All secondary deposits (blue) are added to the primary deposit (black). This clustering is needed to be able to compare the non-vrt data with the vrt-simulation.

After clustering the non-vrt data they are processed through the data-cuts as discussed in section 3.1, with the exception that no energy cut is made yet, to be able to compare with the vrt-simulation. The result is shown in figure 13. All three simulations started with 10 million 1 MeV photons.

0 200 400 600 800 1000

E

dep in fiducial [keV]

10 11 10 8 10 5 10 2 101 104

Counts / 10.0 keV

VRT

_{Geant4 MC}

Python vrt

Figure 13: Comparison of three different photon-simulations: a normal Geant4 simulation (yellow), a vrt Geant4 simulation (blue) and a vrt python simulation by Kesber et al. . All simulations consist of 10 million 1 MeV photons but only photons that scattered once in the detector are shown. The left part of the spectrum, up to roughly 800 keV is the compton-spectrum. The peak at 1 MeV is due to the photo-electric effect. The blue line extends to slightly higher energy deposits than the red line as the blue line includes doppler broadening.

Although the three spectra nicely overlap there are three points worth mentioning. 1. The non-vrt simulation does not seem to produce energy deposits between 750 and 850

(27)

accurate. Because of the weighting of the deposits they have much higher statistics and can have values below 1 per bin of the histogram. If the non-vrt simulation would have used ten times as many particles it would start to fill up this part of the spectrum after appropriate scaling.

2. The photo-electric effect fills two bins for the non-vrt simulation. This is due to the clustering of the events. In about 1% of the cases some energy is missed in the clustering and the total deposit ends up slightly lower.

3. The Compton part of the Geant4 vrt code stretches further than it does in the python code. This is due to Doppler broadening as discussed in section 4.1.2.

4.3 Energy cut

As stated in section 3.1.2 only low energy deposits are interesting for the XENON-project. Implementation of an energy cut for the photo-electric effect is trivial in the sense that photons with an energy above the threshold are not allowed to undergo the photo-electric effect. Achieving the same energy-cut for the Compton part of the spectrum is more complicated and will be the subject of this section. First section 4.3.1 will evaluate the procedure followed by Kesber et al. and why it did not give the correct distribution in this work. Section 4.3.2 will then discuss the final, successful, implementation of the energy cut.

4.3.1 Maximum angle

The initial idea was to follow the approach of Kesber et al. and analytically evaluate equation 12. As follows from equation 15 there is a one-to-one relation between the deposited energy and the scatter angle such that equation 12 can be changed to

w4(Ecut) = 2π σtot Z θ=θ(Ecut) θ=0 Z0dσ dΩsin θdθ. (20)

Here Z0 _{is the energy-dependent form-factor which will be discussed below and can be thought}

of as the effective number of electrons in a xenon-atom. _dΩdσ is given by the Klein-Nishina

formula differential cross-section per electron [47] dσ dΩ = 1 2r 2 ef (θ)2f(θ) + f(θ) −1 − sin2θ , (21) with f (θ) = 1 1 −_mE(θ) ec2(1 − cos θ) .

The scatter angle can consecutively be sampled from the cumulative density function. However, this method breaks down due to Doppler broadening, discussed in section 4.1.2,

(28)

as the spectrum gets elongated past the energy-cut value and thus a different approach is needed.

Form-factor

Even though equation 20 and 21 are not used to calculate the weight for the energy cut, they govern the sampling of Compton-energies in the simulation. Therefore a good

under-standing of these formulae is important for this work and Z0 from equation 20 deserves some

elaboration.

As equation 21 gives the differential cross-section per electron, it seems natural to multiply this term with the number of atoms in a xenon-atom (54). While this is approximately true for high energy photons, this overshoots the cross-section at lower energies. Due to a screening-effect a low energy photon will only see, and thus react with part of the electrons in an atom [48]. As the photon-energy increases the number of electrons per atom it reacts with gets progressively close to the atomic charge as can be seen in figure 14.

10 5 ₁₀ 4 ₁₀ 3 ₁₀ 2 ₁₀ 1 ₁₀0 ₁₀1 ₁₀2

E [MeV]

0 10 20 30 40 50

Form-factor

Figure 14: Energy-dependence of the form-factor for photons in liquid xenon. The form-factor can be thought of as the effective number of electrons in an atom. As the photon-energy increases the form-factor gets closer to the atomic charge.

4.3.2 Look-up table

The second attempt took a different approach. Instead of forcing all scatters to individually obey the energy-cut, energy-deposits were re-sampled until they fell below the energy thresh-old. In this way, the issue of Doppler-broadening is avoided as the re-sampling takes place after the Doppler-effect is calculated. On the other hand it is still more effective than a con-ventional Monte-Carlo simulation as the interaction point does not have to be recalculated. What is left is giving the particles a correct weight.

(29)

discuss the correct weighting. The shape of the spectrum in considerably different from 13 as figure 15 shows the unweighted spectrum. Photons that deposit a large amount of energy will have a lower energy afterwards. As the attenuation length increases with energy these photons are less likely to escape the detector and subsequently have a lower weight. For this reason the bump at the end of the spectrum has disappeared in figure 13.

0 100 200 300 400 500 600 700 800 900

E

100 101 102 103 104 105

Counts / 10.0 keV

VRT, unweighted

Figure 15: Unweighted Compton-spectrum of 1 MeV photons. The coloured area represents an energy cut. The weight for this cut should be equal to the coloured area divided by the total area below the graph. The graph does not continue all the way up to the photon energy as a photon can only deposit a certain amount of energy in one Compton-scatter.

The correct weight for forcing a 1 MeV photon to deposit less then 300 keV is equal to the area of the blue region divided by the total area beneath the curve. It is worth noting that this value has to be calculated only once. Thus, if one simulation is run to determine this factor it can be plugged into the simulation for all 1 MeV with an energy-cut of 300 keV. After the first scatter (or if one needs to simulate different energies), the photon energy will be different and a new value is needed. Therefore, multiple simulations were run with a photon energy spacing of 50 keV ranging from an initial energy of 200 keV up to 1 MeV. For each simulation the cumulative sum was taken with a binning of 2 keV which results in a table with the following structure:

(30)

Initial photon energy

200 keV 250 keV 300 keV ... 1 MeV

Energy dep osit 2 keV ... ... ... ... ... 4 keV ... ... ... ... ... 6 keV ... ... ... ... ... .. . ... ... ... ... ... 1 MeV ... ... ... ... ...

Table 2: Structure of the look-up table. Every column is an individual simulation. The table is filled with the fractional area of the spectra, meaning that the top left will show how often a 200 keV photon deposits an energy below 2 keV in a Compton-scatter in xenon. Note that cells with energy deposits higher than the photon energy, e.g. the bottom right, will have values of 1. 100% of all Compton scatters will have an energy below 1 MeV.

After filling this look-up table it can simply be loaded in at the start of a simulation. The correct weights are found by using a two-dimensional interpolation-method. The filled table is shown in figure 16 where the values are represented by colours.

200 400 600 800 1000 Ephoton[keV] 0 200 400 600 800 1000 Escatter [k eV] 0.2 0.4 0.6 0.8 1.0

Figure 16: Filled in look-up table for photons up to energies of 1 MeV. The values, indicated by the colour-scale, represent the weights corresponding to the different energy cuts. The layout of the table is shown in table 2.

The use of a look-up table yields the correct weights and has two advantages over the other method that was tried.

1. The table has to be generated only once. After the table is filled it only needs to be loaded in and no further calculations are required.

(31)

2. It is straightforward to apply this method for different scatter-processes and particles. One does not rely on an analytical approach as in section 4.3.1,

These two features enable the implementation of a similar scheme for the neutron spectra. Finally, to show that this method gives the correct weighting to the energy-deposits, figure 13 was recreated now including an energy cut. In addition to the three lines already present in figure 13 a green line has been added. This fourth line is a simulation of 10 million 1 MeV photons with the vrt-Geant4 code which includes an energy cut at 200 keV. The right plot of figure 17 shows that all four lines almost overlap, although the python vrt-code gives a slightly higher energy distribution. This is probably caused by Doppler-broadening and atomic relaxation, two processes that are not present in the python code, but should be investigated further.

0 200 400 600 800 1000

E

1011 108 105 102 101 104

Counts / 10.0 keV

VRT

VRTcut

Geant4 MC

Python vrt

0 25 50 75 100 125 150 175 200

E

101 100 101

Counts / 10.0 keV

VRT

VRTcut

Geant4 MC

Python vrt

Figure 17: Simulation of one millions 1 MeV photons. The plots include four lines, out of which three are the same as in figure 13. The fourth line, which is green, shows the vrt-Geant4 code where an energy-cut of 200 keV has been enforced. The plot on the right is a zoomed in version that focuses on the first 200 keV of the spectrum. Do note that the Python vrt-code consistently lies above the Geant4 vrt-codes, this slight deviation might be due to Doppler-broadening and atomic relaxation as these processes were added to the Geant4-code to make it more realistic. However, the difference between the lines should be investigated further.

(32)

5 Neutrons

Although photons form an important background in the XENON-detector, the main focus of this thesis lies on neutrons, which more closely resemble dark matter signals [35]. Section 5.1 will discuss the origin of the neutron background followed by a description of the different ways neutrons can react within the detector in 5.2. This also includes isotope-specific data in section 5.2.4. The final simulated spectrum is shown in 5.3 where it is compared to a conventional Monte-Carlo simulation.

5.1 Origin of neutrons

Neutrons that enter the detector can be split into two branches; cosmogenic neutrons and radiogenic neutrons.

5.1.1 Cosmogenic neutrons

Cosmological events like supernovae and gamma-ray bursts can accelerate particles well be-yond energies reached by man-made accelerators [49]. As mentioned before, the XENON-experiment is located within a mountain that shields the detector from this radiation. How-ever, muons, which can transverse large amounts of rock before they are stopped can still reach the detector.

Upon reacting, a muon will create a shower which can have neutrons with energies in the GeV-range [35]. To be able to remove this background signal, the water tank surrounding the detector is monitored with photo-multiplier tubes [2]. High-energetic muons travelling through the water tank will exceed the local speed of light and emit Cherenkov radiation [50]. This Cherenkov light is detected by the PMTs enabling the tagging of neutron-like signals as originating from a muon. This system, called the muon-veto, is able to tag neutrons originating from muons entering the water tank with an efficiency of >99.5%. If the muon-induced shower instead originates in the surrounding rock the neutrons can be tagged with a >70%-efficiency [51].

(33)

5.1.2 Radiogenic neutrons

Despite an extensive campaign in which radio-pure materials were selected to use in the de-tector, the dominant part of the neutron back-ground is due to radiogenic neutrons [2]. The water tank serves as a passive shield against neutrons originating from radioactive decay in the surrounding rock. Additionally, 120 photo-multiplier tubes are used to monitor the water tank making it an active neutron-shield. However, water shielding does not re-duce the neutrons coming from the detector material itself.

The neutrons are the result of trace amounts of

238_U, 235_U,226_Ra,232_{Th and}228_{Th in the}

dif-ferent parts of the detector [2]. These isotopes are unstable and can undergo spontaneous fis-sion. In the resulting decay chains multiple α-particles are emitted as is shown in figure

18, which is the decay chain of 235_{U. These}

α-particles result in emitted neutrons through (α,n)-reactions.

Figure 18: Decay chain of 235_{U taken from [52].}

The α-particles that are emitted along the chain result in neutrons through (α,n)-reactions.

Table 3 shows the abundance of the important isotopes in the different parts of the de-tector. Due to its relatively high mass the stainless steel cryostat contributes the majority of the neutrons, 36%. Despite its low mass, the photo-multiplier tubes are responsible for 33% of the neutron-background. The large contamination per kg is a result of the ceramic used in the PMTs which is relatively radio-loud [2].

(34)

Component Quantity Contamination [mBq/kg]

[kg] 238_U 235_U 226_Ra 232_Th 228_Th 60_Co 40_K 137_Cs

Cryostat vessels 1120 3.2(9) 0.37(13) 0.37(5) 0.29(7) 0.45(5) 2.5(5) 2.1(3) <0.41 Cryostat flanges 730 1.4(4) 0.06(2) <4 0.21(6) 4.5(6) 14.1(9) <5.6 <1.5 Bell and SS electrodes 190 3.2(7) 0.57(10) 0.62(10) 0.36(14) 0.46(9) 0.78(11) 1.6(6) <0.17

PTFE 128 0.12(5) <0.06 0.10(2) 0.11(5) <0.06 <0.053 2.4(3) <0.038

Cu parts 355 <0.69 <0.28 0.033(5) <0.027 <0.023 0.11(2) <0.29 <0.016

PMTs 94 50(15) 2.1(6) 2.6(10) 3.0(11) 2.4(9) 4.6(9) 72(19) 0.9(3)

Table 3: Table taken from [2] showing the contamination of detector materials by radio-active isotopes. The first column shows the different detector-components with their respective masses in column two. The abundance of the different isotopes is shown in mBq per kg.

5.2 Neutron interactions

As opposed to photons, which couple to electric charge, neutron cross-sections vary signifi-cantly for the different isotopes of xenon [14]. In combination with the higher energy range of neutrons in the detector, as discussed in section 5.1, none of the processes can be neglected. On the other hand, only neutrons that undergo elastic scatters, section 5.2.1, resemble dark matter signals. Neutron capture, inelastic scatters and neutron fission (section 5.2.2, and 5.2.3) can be identified by their reaction products and result in different signatures.

5.2.1 Elastic scatter

Just like wimps can (supposedly) scatter elastically from ordinary matter, so can neutrons. As both particles can scatter through the same process it is hard to distinguish neutron elastic scatter signals from a dark matter signal [35]. It is thus important to have a good understanding of the neutron background.

In elastic scatters no additional particles are created. It is solely a transfer from one particle to the other in which both kinetic energy and momentum are conserved. In principle, scatters in which particles get excited are not elastic scatters, however, as no additional particles are created, they are considered elastic within the Geant4-framework. As the rest mass of neutrons is ≈ 940 MeV they can be described non-relativistically in this work. This means (if no particles get excited)

1 2m1v 2 1 = 1 2m1v 2 2+ 1 2m2v 2 3 (22) and m1v~1 = m1v~2+ m2v~3 (23)

(35)

and the static particle respectively, v1 = | ~v1|2 the speed of the incoming particle pre-collision

and v2 (v3) the speed of the incoming(static) particle post-collision, see figure 19.

The energy deposit (1₂m2v32) will depend on the angle under which the particles

scat-ter, where the maximal deposited energy corresponds to backwards scattering. In this case, equation 23 reduces to

m1v1 = −m1v2+ m2v3. (24)

Solving this equation for v2 and plugging it into equation 22 gives

v2 = m2 m1 v3− v1 v1 = 1 2 m2 m1 + 1 v3 and thus Edep, max= 1 2m1v 2 1 4m1m2 (m1+ m2)2 . (25)

As a xenon-atom is about 131 times heavier than a neutron (depending on the isotope) a non-relativistic neutron can deposit up to 3% of his energy in a single elastic scatter. The full energy-spectrum can be derived by keeping the angular dependence and yields [53]

Edep= Ein 4m1m2 (m1+ m2)2 cos2_θ, ₍₂₆₎

where Ein is the energy of the incoming particle and θ chosen as in figure 19.

Figure 19: Schematic representation of the elastic scattering between a neutron and a xenon-nucleus.

If particles are excited, the energy-distribution also depends on the shell-structure of the xenon, which gives rise to spikes in the cross-section as shown in section 5.2.4.

(36)

5.2.2 Neutron capture

It is also possible that a neutron gets absorbed by the xenon-nucleus instead of scattering away. After absorption, the xenon-nucleus is in an excited state and needs to get rid of its energy surplus.

If the additional energy is released by emitting a photon the process is called neutron capture. Neutron capture raises the mass of a xenon nucleus by one atomic unit as in equation 27.

n +xXe →x+1Xe + γ (27)

Figure 20: Schematic representation of neutron capture by a xenon-nucleus. The excited nucleus relaxes by emitting a photon.

For incoming neutron energies above a few tens of keVs, the neutron-capture cross-section tends to decrease for the different xenon-isotopes, as shown in section 5.2.4. Like for elastic-scatters, the cross-section depends on the shell-structure of the atoms.

5.2.3 Inelastic scatter and neutron induced fission

If the incoming neutron has sufficient energy, the excited nucleus can instead split into multiple fragments as shown in figure 21. Within Geant4 these reactions are treated either as inelastic scatters or neutron induced fission depending on the nature of the post-reaction particles.

Possible inelastic reactions include the emission of an alpha particle leaving a tellurium core:

n +xXe →x−3Te + α. (28)

(37)

For low incident neutron energies, inelastic scatters are subdominant to elastic scatters. However, as the neutron energy increases, inelastic scatters become progressively more im-portant.

For higher energies, X and Y can have comparable mass in which case the reaction is considered a fission-reaction. Due to the high energy-threshold neutron-induced fission is irrelevant for this work. Section 5.2.4 shows the cross-section for neutron-induced fission is zero in the investigated energy range.

Figure 21: Schematic representation of neutron induced fission. The excited xenon-atom releases its energy by decaying into multiple fragments. Depending on the nature of the reaction-products the event is considered either an inelastic scatter (emission of light particles like a helium or tritium core) or a neutron-induced fission (X and Y are similar in size).

5.2.4 Isotope dependent cross-section

As neutrons scatter within the strong-nuclear potential of the xenon nuclei, it is important to distinguish the different isotopes within the detector. The Geant4 package NeutronHP, which stands for neutron high precision, includes cross-sections for different isotopes and is specifically aimed towards low energetic neutrons [14]. Figure 22 shows the cross-sections for elastic neutron scattering, inelastic scattering, capture and induced fission for both natural xenon and the natural occurring isotopes individually. The legend shown in the first graph applies to all plots.

The validity range of the NeutronHP package goes up to 20 MeV. For these energies there is no neutron induced fission in any of the xenon-isotopes.

Xe-132

(39)

103 ₁₀2 ₁₀1 ₁₀0 ₁₀1

E [MeV]

1034 1033 1032 1031 1030 1029 1028 1027 1026

[m

2

]

Xe-134

103 ₁₀2 ₁₀1 ₁₀0 ₁₀1

E [MeV]

1032 1031 1030 1029 1028 1027

[m

2

]

Xe-136

Figure 22: Neutron cross-sections for natural xenon and its different isotopes. Cross-sections were calculated using the Geant4 NeutronHP package up to its validity limit of 20 MeV [14]. All plots share the same colour-code where elastic scatter cross-sections are shown in blue, inelastic in yellow and neutron capture in red. As neutron-fission is not possible for any of the xenon-isotopes for these energies, the cross-section is equal to zero and does not show up in the graphs.

5.3 Validation (neutrons)

As elastic neutron scatters most closely resemble the expected wimp-signal, the vrt Monte-Carlo simulation forces the incoming neutron to scatter elastically.

Furthermore, the cross-sections shown in section 5.2.4 enable the calculation of the atten-uation length λ(E) = 1 σ(E)n = A σ(E)ρNA , (29)

where σ(E) is the total cross-section for neutrons in natural xenon, n is the number density of

xenon atoms, A the atomic mass of natural xenon, ρ the mass density and NA= 6.022 · 1023

g mol-1_{Avogadro’s number. This leads to an attenuation length for 1 MeV neutrons of 10.8}

cm which is in good correspondence with values quoted by [35].

As the python vrt-code only enables the simulation of photons it can not be used to validate the neutron vrt-code. Instead, the neutron simulation will be tested against a conventional non variance reduced Monte-Carlo code. Like in section 4.2 only particles that scattered once in the fiducial volume are taken for comparison with the additional constraint that the scatter was elastic. Figure 23 shows both the vrt and the non-vrt neutron simulations. Both simulations consist of 10 million 1 MeV neutrons. The good correspondence between the vrt (blue) and non-vrt (yellow) code serves to validate the variance reduced neutron simulation.

(40)

[keV]

Figure 23: Single energy deposits by a vrt and a non-vrt neutron simulation. Both simulations consist of 10 million 1 MeV neutrons. From the non-vrt code only those neutrons that scattered exactly once in the detector and where that one scatter was both elastic and within the fiducial volume are shown.

(41)

6 Acceleration

Figure 23 shows that the vrt-neutron code produces the correct neutron spectrum. It also shows a smoother curve for the vrt-code. This section will quantify this effect using the same approach as Kesber et al. as explained in 6.1 [10]. Consecutively, section 6.2 will show the amount of acceleration achieved by using the vrt-neutron code.

6.1 Relative error

As the method is called variance reduced Monte-Carlo simulation instead of accelerated Monte-Carlo simulation the term acceleration might seem misplaced. However, accelera-tion refers to the reducaccelera-tion of the number of simulated particles needed with the vrt-code to reach the same variance as with a standard Monte-Carlo.

Rather than comparing the shape of the energy spectra it is more convenient to compare the area below the curves. This area corresponds to the number of background events in the detector and is an important quantity for the experiment. Therefore, the acceleration will be the reduction of the number of particles that has to be simulated with a vrt-code to reach the same error on the expected post-data cut number of background-events as with a conventional Monte-Carlo simulation.

To calculate the error on expected events it is useful to define a truth value. This truth value is derived using the fact that the chance P that an individual simulated particle survives the data-cuts is independent of the number of simulated particles. Once this chance is derived, one can simply calculate the truth value using

Ntruth= P · Npart. (30)

Here, Ntruth is the expected number of particles that survives the data-cuts and Npart is the

number of simulated particles.

As shown by [10], the vrt-simulation converges faster to the truth value than standard Monte-Carlo simulations and it is therefore most convenient to use a large reference vrt-simulation to calculate P . As can be easily seen from equation 30, P can be calculated with

P = Aref

Nref

, (31)

where Aref is the area beneath the deposited-energy curve and Nref the number of particles

simulated in the reference simulation. The value of P depends heavily on many parameters (e.g. FV-size, neutron-energy). With the detector-dimensions as listed in section 3 we find

P = 3.2 · 10−3 _{for 1 MeV neutrons. This means less than one percent of particles in a}

(42)

As expected, figure 24 shows the vrt-simulation converges rapidly towards Ntruth while

the error-bars on the non-vrt simulation are considerably larger. However, at Npart = 106

there is a discrepancy between the vrt and non-vrt code which seems to indicate that the vrt-simulation overestimates the number of background effects in the detector. First investigations in this issue seem to indicate a small mismatch between the attenuation lengths of both simulations. Further work should pinpoint the problem and resolve it.

truth

Figure 24: Deviation from Ntruth both by a vrt and a non-vrt code. Each point represents a total of

10 simulations. As Ntruth is calculated using a large vrt-simulation it is no surprise that the vrt-code

converges towards this value. In correspondence to what was found by Kesber et al. the error on the vrt-simulation is smaller than on the non-vrt version [10]. The difference between the vrt and non-vrt code at Npart= 106seems to indicate the vrt-code overestimates the number of events in the detector.

First results from investigating this discrepancy seems to hint to a mismatch between the attenuation lengths of the simulations.

To calculate the size of the error-bars in figure 24 we conducted a family of 10 simulations. Simulations within a family have the same number of simulated particles, the same energy and only differ in the random seed which initialised the simulation.

A natural choice for the error seems to be the square root of the variance σ = s P i(Ai− µ)2 Nsim , (32)

where i runs over all simulations within the family, µ is the average number of data-cut

surviving events and Nsim = 10 the number of simulations in a family. However, as σ depends

on the number of particles in the simulation (Npart) it is more useful to use the relative error

RE = σ

µ (33)

(43)

6.2 Acceleration result

To quantify the efficiency of the vrt-simulation figure 25 shows the value of RE for 4 different

simulation-families. The line through the points is fitted by assuming RE depends on Npart

like

RE(Npart) =

q

pNpart

, (34)

where q is the fit-parameter.

103 ₁₀4 ₁₀5 ₁₀6

N

part 10 2 10 1 100

Relative error

VRT

Geant4 MC

Figure 25: Relative error (RE) of different families on simulations. The vrt-simulation has a lower error for the same number of simulated particles (Npart). To reach the same RE as a conventional

Monte-Carlo simulation, a vrt-code needs to simulate 12 times less particles as indicated by the black horizontal line. The fitted curves are shaped like equation 34.

Figure 25 shows the vrt-simulation needs less particles to get to the same relative error. The acceleration α is the horizontal distance between the two fitted lines. This can be calculated by

α = qnon-vrt

qvrt

2

= 12 (35)

The value of 12 is two orders of magnitude smaller than the value quoted by [10] for photons with the same energy. However, this was to be expected as neutrons have a larger penetration depth than photons. Furthermore, this work did not include an energy-cut for neutrons as the neutron-energies considered have energy-deposits well below the dark-matter energy-threshold of 200 keV. For photons of 1 MeV the energy-cut is responsible for a large part of the acceleration so the fact that this was not needed for neutrons of similar energy explains the difference in acceleration achieved for the different particles. Implementation of an energy-cut for higher energetic neutrons will yield higher values of α.

Variance reduced monte-carlo simulation for neutrons in XENON

MSc Physics and Astronomy

GRAPPA

Master Thesis

Variance reduced monte-carlo simulation for neutrons in

XENON

by

Gijs Leguijt

11000279

June 2020

60EC

September 2019 – June 2020

Supervisor/Examiner:

Examiner:

Prof. dr. A.P. Colijn

Dr. M. Vreeswijk

Summary

Contents

1

Introduction

2

The XENON-experiment

3

Monte-Carlo simulations

r

[cm

]

z [

cm

]

Pre-cut

r

[cm

]

z [

cm

]

Post-cut

r

[cm

]

z [

cm

]

Pre-cut

r

[cm

]

z [

cm

]

Post-cut

4

Photon validation

E

Counts / 10.0 keV

VRT

Geant4 MC

Python vrt

E [MeV]

Form-factor

E

Counts / 10.0 keV

VRT, unweighted

E

Counts / 10.0 keV

VRT

VRTcut

Geant4 MC

Python vrt

E

Counts / 10.0 keV

VRT

VRTcut

Geant4 MC

Python vrt

5

Neutrons

E [MeV]

[m

_[cm

_]

_[cm

_]

_[cm

_]

_[cm

_]

_{Geant4 MC}