mini-TimeCube as a Neutron Scatter Camera

TL;DR

This paper demonstrates through Monte Carlo simulations that a compact detector, originally designed for antineutrino detection, can be used as an effective neutron scatter camera with precise spatiotemporal resolution for imaging fast neutrons in the 1-10 MeV range.

Contribution

The study introduces a novel application of the miniTimeCube detector as a neutron imaging device using waveform digitization and Bayesian reconstruction methods.

Findings

Waveform digitization enables ~100 ps timing and ~5 mm spatial resolution.

Bayesian framework effectively reconstructs neutron direction and energy.

Populations of reconstructed neutrons can identify source location and spectrum.

Abstract

We present Monte Carlo (MC) simulation results from a study of a compact plastic-scintillator detector suitable for imaging fast neutrons in the 1 -- 10 MeV energy range: the miniTimeCube (mTC). Originally designed for antineutrino detection, the mTC consists of 24 MultiChannel Plate (MCP) photodetectors surrounding a 13 cm cube of boron-doped plastic scintillator. Our simulation results show that waveform digitization of 1536 optically sensitive channels surrounding the scintillator should allow for spatiotemporal determination of individual neutron-proton scatters in the detector volume to $\thicksim$100 picoseconds and $\thicksim$5 mm. A Bayesian estimation framework is presented for multiple-scatter reconstruction, and is used to estimate the incoming direction and energy of simulated individual neutrons. Finally, we show how populations of reconstructed neutrons can be used to estimate the direction and energy spectrum of nearby simulated neutron sources.

Tables

Table 1. Table 1: r 𝑟 r reflections in a rectangular volume can be modeled by n 𝑛 n virtual points outside the detector.

Reflections	New Points	Total Points
$r$	$p=4r^2+2$	$n=\sum _i^r{p}_i$
0	1	1
1	6	7
2	18	25
3	38	63
4	66	129
5	102	231

Table 2. Table 2: mTC neutron candidate criteria, including time, distance, energy and photon thresholds, and energy ratio thresholds. Statistics shown for 1, 3, and 5 MeV neutrons as well as Pu and background (Bk) spectra. Background spectra is composed of cosmic-ray induced neutrons at sea level Goldhagen2004 .

Criteria	Bk	Pu	1 MeV	3 MeV	5 MeV
$>20$ PE	.19	.38	.43	.33	.27
$>10P_1$ PE	.16	.33	.36	.30	.25
$>5P_2$ PE	.19	.38	.43	.32	.27
$t_2-t_1>1$ns	.14	.31	.38	.22	.15
$P_2-P_1>10$mm	.18	.36	.42	.30	.24
$d{E}_0^1>100$keV	.14	.28	.31	.27	.23
$d{E}_0^1/d{E}_1^2>.20$	.15	.28	.32	.22	.19
	.18	.36	.40	.32	.26
$d{E}_0^1/E_0>.1$	.12	.26	.33	.19	.13
Combined	.07	.16	.22	.12	.07

Table 3. Table 3: P 1 subscript 𝑃 1 P_{1} , P 2 subscript 𝑃 2 P_{2} , t 1 subscript 𝑡 1 t_{1} , t 2 subscript 𝑡 2 t_{2} and d ​ E 0 1 𝑑 superscript subscript 𝐸 0 1 dE_{0}^{1} mTC MC 1 ​ σ 1 𝜎 1\sigma resolutions for 1, 3 and 5 MeV neutrons as well as Pu and background (Bk) spectra. Background spectra is composed of cosmic-ray induced neutrons at sea level Goldhagen2004 .

$E_0$	$P_1$	$P_2$	$t_1$	$t_2$	$d{E}_0^1$
	(mm)	(mm)	(ns)	(ns)	(fraction)
1 MeV	5.4	19	0.13	4.3	0.35
3 MeV	3.2	15	0.04	2.3	0.23
5 MeV	2.8	16	0.03	2.1	0.19
Pu	2.8	16	0.03	2.1	0.19
Bk	2.8	16	0.03	2.1	0.19

Table 4. Table 4: Θ Θ \Theta and E 0 subscript 𝐸 0 E_{0} mTC MC 1 ​ σ 1 𝜎 1\sigma resolutions and efficiencies for 1, 3 and 5 MeV neutrons as well as Pu and background (Bk) spectra. Background spectra is composed of cosmic-ray induced neutrons at sea level Goldhagen2004 .

$E_0$	PEs	efficiency	$\Theta$	$E_0$
		(fraction)	(deg)	(fraction)
1 MeV	76	0.23	37	0.32
3 MeV	272	0.12	26	0.26
5 MeV	430	0.08	24	0.24
Pu	170	0.12	26	0.26
Bk	-	0.08	24	0.24

Equations

$$\sin \Theta =\sqrt{\frac{E_0-E_1}{E_0}}=\sqrt{\frac{d{E}_0^1}{E_0}}$$

$$E_0=d{E}_0^1+E_1$$

$$E_1=\frac{1}{2}m{v}^2$$

$$v=\frac{dx}{dt}=\frac{||P_2-P_1||}{t_2-t_1}$$

$$\displaystyle\theta_{n}=\left[\begin{array}[]{ccc}P_{1}&t_{1}&w_{1}\\ &\vdots&\\ P_{n}&t_{n}&w_{n}\end{array}\right]$$

$$\sum _i^n{w}_i=1$$

$$\displaystyle z=\left[\begin{array}[]{cc}P_{1}^{\gamma}&t_{1}^{\gamma}\\ \vdots&\vdots\\ P_{m}^{\gamma}&t_{m}^{\gamma}\end{array}\right]$$

$$p(\theta |z)=\frac{p(z|\theta )p(\theta )}{p(z)}=\frac{1}{c}p(z|\theta )p(\theta )$$

$$p(\theta |z)=\prod _j\sum _i{w}_{i}p\left(z_j|{\theta }_i\right)p\left({\theta }_i\right)$$

$$p\left(z_j|{\theta }_i\right)={\Lambda }_t{P}_{\Omega }{P}_{\gamma }{P}_{T}Q$$

$$\Lambda_{t}\left(x\right)=\left\{\begin{tabular}[]{ll}$\left(1-e^{-x/t_{r}}\right)e^{-x/t_{f}}$&$\quad x\geq 0$\\ 0&$\quad x<0$\end{tabular}\right.\\$$

$$x=t-({\theta }_t+r/v)$$

$$t_r$$

$$\Omega =2\pi \left(1-\cos \theta \right)$$

$$P_{\Omega }=2\pi \left(1-\frac{r}{\sqrt{r^2+a^2}}\right)\hat{r}\cdot \hat{n}$$

$$f\left(x;\lambda \right)=1-e^{-\lambda x}$$

$$P_{\lambda }=1-f\left(x;\lambda \right)=e^{-\lambda x}$$

$${\theta }_i={\theta }_r$$

$${\theta }_t=\arcsin \left(\frac{n_1}{n_2}\sin {\theta }_i\right)$$

$$R_s$$

$$=$$

$$R_p$$

$$=$$

$$R=\frac{R_s+R_p}{2}$$

$$T=1-R$$

$$p(z|\theta )=\prod _{i}f(k_i;{\lambda }_i)$$

$$f(k;\lambda )=\frac{{\lambda }^k{e}^{-\lambda }}{k!}$$

$${\lambda }_i=E\ast Y\ast p(z_i|\theta )$$

$${\lambda }_i=E\ast Y\ast P_{\Omega }\ast P_{\gamma }\ast P_T\ast {\Lambda }_t\ast Q$$

$$f(k;\lambda )=e^{k\log (\lambda )-\lambda -\log \left(\Gamma (k+1)\right)}$$

$$\Gamma (k+1)=k!$$

Full text

simulation, estimation, optimization, fitting, neutron, gamma, detection, direction, energy, scintillation

††preprint: AIP/123-QED

1. I Introduction

2. I.1 Reactor antineutrinos

3. I.2 Neutron detection

4. II Detector Description

5. III Neutron Imaging Theory

6. IV Neutron Imaging in Practice

7. V Observing Neutrons in the mTC Detector

8. V.1 Individual Neutron Estimation

9. V.2 Mixed MAP Estimator

10. V.2.1 Scintillator Temporal Response

11. V.2.2 Solid Angle

12. V.2.3 Attenuation

13. V.2.4 Reflection

14. V.3 Poisson Energy Estimator

15. V.4 Bayesian Information Criterion (BIC)

16. V.5 Candidate Efficiency

17. V.6 MC Results

18. VI Results of Neutron Imaging Experiments

19. VII Conclusion

20. VII.1 Future detectors

21. VIII Acknowledgments

22. IX Bibliography

I Introduction

Recent advances in methods to reconstruct particle interactions in scintillating mediajgl09 coupled with advanced picosecond (ps) resolution photodetectors and electronics ber13 will allow for new instruments with unprecedented four-dimensional (x, y, z, t) imaging of antineutrinos and neutrons emanating from nuclear reactors and radionuclides, respectively. For antineutrinos this results in improved spatiotemporal resolution of the positron and neutron associated with the Inverse Beta Decay (IBD) reaction. While for neutrons, it allows for observation of “light starved” proton-proton scatters. Such a methodology allows for the reconstruction of particle kinematics (vertex and energy deposition) and subsequent estimation of the particle topology and directionality. This capability has significant scientific applications, in such areas as the study of neutrino oscillations phenomena and the search for sterile neutrinos.

The primary purpose of this paper is to address the issue of fast neutron directional detection in plastic scintillator instrumented with fast electronics and high photocoverage. However, in this introduction we also provide information of why it is important to build multi-purpose detectors, which can be used to monitor nuclear reactor via antineutrino inverse-beta-decay detection and to detect special nuclear materials via fast-neutron re-scattering.

I.1 Reactor antineutrinos

There are significant practical applications of such developments, primarily in the imaging of antineutrinos to geolocate nuclear reactors and detecting neutron scatters to find Special Nuclear Materials (SNR) from increased stand-off distances. Additionally, as part of the development program, study of spatiotemporal resolution will improve algorithmic discrimination of cosmogenic and geophysical backgrounds (which usually plague prior types of neutrino and neutron detectors).

Thus, we present a description of the miniTimeCube (mTC), a proof of concept antineutrino and neutron imager. The detector consists of a $(13\ \mathrm{cm})^{3}$ Boron-doped plastic cube surrounded with four detectors (Photonis 64 pixel multichannel plate photomultipliers) per face, operating in the 850 ps rise-time range ($\thicksim 60$ picoseconds single-PhotoElectron (PE) time resolution with waveform fitting). The detector was demonstrated at the NIST Center for Neutron Research (NCNR) in Gaithersburg, Maryland.

For antineutrinos, the majority of reactor based neutrino experiments have traditionally employed the use of large scale water/scintillator based detectors. The concept of using small neutrino detectors ($<5$ kt) near nuclear reactors to monitor nuclear power and fuel burn up in the core was introduced in 1970s, but was technically impractical then. Subsequently, in Russia the Rovno neutrino spectrometer for the first time observed the distortion of the antineutrino spectrum as a function of uranium burnup and plutonium accumulation within the core kli94 .

The KamLAND kiloton scintillation detector started operation in 2002 and recorded neutrinos from reactors all around Japan, definitively determining the oscillations of electron anti-neutrinos. The similarly-sized Canadian experiment SNO simultaneously settled the solar neutrino anomaly (missing neutrinos) in favor of electron neutrino oscillations; and the smaller Borexino experiment in Italy clinched the issue. In more recent years detectors at 100 m to 2 km from reactors, as at Double CHOOZ, Reno and Daya Bay have been employed to constrain oscillation parameters at various baselines dchoo ; reno ; daya , completing our first order understanding of three-neutrino mixing.

The San Onofre Nuclear Generating Station detector (SONGS1) collaboration proposed ber02 and demonstrated bow07 a compact (about one tonne) antineutrino detector specifically designed for International Atomic Energy Agency (IAEA) nuclear safeguards applications. They demonstrated the possibility of monitoring a (cooperating) nuclear reactor in neutrino output from tens of meter ranges, employing fairly simple modern technology. Note that the special virtue of such observations is that such are completely non-invasive, and that neutrinos cannot be faked or hidden. Several types of detectors have been tested at the now abandoned SONGs station but only one yielded useful results, emphasizing that such development is not completely trivial even with modern equipment.

Recently, a reanalysis of all reactor neutrino experiments since the original Reines and Cowan discovery of neutrinos in the 1950’s, has shown a $\thicksim 6\%$ deficit in the expected antineutrino flux at baselines $<100$ m compared to calculated fluxes, which is referred to as the Reactor Antineutrino Anomaly (RAA) men11 . The RAA coupled with a renewed interest in counter-proliferation have generated an increased focus on small-scale neutrino detectors.

Several small-scale detectors have been built to address non-proliferation, short baseline oscillations, and the recent evidence of the reactor antineutrino anomaly. Given the compact size of the detectors, they are typically deployed near the reactor core, typically within 100 meters. Experiments at such short baselines experience extremely high reactor-induced backgrounds that are otherwise negligible at longer baselines.

For long distances, there is another experiment, currently under construction — WATer CHerenkov Monitor of ANtineutrinos (WATCHMAN Askins:2015bmb ). It will be the first prototype detector to remotely monitor individual reactor operations at a significant distance.

The three main backgrounds of concern are gamma radiation, fast neutrons, and cosmogenic muons. Lead or other heavy metal, and polyurethane or wax shielding is often used to mitigate reactor-induced gamma and neutron radiation. In addition to passive shielding, muon paddles (large flat scintillation counters) are typically used to tag muons. Muons may pose a problem for electronic dead time for larger detectors ($>1$ m3) which are not deeply buried. But muons are also problematic to antineutrino detection even in deep detectors due to the occasional production of long lived isotopes whose decay a second later can fake a neutrino interaction.

I.2 Neutron detection

For neutron imaging, eleven orders of magnitude separate the IBD cross section and the neutron-proton cross section mechanism. Thus, indefinite volumetric expansion of neutron detectors is constrained due to the diminishing returns in rate measurements, since the neutrons will not deeply penetrate the detector volume. The first directional neutron scatter camera measured the double proton recoil in two separate planes of plastic scintillator to determine a crude estimate of the angular resolution of the neutron van07 . Using this method, first presented in 1985 pre85 ; wal96 , neutron energy can be determined through a combination of peak amplitude measured in the primary plane followed by a time-of-flight measurement in the second plane. Leveraging the estimated locations of proton scatters in each plane in conjunction with the energy lost from the pulse amplitude in the 1st plane produces a cone of directionality whose angle is determined by rudimentary kinematic equations for elastic scattering. Similarly, a neutron scatter camera using two spatially separated liquid scintillator cells demonstrated the capability of detecting a 252Cf source at a 30m stand-off distance mas09 . Pinhole apertures and an array of coded pinhole apertures have also successfully imaged fissile neutrons bla11 ; bru09 . However, the pinhole aperture approach to directionality are significantly limited due to the low count rates associated with achieving directionality through restriction of the observing solid angle area of the aperture. Our approach differs from the preceding primarily because we use a single volume to record both scatters.

Similar to cloud chambers or bubble chambers, Time Projection Chambers (TPCs) employ a detector volume filled with a gas or liquid or combination of gases allow for the high resolution imaging of particle scatters. Bubble and cloud chambers permit the photography of the particle induced disturbances. TPCs employ an electric field to drift the ions and electrons to electronic pickups, making a sort of time sequence of the image projected onto a plane. Our concept with the mTC aims at using the now available high precision light detection time to conceptually combine these methods, in a sort of photon-TPC. Key to this method, now being explored for use in higher energy detectors as well as here at nuclear energies, is the realization of Fermat’s Principle, which states that light will always take the shortest path from object to image. In other words, despite the relatively long emission times of scintillating materials, the leading photons will reveal the event topology. In the past we have mainly thought of scintillation detectors as only useful for calorimetry, measuring the total light output and thus energy since the scintillation radiation is isotropic. It has been demonstrated that fast light detection can pinpoint interaction locations (the vertex) and even permit track reconstruction and particle identification Barrett:2013zka .

As applied to neutron detection, such a methodology is advantageous compared to the aforementioned neutron cameras because it allows for the three-dimensional reconstruction and ionization of the recoiling proton. Helium-3 has typically been considered the “gold standard” for neutron detection and TPCs in general due to the high neutron detection efficiency, gamma radiation discrimination kou09 , non-toxicity, and low cost gao11 . However, the recent shortage in global helium supplies she10 has resulted in an increased level of research activity to innovate neutron detectors to mitigate the lack of 3He. Several successful 2H TPC alternatives smi05 ; bow10 ; hun08 have been demonstrated however they have not been widely accepted into the safeguards community due to safety concerns. Gas filled TPCs provide an excellent means of imaging fast neutrons, however they suffer from a low detector cross section due to the $\sim 10^{-3}$ decreased density of gas as compared to plastic. Such detectors also do not promise to allow for robust operation in a field deployment (wire vibration, high voltage, ultra-cleanliness, and high pressure special gases and containers are problematic).

The mTC aims to overcome these shortcomings. Key to this is also the availability of relatively low cost, high density and very fast (ps) digitization electronics, developed at the University of Hawaii. While we originally developed the concept for antineutrino imagingli2016 , the mTC could also potentially image fast neutrons from radionuclides with performance comparable to traditional gas-based TPCs. We leverage a robust GEANT4 and MATLAB Monte Carlo simulation, validated by real world data jocher2013 , which shows the mTC’s ability to image fast neutrons with a high efficiency and significant angular accuracy.

II Detector Description

The mTC is a $13\mathrm{cm}\times 13\mathrm{cm}\times 13\mathrm{cm}$ cube, composed of 2.2 liters of 1% natural boron-doped (0.2% 10B by weight) EJ-254EJ-254 plastic scintillator surrounded by 24 Planacon Micro-Channel-Plate Photomultiplier Tubes (MCP-PMTs). Each MCP-PMT has 64 unique channels, for 256 channels per face and 1536 total across all 6 sides. Figure 1 shows the mTC with surrounding MCP-PMTs and all associated electronics. Each channel is coupled to full-waveform readout electronics sampling at $\sim$2.8 GigaSamples Per Second (GSPS), which are theoretically capable of providing $<$ 100 ps timing resolution of single Photo Electrons (PEs).

During the mTC operation, we investigated different trigger schemes Li:2018dwm . The final version was a 2-level trigger — with preset threshold level (L0) and minimum and maximum numbers of triggered channels (L1 min and L1 max).

III Neutron Imaging Theory

When a neutron in motion is deflected from a straight line trajectory, the angle of the deflection may be calculated given the kinetic energy of the neutron prior to and immediately following the deflection. This simple relationship is expressed by Equation 1:

[TABLE]

where $\Theta$ is the neutron deflection angle, $E_{0}$ is the original neutron energy and $E_{1}$ the neutron energy after deflection; $dE_{0}^{1}$ being the difference. $dE_{0}^{1}$ is directly observable, however $E_{0}$ is not, and must be inferred via Equation 2:

[TABLE]

where the posteriori kinetic energy $E_{1}$ is calculate via the neutron velocity following deflection:

[TABLE]

where

[TABLE]

To solve for $v$, one needs to know $P_{1}$ and $P_{2}$ (shown in Figure 2), the positions of the first and second scatters, and $t_{1}$ and $t_{2}$, the associated scatters times (visible in Figure 3). Neutron mass $m=939.565378$ MeV/c2 and speed of light $c=299.792458$ mm/ns are assumed. Substituting 2, 3 and 4 into 1 allows us to directly solve for $E_{0}$ and $\Theta$ from the available observables.

IV Neutron Imaging in Practice

To successfully apply the neutron direction and energy method described in this paper, 3 observable parameters are required from scatter 1:

Position $P_{1}$ 2. 2.

Time $t_{1}$ 3. 3.

Neutron kinetic energy lost during scatter 1: $dE_{0}^{1}$

and 2 additional parameters from scatter 2:

Position $P_{2}$ 2. 2.

Time $t_{2}$

In practice none of these parameters are perfectly knowable; however, we may attempt to estimate their values from the evidence the neutron leaves behind. This evidence comes in the form of photons, light emitted by the charged particles that neutrons tend to bump into in organic scintillator. Light sensors surrounding the scintillator can record the timing and location of these photon bursts, and this information can in turn be used to estimate the 5 parameters above.

Organic scintillator is composed of long hydrocarbon chains with loosely bonded protons at the perimeter. Neutrons traveling through these chains will tend to hit the loosely bonded protons about 60% of the time, colliding with 12C atoms the rest of the time (A small minority ($\sim$1%) of scatters occur on C13 and other isotopes).

When a neutron hits a proton in the scintillator, it transfers some of its kinetic energy to the proton, which in turn becomes visible light via charged particle ionization of the surrounding atoms the proton passes by. The proton path is usually short, a few mm, and the scintillation yield is usually low, due to the large mass of the proton, yet a high energy ($\sim$ 1 MeV) neutron recoil may still yield significant observable light.

V Observing Neutrons in the mTC Detector

The real world observability, or ‘resolution’, of each of these 5 parameters ultimately dictates neutron $\Theta$ and $E_{0}$ observability, therefore we seek to estimate these values as accurately as possible for every incoming neutron. To this purpose we have created estimators and verified their performance via MC simulation of a real world detector: the UH miniTimeCube (mTC) antineutrino detector. Though intended primarily for antineutrino research, the fast timing of this detector, long recording buffer, and high neutron cross section make it well suited to neutron reconstruction as well, thus we use it as our modeled testbed for simulated neutron reconstruction.

V.1 Individual Neutron Estimation

To fit a single neutron, we apply 5 parallel estimators to the measurements, with each estimator assuming a different number of observable neutron scatters within the detector. The Maximum A Posteriori (MAP) likelihood and parameter count of each estimator are then used to compute the corresponding Bayesian Information Criterion (BIC), and the smallest BIC of the 5 indicates the most likely number of scatters. If 2 or more scatters are determined by the BIC, then the first two scatters are isolated and used to compute neutron $\Theta$ and $E_{0}$ per Equation 1. If the BIC determines only 1 scatter occurred, the event is not considered a candidate and is rejected.

Each estimator fits the MAP values of the 5 parameters ($x$, $y$, $z$, $t$, $E$) which fully describe each recoil. For $n$ recoils the parameters are:

[TABLE]

where $w$ indicates the weight of each scatter. For a finite set of scatters, the non-negative, nonzero weights $w>0$ are constrained to a convex combination summing to unity:

[TABLE]

The measurements available for this task, $z$, are the observed photon arrival positions $P^{\gamma}$ and times $t^{\gamma}$ on the MCPs surrounding the scintillator:

[TABLE]

The measured positions are discretized to the MCP anode centers (each channel is a 5.9 mm $\times$ 5.9 mm square). Measurement time resolutions are determined through MC study, in the case of the mTC we find these to be $<$ 100 ps $1\sigma$ per channel for first arrivals including the effect of the MCP Transit Time Spread (TTS) (we note that empirical Photonis MCP differ substantially from simulated results).

The following steps are used to fit up to 5 possible neutron scatters:

Assume 1 scatter and estimate $\hat{\theta}_{1}^{MAP}$ 2. 2.

Assume 2 scatters and estimate $\hat{\theta}_{1...2}^{MAP}$ 3. 3.

Assume 3 scatters and estimate $\hat{\theta}_{1...3}^{MAP}$ 4. 4.

Assume 4 scatters and estimate $\hat{\theta}_{1...4}^{MAP}$ 5. 5.

Assume 5 scatters and estimate $\hat{\theta}_{1...5}^{MAP}$ 6. 6.

Evaluate Equation 38 BIC for each $\hat{\theta}$ 7. 7.

Scatter count defined by minimum BIC 8. 8.

Assign individual photon probabilities to each scatter 9. 9.

Estimate energy of each scatter ($dE$) using ‘extended Poisson’ Equation 9

V.2 Mixed MAP Estimator

The MAP cost function follows Bayes’ equationBar-Shalom :

[TABLE]

where the likelihood $p(z|\theta)$ forms a mixture distribution with 1 mixture component per neutron recoil. Dropping the normalizing constant $c$:

[TABLE]

where the likelihood $p(z_{j}|\theta_{i})$ of measurement $j$ given source $i$ with weight $w_{i}$ and prior $p\left(\theta_{i}\right)$ is simply an evaluation of the measurement space created by $\theta_{i}$ at $z_{j}$. This measurement space is defined by a point-source position $P_{\theta}$ at time $t_{\theta}$, and is a function of several detector and scintillator characteristics including:

• •

Scintillation spectrum, yield and decay constant(s)

• •

Cherenkov spectrum

• •

Quenching factors for heavy particles

• •

Scintillator attenuation length

• •

Re-emission efficiency of attenuated photons

• •

Refraction indices of the scintillator and PMT glass

• •

PMT quantum efficiency, size and location

• •

Time and energy calibrations

The likelihood of observing a single photo-electron (PE) $z_{j}$ at a single point due to a single point-source $\theta_{i}$ is

[TABLE]

(visualized in Figure 4) where $\Lambda_{t}$ is the temporal likelihood discussed in V.2.1, $P_{\Omega}$ is the solid angle probability discussed in V.2.2, $P_{\gamma}$ is the un-attenuated energy probability discussed in V.2.3, $P_{\mathrm{T}}$ is the transmission probability discussed in V.2.4, and $Q$ is the PMT quantum efficiency.

V.2.1 Scintillator Temporal Response

A scintillator’s temporal response is typically defined by a double exponential distribution. The probability density function (PDF) for the temporal likelihood that we use is shown in Figure 5, and is:

[TABLE]

where

[TABLE]

Here $t$ is the time of the first photon arrival in channel $i$, $\theta_{t}$ is the estimated time of interaction, r is the distance to the MCP channel $i$, v is the velocity of light in the medium, and $t_{f}$ and $t_{r}$ are the scintillator fall time and rise time respectively. The function acts as a prediction of measurement times for photons that propagate from a scintillation vertex.

In the mTC’s EJ-254 scintillator these times are:

[TABLE]

When used in equation 16, the $x$ value in equation 17 is the time difference between the observation $z_{j}$ and the point-source $\theta_{i}$, after compensating for the travel time between the two points. The travel time compensation assumes a scintillation and QE weighted mean photon group speed. This will be the average speed of all observed scintillation PEs. In mTC this speed is 185.8 mm/ns.

V.2.2 Solid Angle

The solid angle of a cone with half-angle $\theta$ is

[TABLE]

from which we can approximate the solid angle of an arbitrarily shaped surface $S$ (the MCP channel centered on $P_{z}$ in our case) with area $A$ as viewed from point $P_{\theta}$ by

[TABLE]

where $\vec{r}=\theta_{\Omega}-P_{z}$ is the vector pointing from $\theta_{\Omega}$ to PMT location $P_{z}$, $\hat{n}$ is the unit normal vector on the surface at the point $P_{z}$, $a$ is the radius of the surface area subtended by the cone, and $r$ is the norm of $\vec{r}$. If a certain PMT (or a group of channels) on the mTC wall observes a large percentage of the emitted photons, maximizing the solid angle probability will involve moving our point-source estimate closer to that area. The approximation comes from substituting the area of a square pixel for the radius-squared of a circular pixel.

An illustration of the mTC’s solid angles from all pixels is shown in Figure 6, and in general for a single source in Figure 7.

V.2.3 Attenuation

The fraction of energy attenuation at distance $x$ from a source is given by the exponential Cumulative Distribution Function (CDF)

[TABLE]

where $\lambda$ equals the inverse of the attenuation length of the medium being traversed. Conversely, the survival fraction is

[TABLE]

V.2.4 Reflection

When a ray of light travels between media with different indices of refraction $n_{1}$ and $n_{2}$, some incident energy will transmit through to the second media while the rest will reflect back into the first. This ratio of transmission and reflection is governed by the Fresnel equations, discovered in the early 1800’s by Augustin-Jean Fresnel, a French engineer and physicist. The law of reflection holds that the reflection angle $\theta_{r}$ will be equal to the incident angle $\theta_{i}$

[TABLE]

and Snell’s law defines the transmission angle $\theta_{t}$ of the refracted light as

[TABLE]

For $s$ polarized light (incident light is polarized with its electric field perpendicular to the incident plane. The light is said to be s-polarized, from the German senkrecht (perpendicular)), the reflection coefficient (which gives the fraction of reflected light), is given by

[TABLE]

and for $p$ polarized light (incident light is polarized with its electric field parallel to the plane containing the incident, reflected, and refracted rays. Such light is described as p-polarized) the reflection coefficient is

[TABLE]

If the incoming light is unpolarized (with equal parts $s$ and $p$ polarization), then the reflection coefficient, shown in Figure 8, is

[TABLE]

Via conservation of energy the transmission fraction is

[TABLE]

In cubical scintillation volumes like the mTC, reflected light can make up a significant part of the overall light collected. The overall light distribution is the sum of the original source plus the reflection source-points. To simulate $r$ reflections, $p$ reflected origins are required, as given by table 1.

V.3 Poisson Energy Estimator

The Poisson Energy Estimator is akin to the traditional photon-counting method, the main difference being that the estimator is capable of assuming multiple point sources spread out in space and assigning expected photon counts to each pixel based on these multiple point sources.

The Poisson energy estimator estimates point-source energies conditional on a given point-source location. It can operate on a single source at a time, or estimate the combined energies of all sources jointly, which may then be apportioned to each source using its likelihood ratio. The estimator maximizes the likelihood of a product of Poisson probabilities, one for each pixel in the detector:

[TABLE]

The Poisson probability mass function (PMF) is

[TABLE]

where $\lambda$ is the expected number of events and $k$ the measured number. An integration through $\Lambda_{t}$ in Eq. 16 (along with scintillator photon yield $Y$ and PMT quantum efficiency $Q$) is sufficient to supply us with an expected number of photons at each channel over a time period

[TABLE]

or,

[TABLE]

The observed numbers of photons $k$, are unfortunately not integers but rather fractions, as are the expected number of photons $\lambda$. A Poisson distribution is only defined for integer values of $k$, thus the problem calls for a different type of distribution capable of non-integer support. This distribution, which we dub the ‘extended Poisson’ probability is shown in Figure 9, and takes on the form

[TABLE]

where the gamma function $\Gamma$ escapes confinement by interpolating the factorial

[TABLE]

between integers as

[TABLE]

V.4 Bayesian Information Criterion (BIC)

The Bayesian Information CriterionBIC (BIC), and its close cousin the Akiake Information CriterionAIC (AIC), can both be used to evaluate competing models on the same measurements. In this application, we find that the BIC is more effective than the AIC.

The BIC and AIC penalize larger numbers of parameters fitting the same number of measurements, and thus prevent over-fitting. They become necessary only when the competing models have differing numbers of parameters, precluding a simple likelihood comparison. A caveat is that since they measure relative performance, not absolute performance, they are unable to warn if all the models fit poorly.

The BIC equation employed is:

[TABLE]

where $p\left(\theta|z\right)$ is evaluated at the ML estimate of $\theta$, and $n_{\theta}$ and $n_{z}$ are the number of parameters and measurements respectively.

The BIC is used exclusively in our work to determine the number of observed recoils a neutron experiences within the mTC, and to assist in distinguishing gammas from neutrons as will be explained in the MC results section. Even though we are only interested in the first two scatters, we consider it a prerequisite to first identify all the scatters before isolating the first two from the rest.

V.5 Candidate Efficiency

Not all neutrons entering the mTC are viable candidates for estimation. Due to the small size of the mTC, many neutrons leave the detector volume before capture, thus depriving the observer of the characteristic delayed bookend which can assist in positive neutron identification. While not specifically required for neutron direction and energy estimates, an observed capture can identify a particle as a neutron with greater confidence than relying solely on the scatter kinematics.

A neutron may stay within the mTC volume yet not capture within the available recording time (12 $\mu$s for the mTC). In other cases, a neutron may be identified correctly, yet may simply not deposit adequate energy for reasonable estimation. To exclude these low-quality neutrons from our measurement pool, certain ‘candidate criteria’ must be established.

As a result of these criteria, only about 6% of 1 MeV neutrons become viable candidates in the mTC, significantly reducing the available measurement population. Table 2 includes a complete list of candidate criteria applied to our simulated neutrons and the resultant candidate efficiencies (the fraction that pass).

Figure 10 shows 1 MeV neutron capture times and locations in the 1% EJ-254 mTC. The mTC electronics are limited to 12 $\mu$s of recording time. Note the number of neutrons captured past the mTC time and space limits. Only 70% of neutrons are captured within $12\mu$s in a 1% EJ-254 mTC, and only 12% of neutrons stay inside, and capture within the mTC.

V.6 MC Results

We use a Monte Carlo simulation of the mTC built with GEANT4GEANT , MATLABMATLAB and Python to determine the observability of each of the 5 required reconstruction parameters. Our results are reported as $1\sigma$ resolutions in Table 3, with the underlying distributions shown in Figure 11.

Table 4 shows $\Theta$ and $E_{0}$ resolutions, along with neutron candidate efficiency over the energies tested, and the number PEs captured in the mTC per neutron event.

Figure 12 shows a histogram of the simulated velocities of the neutron between their first two scatters, and gives us a simple cut with respect to what is a reasonable velocity to expect. We use this to set a velocity rejection cut at $\leq 40$ mm/ns, (note that the speed of light in the scintillator is approximately $190$ mm/ns). This is gives us an extra buffer against accidental fits for double gamma recoils.

A simulation of 50000 gammas with energies of 0-3 MeV and incident upon the top face of the cube yielded only 3 successful fits using the neutron scatter fitter. Those 3 events were then promptly rejected using our secondary criteria, that their velocities were on the order of the speed of light in the medium.

A qualitative explanation for why gammas are rejected in our model is given in Figure 13. The example gamma Compton scatters twice creating two outgoing scintillation spheres. As light takes approximately a nanosecond to traverse the mTC volume, the scatters are separated in time by only a fraction of a nanosecond. The gamma itself is traveling faster than the photon sphere is expanding and our model will likely fit for a single recoil at a weighted average between the two scintillation points. For the neutron recoils the average time between scatters is approximately 2 nanoseconds. The slower moving neutron allows for the early light from the first recoil to reach the PMTs before the first light produced from the second recoil. As shown in figure 3, the light arriving from both neutron recoils will always be mixed due to the scintillator decay time.

The fact that we can have confidence that our model only fits for neutrons is an excellent result and allows us to make bold determinations about neutron directionality. Another avenue for determining neutron directionality is to only consider the average directional vector that the set of vectors $(P_{2}-P_{1})$. In Figure 14 we see a comparison of the fit vs truth for the average directional vector of a set of 10093 neutron events and the agreement is encouraging. $\cos(\theta)$ is .998 +/-.002 corresponding to an angle difference of (-.1,.13) radians or (-5.7,7.4) degrees at a 99.96 percent confidence limit.

Figure 15 shows composite 4pi PDFs over multiple neutron reconstructions. The PDFs for multiple neutrons are simply the summation of individual neutron PDFs. These results include only ‘signal’ neutrons from a point source of interest. A more comprehensive study would include background neutrons from the GoldhagenGoldhagen2004 sea level spectrum, as well as background gammas, which are generally present at two orders of magnitude above the neutron background. The addition of both of these backgrounds to the simulation would provide a more accurate picture of detector performance under non-ideal conditions.

VI Results of Neutron Imaging Experiments

During the neutrino data-collection campaign at NIST the mTC experienced an overheating incident, which resulted in the loss of more than half of the PMTs. The extent of the structural damage to the remaining PMTs (and other components) was unknown. The simple act of turning on the system and performing simple start-up calibrations was nontrivial. The mTC was shipped back to the University of Hawaii at Manoa for diagnosis. Secondary effects were observed including the increased trigger rate and higher noise. Neutron Imaging Tests were conducted in July 2017 using a californium source placed on various positions around the cube and the results were inconclusive.

Tests were preformed on a surviving PMT to analyze its performance. Laser pulses were injected directly into single MCP channels and a large amount of cross-talk was observed in the surrounding pixels. Figure 16 quantifies the timing differences between the target pulse (the pixel in which we inject laser) and the undesired pulses generated in the adjacent pixels. The average delay was $1.76$ ns with a standard deviation of .70 ns is a serious issue for our neutron imaging purposes. Since the average time between the low-light generating first and second recoils is on the order of 2 ns, our ability to track neutrons is not possible with the cross-talk response we observed in our PMTs (unknown to us at the time of detector construction).

VII Conclusion

We presented an efficient model for reconstructing neutron recoils in the mTC which could be adapted for any related design using a single solid volume scintillating target. Our simulations showed that in reconstructing the average directional vector for neutrons entering our scintillating target, the angular agreement between the truth and fit was $\cos(\theta)$ is .998 +/-.002 corresponding to an angle difference of (-.1,.13) radians or (-5.7,7.4) degrees at a 99.96 percent confidence limit.

Unfortunately, before we undertook neutron imaging tests with a source, our detector suffered a catastrophic overheating event that caused it to behave unreliably. Additionally, the PMTs we used did not behave as expected in terms of their timing response. Nonetheless, we are confident that our model is accurate and efficient for neutron recoil reconstruction.

VII.1 Future detectors

The miniTimeCube experiment led to the realization that the single scintillation volume might more ideally be replaced by an array of segmented scintillation volumes to allow for more accurate isolation and reconstruction of individual neutron recoils in future double-scatter detectors. Some of these detector concepts are currently under exploration at the University of Hawaii and various U.S. national labs.

VIII Acknowledgments

We would like to thank the University of Hawaii and Ultralytics LLC for their support in this effort. Lawrence Livermore National Laboratory is operated by Lawrence Livermore National Security, LLC, for the U.S. Department of Energy, National Nuclear Security Administration under Contract DE-AC52-07NA27344. LLNL-JRNL-762419.

IX Bibliography