# Internal mechanical dissipation mechanisms in amorphous silicon

Carl Lévesque<EMAIL_ADDRESS>Sjoerd Roorda
<EMAIL_ADDRESS>François Schiettekatte
<EMAIL_ADDRESS>Normand Mousseau
<EMAIL_ADDRESS>Département de Physique and Regroupement
québécois sur les matériaux de pointe, Université de Montréal, C.P. 6128,
succursale centre-ville, Montréal, Québec, Canada H3C 3J7

###### Abstract

Using the Activation-Relaxation Technique-nouveau, we search for two-level
systems (TLSs) in models of amorphous silicon (_a_ -Si). The TLSs are
mechanisms related to internal mechanical dissipation and represent the main
source of noise in the most sensitive frequency range of the largest
gravitational wave detectors as well as one of the main sources of decoherence
in many quantum computers. We show that in _a_ -Si, the majority of the TLSs
of interest fall into two main categories: bond-defect hopping where neighbors
exchange a topological defect and the Wooten-Winer-Weaire bond exchange. The
distribution of these categories depends heavily on the preparation schedule
of the _a_ -Si. We use our results to compute the mechanical loss in amorphous
silicon, leading to a loss angle of 10-3 at room temperature, decreasing to
10-4 at 150 K in some configurations. Our modeling results indicate that
multiple classes of events can cause experimentally-relevant TLSs in
disordered materials and, therefore, multiple attenuation strategies might be
needed to reduce their impact.

††preprint: APS/123-QED

##  I Introduction

Current gravitational wave detectors (GWD) consist of Michelson
interferometers with arms the length of a few km containing a Fabry-Perot
cavity. Since 2015 [1], GWDs have successfully detected nearly a hundred
events[2], at an accelerating pace thanks to continuous efforts. Among the
targets for improvement are the test masses: massive dielectric mirrors at the
end of each arm. Their reflective surface consist of a stack of alternating
high refractive (HR) and low refractive (LR) index materials. In the current
implementation of LIGO and VIRGO, the LR material for coating is amorphous
silica and the HR material Ti-doped amorphous tantala [3]. Amorphous
materials, especially in the HR layers, present intrinsic fluctuations that
can be directly linked to the mechanical dissipation ($Q^{-1}$) and thermal
noise through the fluctuation-dissipation theorem [4, 5, 6]. Such phenomena
also cause decoherence in some quantum computers [7]. Despite considerable
efforts to reduce these losses in the HR layers [8, 9, 10], low mechanical
loss remains the limiting factor of noise at frequencies around 50 Hz in major
GWDs, at which these detectors are the most sensitive [11].

Here, we investigate the origin of the mechanisms leading to internal
mechanical dissipation through atomistic simulations. In order to simplify the
problem, we investigate amorphous silicon (_a_ -Si), a prototypical model of a
continuous random network. _A_ -Si consists of a single element and its
structure and dynamics have been under investigation for more than 50 years
[12]. Beyond its generic interest, _a_ -Si is directly relevant in the GWD
context as it is considered for future generations of GWD [13]: it features a
high refractive index, reducing the number of layers in the stack, and it can
be synthesized with a ultra-low internal mechanical dissipation [14].

We show that (i) a distribution of two-level systems within the
experimentally-relevant energy and frequency range; (ii) depending on the
relaxation state of the material, relevant events are dominated by a dangling
bond hop or the more elaborated Wooten-Winer-Weaire (WWW) [15] bond exchange
mechanism; and (iii) the loss angle deduced from our models is compatible with
experiments. Overall, these results suggest that since multiple classes of
events can be associated with two-level systems in the same frequency range in
disordered materials, so multiple strategies might be needed to reduce their
number.

##  II Theory and methods

###  II.1 Two-level systems

Due to metastability associated with structural disorder and a distribution in
local strain, amorphous materials present, intrinsically, more possibilities
for the presence of local minima separated by low-energy barriers than their
crystalline counterparts. As they evolve over time, solid-state systems can be
pictured as transiting from one local minimum to the other on a potential
energy landscape [16]. In the context of amorphous solids, a two-level system
happens when two minima are connected by a low-energy single saddle point and
surrounded by much higher energy barriers so that the system is locally
trapped in a two-state basin.

Two-level systems (TLSs) in amorphous solids were first introduced as a way to
explain the behavior of their specific heat and thermal conductivity at low
temperatures [17]. Atoms tunneling through saddle points of TLSs provide more
degrees of freedom to the system, increasing its heat capacity. As illustrated
in Fig. 1, a TLS can be characterized by the energy of its saddle point (S)
relative to the minima (states 1 and 2), called the barrier ($V$), and the
difference in energy between the two minima, called the asymmetry ($\Delta$).

We note that while the rugged nature of the energy landscape of amorphous
solids gives a very large, quasi-continuous global distribution of barriers
and asymmetries [18], events are localized in space and distributed throughout
all the sample. Therefore they can be trapped by important gaps in the local
barrier distributions, making dissipation a local phenomenon.

Figure 1: Potential energy landscape representation of a two-level system. The
minima 1 and 2 are connected by a saddle point S. The mean barrier is V and
the asymmetry is $\Delta$.

Assuming that quantum tunneling is negligible and that the system has the time
to thermalize in a local state after crossing a saddle point, the mean rate at
which it will transition from minimum 1 to minimum 2, $\tau_{12}$ is given by
the Arrhenius law

$\tau_{12}=\tau_{0}e^{\frac{E_{S}-E_{1}}{k_{B}T}},$ (1)

where $\tau_{0}^{-1}$ is the attempt frequency, $E_{s}$ and $E_{1}$ are the
energy levels at the saddle point and the first minimum, respectively. While
activated prefactor can vary considerably [19], calculations of the attempt
frequencies in _a_ -Si with the harmonic transition state theory were found to
be close to $10^{13}$ s-1 with small fluctuations[20]. To shorten the
simulations, we shall use a constant value of $\tau_{0}=10^{-13}$ s.

Because of the finite asymmetry $\Delta$, $\tau_{12}$ is different from
$\tau_{21}$. For the whole TLS, a relaxation time, $\tau$ is defined as a
function of temperature $T,V,\Delta,$ and $\tau_{0}$ [21]

$\tau=\tau_{0}\operatorname{sech}\left(\frac{\Delta}{2k_{B}T}\right)e^{\frac{V}{k_{B}T}}.$
(2)

Mechanical dissipation in the frequency regime of Hz to MHz in amorphous
materials, in response to a strain wave of frequency $\omega$, is widely
thought to originate from the excitation of TLSs [22]. Mechanical energy from
the wave ($V\pm\Delta/2$) can push the TLS to its saddle point. The TLS will
then relax to state 2, transforming the mechanical energy into thermal energy.
Dissipation from a given TLS will be maximized when its relaxation time (Eq.
2) matches the inverse frequency of the dissipated excitation. An
approximation for relevant barriers is given by

$V\sim k_{B}T\ln{\frac{1}{\omega\tau_{0}}}.$ (3)

In the context of gravitational wave detectors, mirrors are kept at
temperatures ranging from room temperature down to cryogenic temperatures of
124 K [23]. At these temperatures and frequencies, quantum tunneling will be
completely negligible in comparison to thermal activation, such that Eq. 1 is
valid. Taking $\omega=50$ Hz we get relevant barriers of 0.28 eV and 0.67 eV
for temperatures of 124 K and 300 K, respectively.

If the entire distribution of TLSs configuration is known, the inverse quality
factor, $Q^{-1}$, also called the loss angle of the bulk material, can be
computed as follows:

$Q^{-1}_{(\omega)}=\frac{1}{E}\sum_{i}\frac{\gamma_{i}^{2}}{k_{B}T}\frac{\omega\tau_{i}}{1+\omega^{2}\tau_{i}^{2}}\operatorname{sech}^{2}\left(\frac{\Delta_{i}}{2k_{B}T}\right).$
(4)

A detailed derivation for this equation can be found in Ref. [22]. The sum
runs over every TLS in the system. $\tau_{i}$ and $\Delta_{i}$ is the
relaxation time and asymmetry of the TLS $i$. $E$ is the elastic modulus and
$\omega$ is the frequency of the applied strain. $\gamma_{i}$ is the strength
of the coupling between the TLS $i$ and the strain called the deformation
potential. This deformation potential can be obtained from the coupling tensor

$\overline{\overline{\gamma}}=\partial\Delta/\partial\overline{\overline{\epsilon}},$
(5)

where $\overline{\overline{\epsilon}}$ is the strain tensor.

Both $\gamma$ and $E$ depend on the nature of the strain. For example, to
compute the attenuation of a longitudinal wave with Eq. 4, $E$ will be the
longitudinal modulus and $\gamma$ will be the longitudinal deformation
potential. A detailed derivation has been carried out by Damart and Rodney in
Ref. [24].

###  II.2 Atomic models

To conduct this study, we consider quenched-melt and hyperuniform network
models.

200 quenched-melt systems of 1000 atoms of _a_ -Si are generated. All 200
models are prepared using the molecular dynamics (MD) simulation software
LAMMPS [25] following the same melt-quench procedure: 1000 atoms are
distributed randomly in a periodic box at a temperature of 3000 K, the system
is then cooled (quenched) at a (relatively) slow rate of $10^{11}$ K/s,
freezing in an amorphous configuration. This method has the advantage of
moderate computing cost and melt quench methods with slow cooling rate have
been shown to generate samples that are in good agreement with well-annealed
experimental _a_ -Si [26, 27].

We compare these models with a nearly hyperuniform network of _a_ -Si
developed by Hejna _et al._ [28]. This system was built using a modified
version of the WWW algorithm [15] developped by Barkema and Mousseau [29]. The
systems are then annealed for further relaxation and meticulously compared to
experimental data.

To correctly simulate the structural and vibrational properties of _a_ -Si, a
modified version of the Stillinger-Weber potential parameter set, developed by
Vink _et al._ [30] was used and applied to the original formulation.

For our topological analysis we define a cutoff value for two atoms to be
connected at the middle point between the first and the second peak of the
radial distribution function our our configurations. This middle point lies
around 3.05 Å. This is the same cutoff definition and value as Ref. [29].

###  II.3 Activation-Relaxation Technique

In this work we search for thermally activated events with characteristic
times of the order of the millisecond to the second. Here we select to used
the Activation-Relaxation Technique nouveau (ARTn) [31, 32], a saddle point
search method that is ideally suited for such tasks as it focuses on finding
high barrier events (with high characteristic times) without having to compute
every thermal atomic vibration.

This technique samples events in the potential energy landscape and finds
their barriers and asymmetries. Characteristic times are then found using Eq
2. MD-based methods can also be used to identify saddle points in such systems
[24, 33]. Because the timescale on which such method operates ( ns), only low-
energy barriers are efficiently identified. Those contribute to dissipation at
low temperature or high frequency according to Eq. 3. Hence, both methods are
complementary.

##  III Results and discussion

###  III.1 Sampling

Using ARTn, we first explore the energy landscape around the final structure
of the 200 quenched-melt 1000-atom configurations. For each sample, 30 ARTn
searches are conducted per local topology. Since, amorphous materials the
number of different environments is much larger than 1000, the number of local
environments centered on the atoms in each of the systems, that means 30 000
event searches per configuration. While events with energy barriers ranging
from 0 to a 5 eV cut-off are generated, only events with an activated barrier
$0.2\leq V\leq 0.7$ eV are considered in the study because, as explained in
section II-A, TLS with lower or higher barriers will not contribute
significantly to internal mechanical dissipation between temperatures of 124
and 300 K (Eq. 3). A cutoff in asymmetry is chosen at $V=\Delta/3$,
represented by the diagonal lines in Fig. 2, because the contribution to
dissipation of events with higher asymmetry is exponentially suppressed by
large $\Delta$ (see Eq. 4). Damart and Rodney have also shown in [24] that
TLSs with larger asymmetry do not contribute to the mechanical loss. The
barrier and asymmetry of the remaining events are plotted in Fig. 2. We note
that the TLS considered here represent only a very small fraction of all
events found in this search.

Figure 2: Barrier and asymmetry of events found by ART. The gray dotted lines
show the asymmetry cutoff.

TLS events found by ARTn can be characterized by the norm of the displacement
of the main atom (i.e. the one moving the most) when the system transitions
between two minima of the potential energy landscape, and the number of atoms
involved in the transition. An atom is considered to be active if its
displacement is greater than 0.1 Å.

Fig. 3 shows the root of the sum of squared displacements for active atoms;
the marker color corresponds to the number of active atoms according to the
scale on the right. We see that the magnitude of the atomic displacements is
correlated with the energy barrier, although with considerable dispersion. In
addition, large displacement for the main atoms is generally associated with a
larger number of active atoms (green dots). Conversely, the blue dots,
representing events with a few active atoms, are at the bottom of the graph
while the green dots, representing events with several active atoms, are at
the top. In general, TLSs events involve between 5 to 30 active atoms, a
number much smaller than the 20 to 150 active atoms found in in oxide glasses
[33, 34]. This indicates that TLSs are much more localized phenomena in a more
locally rigid structure with a higher coordination number, such as amorphous
silicon, than oxides.

Figure 3: Root of squared atomic displacement in each TLS as a function of the
energy barrier for the events plotted in Fig. 2. Symbols are color-coded
according to the number of active atoms during the event (right scale).

###  III.2 Categorizing TLSs

Following the classification of activated events in _a_ -Si developed by
Barkema and Mousseau in [35], we adopt a three-class categorization based the
evolution of the bond network during events. The first category of TLSs (type
1) is associated with a bond hop from one atom to another. These events
usually involve the diffusion of a coordination defect, such as dangling bond
defects, and will be called bond defect hopping TLSs. These jumps are
typically made possible by the movement of a single atom, leading to a
relaxation of the surrounding environment.

In the example of Fig. 4, the main atom (B, in red) is 5-fold coordinated. As
for the yellow atoms, which are the ones that change their bond with the main
atom during the process, atom A (in yellow) initially features coordination of
4, but becomes 3-fold coordinated at the end of the event. The inverse happens
to atom C (also in yellow): it is initially 3-fold coordinated and becomes
4-fold coordinated during the process.

The second category of TLSs (type 2) corresponds to the (WWW) bond exchange
mechanism [15] and is called _bond-exchange TLS_. This mechanism is commonly
observed in _a_ -Si [35, 20] and has been described in other amorphous and
crystalline solids such as graphene[36]. It involves two connected atoms
exchanging their respective bonds between them. An example is provided in Fig.
5. The atoms A and B (in red) stay connected at all times and exchange their
bonds to atoms C and D (in yellow). Interestingly, these events do not depend
on the presence of a defect to occur, as opposed to type 1 events, and can
occur in a material featuring few voids.

Figure 4: Example of bond defect hopping event (type 1). Panel a) and c)
represent the initial and final minima, respectively, while panel b) shows the
system configuration at the saddle point. The central atom is in red. The
atoms in yellow change their bonding status with the main atom. The green
lines show the trajectory of the atoms between the frames.

Figure 5: Example of a bond exchange TLS (type 2). The red atoms are the main
atoms. The yellow atoms break a bond with one red atom and form a new bond
with the other red atom. The green lines show the trajectory of the atoms
between the frames.

The third class of events includes all TLSs that do not fit the first or the
second category. Typically, these TLSs involve three or more atoms. While it
is formally possible to further analyze them [37, 20], their diversity limits
the understanding we can gain from their detailed classification so we will
focus on type 1 and 2 TLSs.

Fig. 6 shows that bond defect hopping and bond exchange TLSs exhibit different
barrier distributions. Bond defect hopping TLSs are associated with lower
barriers than bond exchange, due to the presence of a bonding defect on the
central atoms. The strain associated these defects increases the potential
energy of the two minima and lowers the saddle point barrier. Asymmetries, on
the other hand, are distributed relatively evenly for all types of TLSs.

Figure 6: Barrier distribution of bond defect hopping (cyan, full bars) and
bond exchange TLSs (empty bars).

###  III.3 Mean trends

Local environments that support TLSs display some structural characteristics
that separate them from the rest of the sample.

These can be summarized by looking at the local density obtained from
computing the Voronoi volume surrounding each atom. While the average density,
as measured over all atoms, is 2.20 g/cm3, the local value surrounding the
dominant atom associated with TLSs is only 2.08 g/cm3. These zones of low
density are associated with strained or under-coordinated atoms, creating
local instabilities that favor the formation of TLSs. Figure 7 confirms this
by showing that the dominant atom associated with TLSs have a longer bond (top
graph, dashed curves) and a wider bond angle distributions (bottom graph,
dashed curves) compared to all the atoms of the system (black solid curve).

Figure 7: Smoothed distribution of the bond lengths (top) and the bond angles
(bottom) of the main active atom in the bond defect hopping (blue dashed
curves) and bond exchange TLSs (red dotted curves). The black dotted lines
represent the first peak and bond angle distribution of the RDF of our _a_ -Si
configurations. A raw histogram (blue, top) is also shown to illustrate the
effect of the smoothing.

###  III.4 Relation between structure and mechanical loss

Experiments have shown that thermal annealing or deposition at high substrate
temperature reduces significantly the mechanical loss of amorphous materials
[38], especially silicon [14]. It is generally argued that this reduction is
due to thermal activation that allows the material to reach more relaxed
states.

To assess the importance of this effect here, we compare configurations
obtained by different means. More specifically, we compare the systems
discussed until now, obtained by melt-quench models, with a nearly
hyperuniform network (NHN) model built by Hejna _et al._ [28] and presented in
section II.2.

Table 1: Characteristics of _a_ -Si configurations prepared by melt-quench (this work) and by bond-exhange (NHN) [28, 15]. Uncertainties in the first column correspond to the standard deviation computed over 200 independent samples.  Sample | Melt-Quench | bond-exchange
---|---|---
|  | (NHN)
Energy | $-3.078\pm 0.004$ | -3.089
(eV per atom)
Density | $2.200\pm 0.006$ | 2.21
(g/cm3)
Over-coordination defect | $17\pm 4$ | 0.4
(per 1000 atoms)
Under-coordination defect | $7\pm 2$ | 0.45
(per 1000 atoms)
TLSs found |  |
Nominal | 390 | 33
Per 1000 atoms | 1.95 | 1.65

We present some physical and structural characteristics of samples prepared
using the two methods in Tab. 1. Despite similar atomic densities and
potential energy per atom, the configuration obtained with a bond-exchange
approach shows a significantly lower density of point-defects compared to
those obtained by melt-quench.

Figure 8: Energy barrier distribution for the TLSs in models obtained by the
melt-quench and NHN preparation methods (described in text).

The lower density of the NHN configuration is reflected in its TLSs
distribution generated with ARTn, as before. Figure 8 shows a depleted TLS
distribution in the NHN configuration at low barriers (between 0.1 to 0.4 eV)
as compared with the melt-quench models. This depletion is associated with a
reduced ratio of bond defect hopping (type 1) to bond exchange (type 2) TLSs.
This is not too surprising considering that bond hopping events have in
general lower barriers (Fig. 6) and often involve coordination defects, which
are much rarer in the NHN model as seen from Tab. 1. Figure 9 shows histograms
of the number of events by type in the melt-quench and NHM models normalized
by the number of atoms in each model. This confirms the much lower density of
bond defect hopping event in NHN models, but it also underlines that the
similar strain leads to a roughly similar total density of TLSs on both the
NHN and the melt-quenched systems (also see Tab 1).

Figure 9: Density of the three types of TLSs in quenched _a_ -Si (blue) and
NHN _a_ -Si (orange).

###  III.5 Loss angle calculations

Knowing the microscopic details of TLSs, it is possible to compute directly
the mechanical loss of _a_ -Si using Eq. 4. The details and approximations
used are described in section II.1.

The first step is to compute the strain-asymmetry coupling tensor in Eq 5. We
obtain this quantity by applying a small affine deformation to both minima of
every TLS. The potential energy of each minimum is then computed. This is done
for both positive and negative strain and results are averaged over both. We
verify that varying the amplitude of the deformation does not change the
results.

The deformation potentials for longitudinal and traverse strain waves are then
computed using the formula derived by Damart and Rodney [24]. For the
longitudinal deformation potentials, we get values between 0.5 and 9.1 eV and
a mean value of 4.1 eV and for the transverse case we get energies ranging
from 0.4 to 7.8 eV with a mean value of 3.4 eV. Deformation potentials are
found to be uncorrelated to the barrier or the asymmetry. Experiments have
reported average transverse deformation potentials of 1 eV [39], however those
were conducted below 1 K, so barriers involved were much lower than those
studied in this work, so we suppose the minima of active TLSs to be closer in
geometry and their deformation potential to be lower.

Mechanical loss predictions are presented in Fig. 10 for both the melt-quench
(black) and NHN configurations (red). The thin curves represent each term in
the sum (Eq. 4), while the thick curve is the total. The melt-quench
configurations lead to an almost constant $Q^{-1}$ value close to $10^{-3}$ at
temperatures between 100 and 300 K.

Figure 10: Internal mechanical dissipation computed for quenched _a_ -Si
(black) and NHN _a_ -Si (red). The pale curves show the contributions from
individual TLSs.

Because we consider only one sample, the NHN system presents a much lower
absolute number of TLSs — 33 vs 390, leading to more noise: a few TLSs are
causing peaks in the dissipation for this system, namely the ones centered at
100, 225 and 260 K. These TLSs have such high contributions because of their
very low asymmetry (0.009, 0.021 and 0.034 eV, respectively). We expect that
these peaks would flatten with more data obtained on a larger system or better
statistics.

Nevertheless, it is remarkable that the NHN configuration yields a similar
mechanical loss to that of the melt-quench around room temperature. Dropping
down to cryogenic temperatures however, the mechanical loss decreases
significantly, getting as low as $10^{-5}$ (although this result is strongly
influence by single TLSs, such as the one producing a peak near 100 K). In the
previous section we showed that the bond defect hopping TLSs, which compose
most of the low barrier TLSs in _a_ -Si (Fig. 9), are much less frequent in
the NHN configuration, causing the rarefaction of low barriers in this
configuration (Fig 8.). The same explanation applies here as well, as low
barriers are active at these low temperatures. Despite the smaller sample size
of the NHN causing sharp peaks in the mechanical loss calculations, the
reduction of low barrier TLSs is very clear (Table. 1), therefore this
decrease of the loss angle with temperature should be significant.

This behavior of $Q^{-1}$ in the NHN _a_ -Si is similar to that of
experimental hydrogenated _a_ -Si, where the mechanical loss is high at room
temperature but decreases significantly between 300 and 10K [40].

##  IV Discussion and conclusion

This study aims to identify the structural origin of internal mechanical
dissipation in amorphous solids. To do so, we focus on _a_ -Si, a classical
reference for covalent disordered materials. More precisely, we characterize
the two-level state mechanisms (TLSs) found in 200 sets of 1000 aotms of melt-
quench generated configurations and a single 20 000 atoms nearly hyperuniform
network (NHN) built built by Hejna, Steinhardt and Torquato [28] using the WWW
algorithm. Both configurations feature similar energies per atoms, however
melt-quench generated _a_ -Si has around 2.5% topological defects (3 and 5
fold cooridinated atoms), while the NHN shows almost none. The potentiel-
energy landscape is explored with ARTn [31, 32]. We keep only TLSs with small
barriers (0.2-0.7 eV) and asymmetries (< $V/3$) because they correspond to the
experimental observation window. 390 TLSs in the melt-quench _a_ -Si and 33 in
the NHN are considered relevant to our study and further analyzed.

Loss angle calculations on these systems lead to a high loss angle close to
$10^{-3}$ for both NHN and melt-quench _a_ -Si at room temperature. While the
loss angle in the melt-quench configurations stays relatively constant with
temperature, NHN _a_ -Si shows an important decrease in loss angle when the
temperature drops from 300 to 100 K, a similar behavior to that observed
experimentally in well-relaxed _a_ -Si [14].

With the detailed information obtained through ARTn, TLSs can be classified by
the bond-network change associated with the two-level system. Two-thirds of
events fall into two categories : bond defect hooping and bond exchange. These
two types of TLSs are associated with different local configurations: bond
defect hopping happen around coordination defects and where bonds are
stretched; and bond exchange TLSs are found in regions where all atoms are
4-fold coordinated but present small angle defects.

Our simulations demonstrate that various classes of TLSs can occur with
different concentrations according to the preparation schedule. Quenching from
a melt result in bond defect hopping being the dominant TLS type, because of
the high concentration of trapped point defects in these configurations. NHN
of _a_ -Si has next to no point defects. This drastically reduces the
concentration of bond hopping TLSs but has has little effect on the
concentration of bond exchange TLSs. From this we learn that TLSs are
independent: the presence of one type of TLSs do not depend on the existence
of other types. We also show that types of TLSs we find in different
configurations are also specific to the preparation, or the relaxation path of
said configuration.

Investigating other systems is necessary to improve our understanding of the
atomistic origin of TLSs and internal mechanical dissipation. Work on oxide
glasses was done in Ref. [33, 34, 24]. In Ref. [24], Damart and Rodney
analyzed TLSs in SiO2, and found different archetypes of TLSs than in _a_ -Si,
such as rotations of Si-O-Si chains. This suggests that TLSs are also specific
to the the nature of the material, for instance between a-Si and and oxide
glasses such as SiO2 and Ta2O5.

Decreasing TLSs in amorphous materials, means therefore understanding the
nature of events for specific systems, but also addressing separately each of
the potential classes of TLSs.

Because of the importance of specificity, further study would benefit from
better physical description to ensure that the atomistic details of TLSs are
accurate. Since larger amorphous configurations tend to better reproduce
experimental observations [29, 41], ab initio approaches are probably not
appropriate. Recent development in machine-learning forcefields [27] offer
here, new opportunity to deepen our understanding of these fascinating
problems.

##  V Code and data availability

The ARTn packages as well as the data reported here are distributed freely.
Please contact Normand Mousseau (normand.mousseau@umontreal.ca).

###### Acknowledgements.

The authors would like to thank K. Prasai and the other members of the Optical
Working Group of the LIGO Scientific Collaboration for fruitful discussions.
They also thank M. Hejna _et al._ for sharing their NHN model. This project is
supported by a team grant from the Fonds de recherche du Québec - Nature et
technologie. NM and FS acknowledges partial support through Discovery grants
from the Natural Science and Engineering Research Council of Canada (NSERC).
CL holds a NSERC graduate scholarship. We are grateful to Calcul Québec and
Compute Canada for generous allocation of computational resources.

## References

  * Abbott  _et al._ [2016] B. P. Abbott _et al._ (LIGO Scientific Collaboration and Virgo Collaboration), Phys. Rev. Lett. 116, 061102 (2016).
  * Abbott  _et al._ [2021] R. Abbott _et al._ , Gwtc-3: Compact binary coalescences observed by ligo and virgo during the second part of the third observing run (2021).
  * Harry  _et al._ [2006] G. M. Harry, M. R. Abernathy, A. E. Becerra-Toledo, H. Armandula, E. Black, K. Dooley, M. Eichenfield, C. Nwabugwu, A. Villar, D. R. M. Crooks, G. Cagnoli, J. Hough, C. R. How, I. MacLaren, P. Murray, S. Reid, S. Rowan, P. H. Sneddon, M. M. Fejer, R. Route, S. D. Penn, P. Ganau, J.-M. Mackowski, C. Michel, L. Pinard, and A. Remillieux, Classical and Quantum Gravity 24, 405 (2006).
  * Kubo [1966] R. Kubo, Reports on Progress in Physics 29, 255 (1966), publisher: IOP Publishing.
  * Saulson [1990] P. R. Saulson, Physical Review D 42, 2437 (1990), publisher: American Physical Society.
  * Levin [1998] Y. Levin, Phys. Rev. D 57, 659 (1998).
  * Palma  _et al._ [1996] G. M. Palma, K.-A. Suominen, and A. K. Ekert, Proc. Math. Phys. Eng. Sci. 452, 567 (1996).
  * Abernathy  _et al._ [2021] M. Abernathy, A. Amato, A. Ananyeva, S. Angelova, B. Baloukas, R. Bassiri, G. Billingsley, R. Birney, G. Cagnoli, M. Canepa, M. Coulon, J. Degallaix, A. D. Michele, M. A. Fazio, M. M. Fejer, D. Forest, C. Gier, M. Granata, A. M. Gretarsson, E. M. Gretarsson, E. Gustafson, E. J. Hough, M. Irving, É. Lalande, C. Lévesque, A. W. Lussier, A. Markosyan, I. W. Martin, L. Martinu, B. Maynard, C. S. Menoni, C. Michel, P. G. Murray, C. Osthelder, S. Penn, L. Pinard, K. Prasai, S. Reid, R. Robie, S. Rowan, B. Sassolas, F. Schiettekatte, R. Shink, S. Tait, J. Teillon, G. Vajente, M. Ward, and L. Yang, Classical and Quantum Gravity 38, 195021 (2021).
  * Granata  _et al._ [2020] M. Granata, A. Amato, G. Cagnoli, M. Coulon, J. Degallaix, D. Forest, L. Mereni, C. Michel, L. Pinard, B. Sassolas, and J. Teillon, Applied Optics 59, A229 (2020), publisher: Optical Society of America.
  * Lalande  _et al._ [2021] E. Lalande, A. W. Lussier, C. Lévesque, M. Ward, B. Baloukas, L. Martinu, G. Vajente, G. Billingsley, A. Ananyeva, R. Bassiri, M. M. Fejer, and F. Schiettekatte, Journal of Vacuum Science & Technology A 39, 043416 (2021), https://doi.org/10.1116/6.0001074 .
  * Nawrodt  _et al._ [2011] R. Nawrodt, S. Rowan, J. Hough, M. Punturo, F. Ricci, and J.-Y. Vinet, General Relativity and Gravitation 43, 593 (2011).
  * Lewis [2022] L. J. Lewis, Journal of Non-Crystalline Solids 580, 121383 (2022).
  * Birney  _et al._ [2018a] R. Birney, J. Steinlechner, Z. Tornasi, S. MacFoy, D. Vine, A. S. Bell, D. Gibson, J. Hough, S. Rowan, P. Sortais, S. Sproules, S. Tait, I. W. Martin, and S. Reid, Phys. Rev. Lett. 121, 191101 (2018a).
  * Liu  _et al._ [2014] X. Liu, D. R. Queen, T. H. Metcalf, J. E. Karel, and F. Hellman, Physical Review Letters 113, 025503 (2014).
  * Wooten  _et al._ [1985] F. Wooten, K. Winer, and D. Weaire, Physical Review Letters 54, 1392 (1985), publisher: American Physical Society.
  * Wales [2004] D. Wales,  _Energy Landscapes: Applications to Clusters, Biomolecules and Glasses_, Cambridge Molecular Science (Cambridge University Press, Cambridge, 2004).
  * Anderson  _et al._ [1972] P. W. Anderson, B. I. Halperin, and c. M. Varma, The Philosophical Magazine: A Journal of Theoretical Experimental and Applied Physics 25, 1 (1972), https://doi.org/10.1080/14786437208229210 .
  * Kallel  _et al._ [2010] H. Kallel, N. Mousseau, and F. Schiettekatte, Phys. Rev. Lett. 105, 045503 (2010).
  * Gelin  _et al._ [2020] S. Gelin, A. Champagne-Ruel, and N. Mousseau, Nature Communications 11, 3977 (2020), bandiera_abtest: a Cc_license_type: cc_by Cg_type: Nature Research Journals Number: 1 Primary_atype: Research Publisher: Nature Publishing Group Subject_term: Atomistic models;Statistical physics Subject_term_id: atomistic-models;statistical-physics.
  * Valiquette and Mousseau [2003] F. Valiquette and N. Mousseau, Physical Review B 68, 125209 (2003), publisher: American Physical Society.
  * Jäckle  _et al._ [1976] J. Jäckle, L. Piché, W. Arnold, and S. Hunklinger, Journal of Non-Crystalline Solids 20, 365 (1976).
  * Phillips [1987] W. A. Phillips, Reports on Progress in Physics 50, 1657 (1987).
  * Birney  _et al._ [2018b] R. Birney, J. Steinlechner, Z. Tornasi, S. MacFoy, D. Vine, A. Bell, D. Gibson, J. Hough, S. Rowan, P. Sortais, S. Sproules, S. Tait, I. Martin, and S. Reid, Physical Review Letters 121, 191101 (2018b).
  * Damart and Rodney [2018] T. Damart and D. Rodney, Physical Review B 97, 014201 (2018), publisher: American Physical Society.
  * Thompson  _et al._ [2022] A. P. Thompson, H. M. Aktulga, R. Berger, D. S. Bolintineanu, W. M. Brown, P. S. Crozier, P. J. in ’t Veld, A. Kohlmeyer, S. G. Moore, T. D. Nguyen, R. Shan, M. J. Stevens, J. Tranchida, C. Trott, and S. J. Plimpton, Computer Physics Communications 271, 108171 (2022).
  * Ishimaru  _et al._ [1997] M. Ishimaru, S. Munetoh, and T. Motooka, Physical Review B 56, 15133 (1997).
  * Deringer  _et al._ [2018] V. L. Deringer, N. Bernstein, A. P. Bartók, M. J. Cliffe, R. N. Kerber, L. E. Marbella, C. P. Grey, S. R. Elliott, and G. Csányi, The Journal of Physical Chemistry Letters 9, 2879 (2018).
  * Hejna  _et al._ [2013] M. Hejna, P. J. Steinhardt, and S. Torquato, Physical Review B 87, 245204 (2013).
  * Barkema and Mousseau [2000] G. T. Barkema and N. Mousseau, Physical Review B 62, 4985 (2000), publisher: American Physical Society.
  * Vink  _et al._ [2001] R. L. C. Vink, G. T. Barkema, W. F. van der Weg, and N. Mousseau, Journal of Non-Crystalline Solids 282, 248 (2001).
  * Barkema and Mousseau [1996] G. T. Barkema and N. Mousseau, Physical Review Letters 77, 4358 (1996), publisher: American Physical Society.
  * Malek and Mousseau [2000] R. Malek and N. Mousseau, Physical Review E 62, 7723 (2000), publisher: American Physical Society.
  * Hamdan  _et al._ [2014] R. Hamdan, J. P. Trinastic, and H. P. Cheng, The Journal of Chemical Physics 141, 054501 (2014), publisher: American Institute of Physics.
  * Trinastic  _et al._ [2016] J. P. Trinastic, R. Hamdan, C. Billman, and H.-P. Cheng, Physical Review B 93, 014105 (2016), publisher: American Physical Society.
  * Barkema and Mousseau [1998] G. T. Barkema and N. Mousseau, Physical Review Letters 81, 1865 (1998), publisher: American Physical Society.
  * Stone and Wales [1986] A. Stone and D. Wales, Chemical Physics Letters 128, 501 (1986).
  * Mousseau and Barkema [2000] N. Mousseau and G. T. Barkema, Physical Review B 61, 1898 (2000).
  * Vajente  _et al._ [2018] G. Vajente, R. Birney, A. Ananyeva, S. Angelova, R. Asselin, B. Baloukas, R. Bassiri, G. Billingsley, M. M. Fejer, D. Gibson, L. J. Godbout, E. Gustafson, A. Heptonstall, J. Hough, S. MacFoy, A. Markosyan, I. W. Martin, L. Martinu, P. G. Murray, S. Penn, S. Roorda, S. Rowan, F. Schiettekatte, R. Shink, C. Torrie, D. Vine, S. Reid, and R. X. Adhikari, Classical and Quantum Gravity 35, 075001 (2018).
  * Fefferman  _et al._ [2017] A. Fefferman, A. Maldonado, E. Collin, X. Liu, T. Metcalf, and G. Jernigan, Journal of Low Temperature Physics 187, 10.1007/s10909-016-1686-6 (2017).
  * Liu  _et al._ [1997] X. Liu, B. E. White, Jr., R. O. Pohl, E. Iwanizcko, K. M. Jones, A. H. Mahan, B. N. Nelson, R. S. Crandall, and S. Veprek, Physical Review Letters 78, 4418 (1997).
  * Igram  _et al._ [2018] D. Igram, B. Bhattarai, P. Biswas, and D. A. Drabold, Journal of Non-Crystalline Solids 492, 27 (2018).


# WAVE: Machine Learning for Full-Waveform Time-Of-Flight Detectors

Glenn R. Jocher<EMAIL_ADDRESS>Kurtis Nishimura
<EMAIL_ADDRESS>John Koblanski<EMAIL_ADDRESS>Viacheslav
A. Li<EMAIL_ADDRESS>Ultralytics LLC, Arlington, VA, USA  Department of
Physics and Astronomy, University of Hawaii, Honolulu, HI 96822, USA  Lawrence
Livermore National Laboratory, CA 94550

###### Abstract

We propose a WAveform Vector Exploitation (WAVE) deep neural network for full-
waveform Time-Of-Flight (TOF) physics detectors, and evaluate its performance
against traditional reconstruction techniques via Monte Carlo study of a small
plastic-scintillator scatter camera. Ultralytics LLC
(https://www.ultralytics.com) provides WAVE freely under the open source
GPL-3.0 license at https://github.com/ultralytics/wave.

###### keywords:

Machine Learning , Deep Neural Networks , Time-Of-Flight , Scintillation ,
SiPM

###### Contents

  1. 1 Introduction
  2. 2 Time-Of-Flight Methods
    1. 2.1 Constant Amplitude (CA) thresholding
    2. 2.2 Constant Fraction (CF) thresholding
    3. 2.3 Machine Learning
  3. 3 Results
  4. 4 Conclusion
  5. 5 Acknowledgements
  6. A

##  1 Introduction

Figure 1: TOF-based scintillating-fiber scatter camera. Two SiPM arrays (SensL
ArrayJ-60035-64P[7]) sandwich a vertical bundle of 5x5 mm square, 20 cm long
plastic fibers (one fiber shown). Matching upper and lower SiPM waveform pairs
(shown in Figure 2) are input to WAVE for scatter reconstruction.

This paper introduces a Machine Learning (ML) Time-Of-Flight (TOF) solution
called WAveform Vector Exploitation (WAVE), and evaluates its performance in a
common physics setting against traditional TOF techniques. TOF techniques are
widely used in the physics community for particle tracking and reconstruction
(e.g., [1, 2, 3]), in stripline anode electronics readouts[4], and in the
medical industry in Positron Emission Tomography (PET) scanners. Comparisons
of different timing techniques relevant to TOF have been performed in the
past, in particular by Genat et al.[6], though for only single waveforms, and
ML solutions were not included.

Traditional TOF relies on a two step process: first timing waveform arrivals
at separate locations, and subsequently reconstructing the source of the
waveforms given the measured times. This technique is heavily reliant on
accurate timing of the waveform arrivals. In contrast, WAVE directly solves
for the waveform’s source in a single pass, without the need for any
explicitly calculated features. This means that WAVE does not need to solve
separately for the arrival time of either waveform, a novel concept which
appears to provide performance benefits in our study.

##  2 Time-Of-Flight Methods

Assuming known propagation speeds $v$ and similar waveform shapes, one may
solve for the waveform origin time $t$ and position $x$ via Equations 1 and 2,
respectively.

$\displaystyle t$ $\displaystyle=$ $\displaystyle\frac{t_{1}+t_{2}-L/v}{2}$
(1) $\displaystyle x$ $\displaystyle=$ $\displaystyle
v\left(\frac{t_{1}-t_{2}}{2}\right)$ (2)

where $t_{1}$ and $t_{2}$ are the measured waveform arrival times at sides 1
and 2, and $L$ is the distance between measurement locations. $x_{0}$ is at
the center of the fiber, and $-x$ trend towards side 1.

In our example application, Equations 1 and 2 are applied to the particle
scatter camera shown Figure 1 to reconstruct the scatters within the detector
volume. The detector consists of one hundred $5\times 5\times 200$ mm3 plastic
scintillating fibers arranged in a $10\times 10$ grid. Gamma or neutron
scatters in these fibers create isotropic scintillation light. Due to a high
index of refraction mismatch between the fiber cores (n = 1.59) and the small
air gaps between the fibers (n = 1.00), much of this light is captured in the
fiber cores and reflected to photodetectors at each end, modeled as Silicon
PhotoMultipliers (SiPMs)[7] for this study. About 25% of the light produced in
the fiber is retained in this manner. The other 75% attenuates or leaves the
detector.

When a photon arrives at one of the SiPM microcells it is converted into an
observable current with about 50% Photon Detection Efficiency (PDE)[15],
resulting in a voltage response that is converted to a digital waveform. For
our proof of concept simulations, we assume a response consistent with the
‘fast‘’ output of commercially available SiPMs, and assume that these voltage
waveforms are digitized by DRS4 integrated circuits[8], which are
differentially-driven switched capacitor array ASICs developed by the Paul
Scherrer Institut (PSI) in Switzerland. The DRS4 can sample a waveform at up
to 5 GigaSamples Per Second (GSPS) with a buffer depth of 1024 samples.
Similar switched capacitor array chips include LABRADOR[9] and PSEC4[10].

Figure 2: Simulated waveform pairs produced by the detector in Figure 1. The
arrival time difference between the side 1 and side 2 waveforms is used in
Equations 1 and 2 to locate the source time and position along the fiber. WAVE
accepts these joined waveforms as inputs of length 256 + 256 = 512.

We use a validated GEANT[11] \+ MATLAB[12] MC engine to model our waveforms.
This MC package was created by Ultralytics for the mini-TimeCube[13] and Large
Area Picosecond Photodetector (LAPPD)[14] collaborations, and allows quick
prototyping of new detector designs. Waveforms modeled by this package, along
with known truth information are used to analyze the TOF performance of 3
methods:

  1. 1.

Constant Amplitude (CA) TOF

  2. 2.

Constant Fraction (CF) TOF

  3. 3.

Machine Learning (WAVE)

The MC engine places ‘perfect’ point sources randomly within the scintillating
fibers, randomized in time and energy, to create a dataset of 100,000 events,
which form the basis for all the results published in this paper. The timing
of these point sources was randomized to start anywhere within a 3 ns window,
which is a similar trigger jitter to what the DRS4 evaluation board
experiences on internally triggered events. The energy of the points was
randomized per an exponential distribution with a 250 keVee (keV electron
equivalent) exponential constant, to mimic the low energy events a plastic
$\sim$MeV neutron scatter camera would likely target.

Digital readout noise was modeled as white and Gaussian with a standard
deviation of 0.5 mV. Electronics saturation was modeled at $\pm$ 1.1 V. Single
Photo-Electron (PE) waveform amplitudes in our model averaged 0.4 mV $\pm$ 10%
at +5 V SiPM overvoltage. All SiPM characteristics including Photo-Detection
Efficiency (PDE), spectral response, dark count rates (80 kHz/mm2 at 20∘ C)
and microcell crosstalk (though not afterpulsing) were modeled per the SensL
J-60035 datasheet[15]. A Gaussian Transit Time Spread (TTS) of 120 ps
$1\sigma$ was modeled for SensL J-60035 PEs as well.

Optical properties of the scintillating fibers, including scintillation
spectrum, index of refraction and yield of 8000 photons/MeV were modeled using
a Saint Gobain BCF-10 datasheet[5]. Optical surfaces are assumed to be
perfectly flat.

Although we have modeled the full detector, the relevant quantities of
comparison for this study are the capability to localize the interaction
position along the longitudinal axis of the fiber, as well as the absolute
time of the interaction. To that end, the comparisons between the following
techniques are done on their performance for a single bar at a time, not
including effects of optical crosstalk that might be present in the full
detector.

###  2.1 Constant Amplitude (CA) thresholding

CA thresholding supplies Equation 1 with crossing times $t_{1}$ and $t_{2}$
extracted from the side 1 and 2 waveforms at the same amplitude. This
amplitude is typically set low, in order to extract the very beginning of a
waveform risetime, yet high enough to avoid false positives from the noise
(which in our model is 0.5 mV 1$\sigma$). Figure 3 shows a study we performed
on this threshold amplitude to determine the optimum value in our application.
This search yielded a best amplitude of 6 mV, which we use in reporting our CA
results.

Figure 3: CA TOF performance vs amplitude threshold. Amplitudes from 6 - 200
mV are analyzed. Typical waveforms amplitudes are 10-20 mV, as seen in Figure
2, but can exceed 500 mV for bright events ($\sim$1 MeVee). The best
performance was observed at the lowest threshold (while staying above the
noise). Candidate efficiency trends downwards at higher thresholds.

###  2.2 Constant Fraction (CF) thresholding

CF thresholding supplies Equation 1 with crossing times $t_{1}$ and $t_{2}$
extracted from the waveforms at the same fraction of max amplitude. This
fraction is typically set at 1/2. Figure 4 shows a study we performed to
determine an optimum value for this parameter. The best results are observed
around 0.3 - 0.5 fraction. 0.5 is a commonly used, so we will also use 0.5 for
reporting CF results.

###  2.3 Machine Learning

Figure 4: CF TOF performance vs fraction threshold. Fractions from 0.0 - 0.99
were analyzed. The best performance was observed around $\sim$0.4.

WAVE is a 5-layer function-fitting feedforward network, shown in Figure 5. We
constructed WAVE in Python using the PyTorch 0.4.0[16] package, and made it
freely available under the open source GPL-3.0 license at
https://github.com/ultralytics/wave. A 512-length input layer, composed of two
joined waveforms, connects to 3 hidden layers of size 76, 23 and 7, with an
output layer of size 2: 1 position $x$ and 1 time $t$. Side 1 and Side 2
waveforms shown in Figure 2 of length 256 each are joined sequentially to
create each 512-length input vector.

All hidden layer neurons are linear (weights + bias), with $tanh()$
activations, except layer 4, which lacks activation. This structure was
settled on using trial and error, though the exact geometry is not critical,
as we observed that our Section 3 results were robust to moderate changes in
network geometry. Our 100,000 event dataset was split into 3 groups:

  1. 1.

Train (70%). Used during WAVE training to update network parameters (weights
and biases).

  2. 2.

Validate (15%). Used to stop WAVE training to prevent over-training. Not used
to update network parameters.

  3. 3.

Test (15%). Used to evaluate performance of CA, CF and WAVE. Not used during
WAVE training.

Figure 5: WAVE neural network: a 5 layer function-fitting feedforward deep
neural network used to reconstruct scatters in the Figure 1 detector. 512
inputs (the joined waveforms shown in Figure 2) pass through 3 fully-connected
hidden layers (weights + bias + $tanh()$ activation). The 2 outputs are
scatter position $x$ (mm) and time $t$ (ps) in the detector.

We used an Adam optimizer with near-default hyperparameters ($lr=0.001$,
$\epsilon=1\times 10^{-8}$, $\beta_{1}=0.9$, $\beta_{2}=0.999$) and a Mean
Squared Error (MSE) loss criterion. Supervised learning (on the Train and
Validate datasets) was performed on a Google Cloud Platform (GCP) Compute
Engine Virtual Machine (VM) instance with one Nvidia Tesla P100 GPU, taking
$\sim$ 5 minutes of wall-clock time for 50,000 epochs, and yielding a minimum
validation loss of 0.023 after 14,278 epochs. The epoch 14,278 model was
frozen for inference on the Test dataset to produce the Section 3 results.

##  3 Results

Results are summarized in Table 1 and detailed in Figure 6 as a function of
scatter energy. WAVE performance exceeds CA and CF in our study by a wide
margin in both position and time reconstruction. The largest performance
increases are observed at the lowest energy events, the most challenging
region of the parameter space to reconstruct.

method | candidates | position $x$ | time $t$
---|---|---|---
| fraction | mm 1$\sigma$ | ps 1$\sigma$
CA - 6 mV | 0.95 | 38 | 627
CF - 0.5 | 1.00 | 26 | 150
WAVE | 1.00 | 11 | 90
Table 1: MC study results (lower 1$\sigma$ is better), showing WAVE exceeding
CA and CF performance in both spatial and temporal reconstruction. This table
summarizes Figure 6, which plots the above metrics vs scatter energy. Figure
6: Results for CA, CF, and WAVE vs scatter energy E. TOP: Reconstruction
errors. BOTTOM: Reconstruction error 1$\sigma$, averaged over energy in Table
1 (biases are ignored in our results, though we note WAVE automatically learns
and eliminates these). WAVE performance exceeds CA and CF for both metrics
across all energies in our dataset.

##  4 Conclusion

Table 1 and Figure 6 show clear performance increases for WAVE compared to
traditional techniques across all regions of the parameter space evaluated in
our study.

The next steps are the application of this technique to real-world data.
Challenges to real-world implementation include the need to accumulate
training data of sufficient quantity and variation to generalize well in real-
world scenarios. For applications where such training data it attainable, it
appears WAVE may offer significant performance benefits over traditional
solutions.

##  5 Acknowledgements

This effort was initiated and funded by Ultralytics LLC
(https://www.ultralytics.com) in collaboration with the University of Hawaii
and Lawrence Livermore National Laboratory. 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. Release LLNL-JRNL-761634.

## References

## References

  * [1] B.W. Adams, A. Elagin, H.J. Frisch, R. Obaid, E. Oberla, A. Vostrikov, R.G. Wagner, J. Wang, M. Wetstein, Timing characteristics of Large Area Picosecond Photodetectors, Nuclear Instruments and Methods in Physics Research Section A, Volume 795, 21 September 2015, Pages 1-11, ISSN 0168-9002, http://dx.doi.org/10.1016/j.nima.2015.05.027
  * [2] A. Abashian, et al., The Belle detector, Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment, Volume 479, Issue 1, 21 February 2002, Pages 117-232, ISSN 0168-9002. http://dx.doi.org/10.1016/S0168-9002(01)02013-7
  * [3] Wei-Dong Li, Ya-Jun Mao, and Yi-Fang Wang, Int. J. Mod. Phys. A 24, 9 (2009). http://dx.doi.org/10.1142/S0217751X09046424
  * [4] Grabas, H. and Obaid, R. and Oberla, E. and Frisch, H. and Genat, J.-F. and Northrop, R. and Tang, F. and McGinnis, D. and Adams, B. and Wetstein, M., “RF strip-line anodes for psec large-area mcp-based photodetectors” Nuclear Instruments and Methods, A 711 (0) (2013) 124-131 http://dx.doi.org/10.1016/j.nima.2013.01.055
  * [5] Saint Gobain Scintillating Fibers website
  * [6] Jean-Francois Genat, Gary Varner, Fukun Tang, Henry Frisch, Signal processing for picosecond resolution timing measurements, Nucl.Instrum.Meth.A, Volume 607, Issue 2, 11 August 2009, Pages 387-393, ISSN 0168-9002, http://dx.doi.org/10.1016/j.nima.2009.05.193
  * [7] SensL ArrayJ-60035-64P-PCB SiPM Array
  * [8] “DRS4 Documentation,” Paul Scherrer Institut Website.
  * [9] G.S. Varner, L.L. Ruckman, J.W. Nam, R.J. Nichol, J. Cao, P.W. Gorham, M. Wilcox, The large analog bandwidth recorder and digitizer with ordered readout (LABRADOR) ASIC, Nucl.Instrum.Meth.A, Volume 583, Issues 2–3, 21 December 2007, Pages 447-460, ISSN 0168-9002, http://dx.doi.org/10.1016/j.nima.2007.09.013
  * [10] Eric Oberla, Jean-Francois Genat, Hervé Grabas, Henry Frisch, Kurtis Nishimura, Gary Varner, A 15 GSa/s, 1.5 GHz bandwidth waveform digitizing ASIC, Nucl.Instrum.Meth.A, Volume 735, 21 January 2014, Pages 452-461, ISSN 0168-9002, http://dx.doi.org/10.1016/j.nima.2013.09.042
  * [11] GEANT4 - A Simulation Toolkit, S. Agostinelli et al., Nuclear Instruments and Methods A 506 (2003) 250-303. http://geant4.cern.ch/
  * [12] MATLAB R2018a, a computer program. http://www.mathworks.com
  * [13] “Invited Article: miniTimeCube,” Review of Scientific Instruments 87, 021301 (2016).
  * [14] G. Jocher et al., ”Multi-Photon Disambiguation on Stripling-Anode Multi-Channel Plates,” Nucl. Instr. and Meth. A Volume 822, (2016), Pages 25-33, ISSN 0168-9002
  * [15] SensL MicroFJ-60035-TSV SiPM
  * [16] PyTorch, a python package. https://pytorch.org/

## Appendix A

Figure 7: Results vs scatter position (in fiber). Edge affects appear in all 3
fitting techniques, though interestingly WAVE is the only technique to show
improved performance near the fiber edges. WAVE outperform CA and CF across
all regions of this parameter space (scatter position space), as it does in
Figures 6 and 8. Figure 8: Results vs scatter time. No significant edge
effects are observed. WAVE outperform CA and CF across all regions of this
parameter space (scatter time space), as it does in Figures 6 and 7. Figure 9:
WAVE training loss vs number of epochs trained. Validation and Test losses
reach a minimum around 20,000 epochs, after which they begin to increase. Test
losses continue to decrease. Training is terminated after 5,000 validation
epochs fail to register a new absolute minimum. Figure 10: Hyperparameter
study: effect of learning rate on best validation loss and epochs required to
reach best validation loss. The default learning rate of 0.001 appears to be
well placed, as it allows for a near optimal final loss in a relatively short
training time. Figure 11: Activation function study. 5 activation function
candidates are evaluated for use in WAVE. ELU and Tanh both appear well suited
to the problem. Softsign produces the worst results. Figure 12: Hidden layer
structure study. Different hidden layer structures were analyzed, starting at
1 hidden layer with 23 neurons and moving to 5 hidden layers with 181 neurons
in the first layer, 64 in the second layer, etc. The largest performance
increase was observed when moving from 1 to 2 layers, with corresponding
training speed improvements as well, yet with diminishing returns beyond 3-4
layers. Figure 13: WAVE pytorch training speed on various hardware, measured
in epochs per second. About 20,000 epochs are required for training,
corresponding to $\sim$ 100 seconds on a Nvidia Tesla P100 GPU.