DEVELOPMENT AND EXPERIMENTAL VALIDATION OF RESPONSE MODELLING FOR TIME-OF-FLIGHT NEUTRON DETECTION AND IMAGING SYSTEMS

The application and feasibility of a time-of-flight neutron detection system is explored for sources with time correlated gamma-ray and neutron emissions, such as the spontaneous fission emitter, Cf-252. For the emission of multiple gamma rays and neutrons from a single spontaneous fission event, a near instantaneous gamma-ray detection followed by a later neutron detection on a multi-detector array allows for an associated time-offlight to be determined for a neutron arising from that event. Using a suite of purpose developed analysis tools, Monte-Carlo simulation and experimental data are compared for the Cf-252 water tank source facility at Lancaster University. Applying a bespoke timeof-flight imaging algorithm, vector-based optimisation (VBO), the true source location is determined within 21 cm by this approach.


INTRODUCTION
The implementation of a time of flight (ToF) approach to fast neutron detection offers a novel alternative to traditional, indirect fast neutron detection systems. The time-of flight approach relies on the detection of correlated gamma and neutron emissions, for example, ones originating from a unique fission event. In the case of Cf-252 the average number of emissions following a spontaneous fission event is 3.757 [1] for neutrons (ν) and 8.32 for gammas [2]. The time difference measured between the detection of gamma and neutron events, here assumed to be in spatially separated detector elements sensitive to both gammas and neutrons, can be used to provide coupled information relating to the energy and source location of the detected neutron.
Nuclear applications of the ToF imaging technique is a rapidly growing area of research interest. Prior to technological developments which have enabled this growth, detection of fast neutrons relied on moderation and the resulting detection of thermalised neutrons, suitable only for counting applications. An alternative approach involves the direct detection of fast neutrons using scintillation detectors, which have the advantage of being responsive to both fast neutrons and gammas. The development of EJ-309, a non-flammable, non-carcinogenic xylene based liquid scintillant, has allowed for the pioneering use of scintilation detectors for fast neutron detection without the potential harmful consequences of early variants. Newer plastic scintillants, such as EJ-299 and its variants, are also promising and provide comparable pulse shape discrimination (PSD) to the EJ-309 liquid scintillant. Coupling this to the parallel development of fast digitizer technology, commercially available systems are now capable of performing on board pulse shape discrimination [3] between gammas and neutrons, whilst also achieving resolving times on nanosecond scales.
Imaging nuclear material by ToF methods has a wide range of potential applications in industry and medicine [4]. Fissile material exists in many forms, including nuclear reactor fuel, plutonium stockpiles, metals, oxides, fluorides and other forms, all of which are particularly abundant in the UK. Radiation monitoring of these materials is of significant national and global importance while civil nuclear power and weapons programs are ongoing. Another potential application may arise in the imaging of shipping containers which may contain illicitly transported fissile material; this would be a significant step forward in the technologies available to security and customs officials. In medicine, hadron therapy used for cancer treatment is on the ascendancy; ToF imaging could provide an alternative method of monitoring neutron dose induced in a patient or staff by proton beam interactions, allowing for the optimisation of beam line intensities and energies.
In this paper, the feasibility of using ToF has been explored using a simulation approach; a bespoke suite of analysis tools, 'pytrac', has been developed to interpret particle track data from Monte-Carlo simulations, facilitating the visualisation of correlated gamma-neutron emissions, perform detector response characterisation, and generate neutron source emission maps using time-of-flight data. Experimental data collected from the Lancaster University water moderated Cf-252 source facility has been utilised to tailor validation of the modelling approach and algorithm development. For time-of-flight analysis, correlated timing information is required in order to determine prompt gamma and neutron arrival times on an event by event basis. The standard MCNP method involves sampling a Watt probability distribution, invoking a new source history for each particle emission from the source. In order to model particle histories on a fission-by-fission basis, the Lawrence Livermore National Laboratory fission model is used [5]. The algorithm used provides a more realistic representation of each individual fission, generating multiple prompt gammas and neutrons and their correlations in emission, and in some cases, subsequent neutron induced fissions occurring later in time. Other physics card options included are energy cut offs for photons, neutrons and protons. Appropriate reduction of minimum energy cut off levels for each of these particles allows the tracking of low energy interactions at a greater computational expense. This will also add greater resolution to the generated pulse profile specifically for gamma and neutron interactions in the detector elements, allowing for a more accurate determination of pulse shape parameters such as pulse rise time. A final inclusion in the model allows for light ion recoil induced by elastic scattering of neutrons; this is the primary mechanism by which protons are generated in scintillator detector elements.

Geometry and materials
Detailed modelling of the detector elements has been carried out using specifications presented in [6]. The elements are each 10 x 10 x 10 cm, with a 3 mm aluminium casing added around the edge of the liquid scintillator element, containing EJ-309 as specified from [7]. Additional detail has been added to the detector's external structures; this includes the addition of the photomultiplier tube and magnetic shield. An X1 source capsule model has also been used. Labelled diagrams of the source and detector models are shown in Figure 1. Both inner tanks have a thickness of 10 mm, and are contained within one another, in a "Chinese doll" type configuration. Detailed CAD data has been provided for the water tanks, which has been simplified for subsequent conversion into CSG geometry for use in MCNP, using the SuperMC software [8]. Detailed models of the EJ-309 detectors have been integrated into the water tank model, positioned to replicate experiments described in Section 3.2. A 3D representation of the experimental set up is shown below in Figure 2.

PYTRAC development
In MCNP6.1, physical quantities such as flux, fluence and energy deposition are extracted from a simulation using user defined tally definitions, set up in the input deck. Whilst it is possible to record the energy deposition of gamma and neutrons in MCNP with standard tallies, as well as the time profile of these depositions in discrete time bins, one cannot isolate true coincidences in this way. An alternative method to extract data from MCNP is by using the particle-track (PTRAC) output file. This contains detailed event-by-event data, specifying quantities such as position, direction, interaction type, and energy deposition. This data can then be interrogated as a postprocess; not only does this enable the user to extract data of interest after the simulation (rather than defining tallies as a pre-requisite), correlated pairs of gammas and neutrons can be identified.
The PYTRAC module is a suite of python based tools developed for a number of time-of-flight applications. At the core of the toolkit is a PTRAC file parser, which has been set up to interrogate the PTRAC output file and identify gamma-neutron coincident events which deposit energy in separate detector elements. Prompt gammas and neutrons induced by fissions are neutral particles and therefore do not directly produce signal carriers. Only the energy deposition of induced secondary charged particles are registered; in the case of gammas, we see electrons being produced by photoelectric absorption, Compton scattering and pair production, whereas neutrons produce recoil protons in elastic collisions with the scintillator material, typically hydrogen in organic scintillators. In order to discriminate against low energy, partial depositions induced by scattered photons or neutrons, a discriminator level at 500 keV has been introduced. This is consistent with experi-mental acquisition systems like those used to acquire the data in Section 3.2.
A visualisation feature has been built into PYTRAC. Event histories of interest, for example, cases where signal carrier energy depositions induced by correlated photon and neutrons both exceed 500 keV in separate detector elements, can be plotted in order to observe characteristic features. This is useful for diagnostic purposes, and also predicting optimal locations to position detectors in arrays. Figure 3 shows an example of the event visualiser applied to a simple test model scenario, containing a central sphere of Cf-252 and two EJ-309 detector elements either side.

Imaging Algorithm Development
Source imaging is a fundamental requirement of the time-of-flight technique for it to be a tangible application in fields such as medical imaging or portal monitoring. A bespoke imaging algorithm has been tailored for use on both experimental and simulation data. For the experiments described in Section 3.2, pulse timing electronics have been used to record the time-of-flight between gamma and neutron signals, determined by on board real time pulse shape discrimination. As a result, only time-of-flight data can be obtained from this setup, so imaging algorithms reliant on pulse height information cannot be applied to this data. We therefore focus solely on the time-of-flight only variant of the vector optimisation method, described below in Section 2.2.1.

Vector-based least squares mesh optimisation
Vector based optimisation (VBO) on a global imaging volume involves iterative optimisation of a single fitting parameter,v, with respect to a measured or simulated result, by minimising the least squared difference between these two quantities. In the case of time-of-flight, a fundamental assumption is made that the mode time-of-flight,t measured by a single detector element is in proportion to its distance from the true source location, i.e. d =vt.
An assumption is made that by taking a mode value for all detectors, these effects should cancel each other given adequate statistics. The time-of-flight optimisation for a single imaging mesh element is expressed below in equation 1, where χ 2 is the quantity to be minimised by manipulation ofv, N is the number of detectors in the array, d n is the distance between detector n and the centre of the mesh element being evaluated, andt n is the mode time-of-flight determined for detector n.
At first instance, the global image volume is split into a 2 x 2 x 2 mesh, with the optimisation applied to each element in the mesh, using an inbuilt numerical solver within the Python SciPy package. The element with the lowest χ 2 is subsequently selected, and split into a 2 x 2 x 2 mesh; the process is then repeated for a user defined number of loops. Following this optimisation, a best fit parameter,v, minimised χ 2 value, and most likely source location is returned. Whilst this information is useful for a quick estimation of the source position, it does not constitute an image as the imaging mesh becomes non-uniform during the iteration procedure. We therefore redefine a uniform, fine resolution image mesh, calculating χ 2 , and corresponding p-value for each element fine mesh using the optimisedv parameter determined by iteration. In this application the p-value, derived for statistical hypothesis testing, evaluates the statistical probability that the null hypothesis has been satisfied, i.e. the probability that equation 1 is a valid model for a given mesh element. P-value arrays are therefore plotted for image visualisation in Section 3. A visualisation of the method can be seen in Figure 4, showing the convergence from an array of neutron detectors (shown in orange).

Mode time-of-flight determination
A requirement of the VBO approach in Section 2.2.1 is that a mode time-of-flight must be known.
In experimental practice, a distribution of time-of-flight values will be measured in discrete histogram bins. For consistency, discrete time-of-flight data generated by simulation, and processed by PYTRAC, is binned into the same histogram time bins. The mode could be approximated by identifying the maximum bin value and using the associated time-of-flight, however this would be limited by the time bin resolution, and statistical fluctuation in the binned data. A more rigorous solution involves applying a fitting function, extracting the mode from parameter estimates. A background corrected Gaussian has been determined as the most promising candidate function, shown in equation 2, where A is the pulse height term, σ is the variance term, and C is the  Figure 4: Convergence of iterative VBO to determine the optimal best fit parameter for the system,v.
background term. A separate fit should be performed for each detector, extracting the mode timeof-fight,t n from each fit. An example fit, based on simulation data described in Section 3.2, is shown in Figure 5.

RESULTS AND APPLICATION
In order to validate the methods outlined in Section 2.2, the VBO method has been applied to two hypothesis cases, shown in Section 3.1, and later to experimental and simulated data sets, shown in Section 3.2.

Hypothesis Cases
The results of two hypothesis cases are shown in Figure 6; an idealised case, where mode timeof-flight values have been exactly generated to ensure equation 1 = 0 for all detectors, and also a random case, where random mode time-of-flights have been generated for each detector. For valid conclusions to be drawn from the VBO method, one should expect p-values very close to 1 at the true location of the source in the idealised case, and close to 0 in the random case.

Simulation vs Experiment
Experiments were performed at the water moderated Cf-252 source facility at Lancaster University, using EJ-309 liquid scintillators and a 4 channel mixed field analyser (MFA). A single detector was used as a gamma trigger, which remained in the same position for all measurements. An energy threshold was set to 500 keV for all measurements, and a pulse shape discrimination factor was set such that neutron attribution was conservative. TTL output signals were processed by an FPGA-SoC board, which returns time-of-flight histograms between gamma and neutron signals in 5 ns bins. A location map for all detector positions is shown in Figure 2, where the fixed gamma detector is indicated in orange. Observation times varied from 88 minutes to 202 minutes, depending on the incident count rate at respective positions. Figure 7 shows VBO images generated for both of these data sets.

Statistical Analysis
Two simple statistical checks have been applied to both hypothesis cases, and the experimental and simulation data. As the true source position is known in all of these cases, the first test calculates  Comparing results for the random and idealised data, as expected, we can assume we do not see any true statistical convergence with a p-value of 2.48 × 10 −8 . However, for the idealised data, the calculated source position is within 0.006 cm, with a p-value of 1.0. These findings suggest the outcome of VBO on both hypothesis cases yields expected results, and thus we can apply VBO to real data sets with some confidence. In the simulation and experimental results tested here, both data sets have an image location accuracy, ∆, of 21 cm when rounded to the nearest integer, with p-values of 0.975 and 0.996 for experimental and simulation data respectively. With a 0.95 confidence limit typically set in null hypothesis testing, we can assume that in both cases, we see true statistical convergence.

CONCLUSIONS
In light of technological developments in nuclear instrumentation, facilitating real time time-offlight detection systems, a suite of analysis tools, PYTRAC, has been developed to interpret particle track data from Monte-Carlo simulations in MCNP. Using PYTAC to extract time-of-flight data, the response of a detection array can be predicted and subsequently optimised by iteration. An experimental validation of the MCNP PYTRAC approach has been undertaken by modelling the Cf-252 source facility at Lancaster University. Our bespoke imaging algorithm, vector-based optimisation (VBO), has been applied to both experimental and simulated data; when rounded to the nearest integer, the true source location is imaged to within 21 cm.
Vector-based optimisation appears very promising for future time-of-flight imaging systems, and work is ongoing to integrate the VBO methodology into real-time instrumentation, however there are a number of limitations in the approach which would benefit from further research. A limiting factor of the experimental configurations assessed in this paper is that pulse height information, which could be unfolded to resolve neutron energy deposition within the detector, is unavailable. Another assumption that has been made throughout is that the measured (or calculated) time-of-flight is always in direct proportion to the distance between the true source and detector. Whilst both of these limitations would have little effect for mono-energetic, un-moderated neutron source imaging applications, facilities such as the Cf-252 source at Lancaster University are water moderated, thus the average neutron velocity will depend on the detector position and how much moderation it has undergone along its trajectory. In summary, improving VBO to take account of non-uniform moderation effects, and interpreting pulse height data, should reduce the uncertainty in source position measurements, and therefore improve the accuracy and versatility of the approach in general.