Paper The following article is Open access

Viterbi decoding of CRES signals in Project 8

, , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , and

Published 6 May 2022 © 2022 The Author(s). Published by IOP Publishing Ltd on behalf of the Institute of Physics and Deutsche Physikalische Gesellschaft
, , Citation A Ashtari Esfahani et al 2022 New J. Phys. 24 053013 DOI 10.1088/1367-2630/ac66f6

Download Article PDF
DownloadArticle ePub

You need an eReader or compatible software to experience the benefits of the ePub3 file format.

1367-2630/24/5/053013

Abstract

Cyclotron radiation emission spectroscopy (CRES) is a modern approach for determining charged particle energies via high-precision frequency measurements of the emitted cyclotron radiation. For CRES experiments with gas within the fiducial volume, signal and noise dynamics can be modelled by a hidden Markov model. We introduce a novel application of the Viterbi algorithm in order to derive informational limits on the optimal detection of cyclotron radiation signals in this class of gas-filled CRES experiments, thereby providing concrete limits from which future reconstruction algorithms, as well as detector designs, can be constrained. The validity of the resultant decision rules is confirmed using both Monte Carlo and Project 8 data.

Export citation and abstract BibTeX RIS

1. Introduction

Cyclotron radiation emission spectroscopy (CRES) is a recently-demonstrated technique for high-resolution energy measurements of individual charged particles [1, 2]. When placed within a magnetic field, charged particles travel in helical trajectories, thereby emitting cyclotron radiation. This radiation can serve as a sensitive, non-destructive probe into the kinematic properties of the radiating particle.

In a uniform magnetic field B, a charged particle with kinetic energy Ke will emit cyclotron radiation with angular frequency

Equation (1)

where γ is the Lorentz factor, and q and m are the charge magnitude and mass of the particle, respectively. Therefore, for known magnetic fields, accurate frequency measurements of the emitted radiation map into a high-resolution measurement of the kinetic energy.

While CRES can have many potential use cases [3], the discussion presented here will be in the context of Project 8 [4], as this is the current largest application of the technique. Project 8 will attempt to determine the absolute mass scale of the neutrino via spectroscopic measurement of electrons near the endpoint of the beta decay spectrum of gaseous tritium. All currently published implementations of CRES experiments contain gas within the fiducial volume as a source of radiating electrons.

CRES data are time series of voltages. These voltages are induced by the electromagnetic radiation incident to receiver elements, and then amplified, downmixed, and digitized. The data consists of both Johnson–Nyquist noise, resulting from thermal motion of charge carriers in the amplification chain, and if an electron is present, a coherent cyclotron radiation signal. CRES signals are modelled by linear frequency-modulated chirps [5] of the form: s(t) = A exp i(ϕ + ω0 t + αt2/2), as the radiation of electron energy yields a linearly increasing instantaneous cyclotron frequency (equation (1): α = ωc Pe/()). Assuming the power released via the electron's radiation, Pe, is known, application of the dechirping operator $({\text{e}}^{-\text{i}\alpha {t}^{2}/2})$ to the recorded CRES signal converts it to a time-limited monofrequency sinusoid. Data is typically represented in time-frequency space via the short-time Fourier transform.

One of the primary challenges in CRES experiments is the reliable detection of signals over thermal noise, given that the power radiated via cyclotron emission is relatively weak. In a tritium beta decay experiment operating at a magnetic field of 1 T, only ≈1 fW is radiated by a given 18.6 keV endpoint electron. While efficiently detecting such events with minimal latency is critical to maximizing the sensitivity reach of the experiment [6], one must simultaneously minimize the number of false positives, which occur when noise fluctuations are incorrectly reconstructed as genuine electron events. In the context of Project 8, false or mis-reconstructed events constitute the dominant background, with a single false positive near the tritium endpoint energy having the potential of jeopardizing the experimental sensitivity of a multi-year discovery run.

In this document, we apply the framework of detection theory to CRES signals in gas-filled detectors. While the signal phenomenology can have additional spectral complexity within specific detector and magnetic field designs [7], the basic requirements of a CRES experiment impose several key features onto the resultant data, which lead to unique analytic decision rules and detection properties.

2. Viterbi algorithm & hidden Markov models

In this section, we discuss a novel application of the Viterbi algorithm [8, 9] for generic CRES signal detection. The Viterbi algorithm is a well-established dynamic programming method, with applications in diverse fields such as natural language processing, bioinformatics, and telecommunications [1012]. Fundamentally, the algorithm allows one to reconstruct the most probable underlying state sequence given time series data.

The application of the Viterbi algorithm to CRES data is motivated by re-expressing the data sequence as a so-called hidden Markov model (HMM) [13, 14], as illustrated in figure 1. By definition, in a Markov model, the state of the system at a given time depends only on the state immediately preceding it. In a HMM, there is an underlying stochastic Markov model controlling the sequence of states, which are hidden from the observer. Each hidden state (upper) then emits an observable symbol (lower), drawn from a known probability distribution function. Given the observation of a sequence of emitted symbols as a result of a HMM, the Viterbi algorithm reconstructs the most probable sequence of hidden states. The algorithm has asymptotic time complexity $\mathcal{O}({S}^{2}\times D)$, where S is the number of hidden states and D is the number of observed data samples.

Figure 1.

Figure 1. State diagram of HMM characteristic of CRES experiments. Hidden states (upper) produce observed data (lower) with different probabilities.

Standard image High-resolution image

In CRES experiments, we claim that the data can generally be represented in terms of a HMM, between hidden states ${\mathcal{H}}_{0}$ and ${\mathcal{H}}_{1}$, corresponding to the noise-only and signal-plus-noise hypotheses, respectively. The transition probability Tij denotes the probability that the system switches from hidden state i to j at each discrete time sample.

As one does not have direct knowledge about whether ${\mathcal{H}}_{0}$ or ${\mathcal{H}}_{1}$ is true at a given time, one must instead use the data, in the form of Fourier voltage amplitudes. For simplicity, we initially consider binary-quantized data, referred to as the sparse spectrogram, in which the Fourier amplitudes are compared to a quantization threshold and assigned either a 0 or a 1 (figure 2: lower). In section 2.2, we derive the detector performance for fully continuous, or raw, data.

Figure 2.

Figure 2. Example of a three-track event generated by a HMM in short-time Fourier transform time-frequency representation, illustrated in both the raw (top) and sparse (bottom) spectrogram data formats. (Red line: MC truth; blue line: Viterbi reconstruction.) Bins have dimensions 40.96 μs × 24.4 kHz.

Standard image High-resolution image

Application of the Viterbi algorithm requires that data be describable by a HMM such that signal and noise state transitions satisfy Markovian time-independence. CRES data satisfy this condition due to a few key characteristics of the experimental setup, which are uncommon with respect to the majority of generic signal detection problems, such as those that arise in radar or sonar applications.

Firstly, so long as the radioactive source gas (e.g. tritium) is properly replenished, the state transition from noise-only to signal-plus-noise $({\mathcal{H}}_{0}\to {\mathcal{H}}_{1})$ will be constant, as decays of the gaseous source occur spontaneously and with constant rate. Secondly, in Project 8 and in CRES experiments in general, the inverse process ${\mathcal{H}}_{1}\to {\mathcal{H}}_{0}$, corresponding to the end of an electron event, also occurs with a time-independent rate, since the time-duration of events is overwhelmingly limited by scattering between the electrons and the residual gas within the fiducial volume [15]. For electrons within a volume of gas with uniform density, traveling at nearly constant velocity, the probability of scattering with a gas molecule so as to expel the electron from the magnetic trap is also nearly independent of time. We can therefore see how the unique case of CRES electrons embedded within a gaseous volume allows for this novel application of the Viterbi algorithm by equating electron appearances and disappearances to state transitions of a HMM.

The Viterbi algorithm is a maximal likelihood algorithm that explicitly evaluates the probability of all possible hidden state sequences, given the observed data. As a result, it is capable of explicitly determining the exact conditions such that a data segment is classified as being caused due to a genuine signal event, as opposed to being caused by consecutive noise fluctuations. We denote the probability of transitioning from hidden state i to j with Tij , and the probability of emitting a high-amplitude bin in hidden state i with pi . For a data segment with M of N bins exceeding the quantization threshold, the likelihoods for the ${\mathcal{H}}_{0}$ and ${\mathcal{H}}_{1}$ hypotheses are respectively:

Equation (2)

Equation (3)

where y is a vector containing the observed (noisy) data. The maximal likelihood detection criterion $p(\boldsymbol{y}\vert {\mathcal{H}}_{1}) > p(\boldsymbol{y}\vert {\mathcal{H}}_{0})$ can therefore be rearranged to constrain the minimum number of high-amplitude bins (M) needed for CRES signal detection, as a linear function of the total number of bins in the data segment (N):

Equation (4)

Transition matrix elements are derived from the properties of the gas in the trap volume. In particular, T01, the probability of state transition ${\mathcal{H}}_{0}\to {\mathcal{H}}_{1}$ during a given time step is Γtb, where tb is the time bin length and Γ is the underlying event rate of electrons with kinetic energy within the Fourier bin. Likewise, the matrix element T10 corresponds to the probability that an event will end in a given time bin, which depends on the electron mean free time τ as $1-{\mathrm{e}}^{-{t}_{\mathrm{b}}/\tau }$. Other matrix elements are fixed by the normalization constraint ∑j Tij = 1.

The quantization threshold itself is optimized numerically. An excessively large quantization threshold would result in p0, p1 → 0, such that the sparse spectrogram in figure 2 (lower) would be completely empty. Likewise, an arbitrarily small quantization threshold would yield a completely filled (black) sparse spectrogram. The optimal quantization threshold minimizes the expected number of above-threshold bins required for detection (equation (4)).

2.1. Incorporating event structure

The first extension to this simplified application of the Viterbi algorithm to CRES signal reconstruction is to include information using the signal event topology to increase sensitivity to temporally short tracks. In particular, CRES electron tracks are usually not isolated; due to the inelastic scattering of electrons with the residual gas, a single electron can produce multiple continuous-phase tracks, punctuated by abrupt frequency discontinuities, which together comprise a single electron event (figure 2). Eventually, a scatter will sufficiently redirect the electron velocity vector to result in the escape of the electron from the magnetic trap, ending the event. Only the start frequency of the first track in a given event is relevant for energy spectroscopy.

A brief high-amplitude fluctuation is more likely to be a result of a genuine CRES event if there is a definitive, long-duration signal nearby in time-frequency space than if the same transient signal was found isolated within noise. Fully incorporating this additional information from the event topology would yield additional sensitivity to shorter tracks than otherwise possible.

Our formulation of the Viterbi algorithm can be extended to incorporate this scattering information by simply adding additional signal states, labeled ${\mathcal{H}}_{i}$, corresponding to signals within the ith frequency bin of the short-time Fourier transform, with indices starting from 1, and with ${\mathcal{H}}_{0}$ still representing the noise-only hypothesis. The transition matrix is then extended to (Nb + 1) × (Nb + 1), where Nb is the number of frequency bins, and contains the cross-section weighted energy loss distribution [16] of the gas composition within the detector volume. These scattering matrix elements Tij  (i, j ≠ 0) correspond to the probability that a track in frequency bin i will scatter to bin j within a given time step, which is a product of the scattering rate and the energy loss probability distribution function, for the relevant gas particle. Notably, the scattering matrix element Tij is to first-order only dependent on the difference between i and j, corresponding to the energy loss of the scatter. As a result, the block matrix with elements Tij  (i, j ≠ 0) is approximately upper triangular, since inelastic scattering in the lab frame overwhelmingly results in decreases in the kinetic energy of the incident electron (i > jTij = 0).

Analogous to the decision criteria given in equations (2)–(4), we can explicitly derive the Viterbi detection criterion for a short track candidate in bin i with Mi of Ni high-amplitude bins, given that there is a track in bin j immediately succeeding it. Provided that the second track is sufficiently long (equation (4)) so as to exclude the noise-only hypothesis, the detection criterion is found by comparing the relative probabilities of the hidden state sequences ${\mathcal{H}}_{0}\to {\mathcal{H}}_{i}\to {\mathcal{H}}_{j}\to {\mathcal{H}}_{0}$ and ${\mathcal{H}}_{0}\to {\mathcal{H}}_{j}\to {\mathcal{H}}_{0}$, that is, with and without the initial track candidate:

Equation (5)

where pi is the probability of emitting a high-power Fourier bin given hidden state i. The resulting detection criterion is of the same linear form as equation (4), though with an altered final term, depending on the transition matrix elements of the proposed scatter.

In most cases, the direct event rate to a given frequency bin is significantly smaller than the rate of scattering to that bin, given an already present electron event, yielding transition matrix elements more sensitive to scattered tracks. Detection in the vicinity of another long track therefore has a reduced threshold, yielding greater sensitivity to short tracks. This defines the detection decision regions illustrated in figure 3, in which the incorporation of event topology information expands the region of detectable tracks. Short duration or low signal-to-noise ratio (SNR) tracks are classified by the Viterbi algorithm as being due to noise-only.

Figure 3.

Figure 3. Illustration of detection regions for Viterbi decoding of CRES signals, as a function of the number of bins above threshold (M) and the total number of bins (N) in a track candidate. For fixed track length, the presence of a nearby event lowers the number of high-power bins M required for detection, defining the conditional region (light grey).

Standard image High-resolution image

Equation (5) also encapsulates the impact of different gas compositions, and their associated energy loss spectra, on signal detection via the transition matrix elements Tij  (i, j ≠ 0). Induction of the above derivation implies that equation (5) can be applied to any subsequent track in the event, consequently describing the Viterbi track detection condition for arbitrary event topologies.

This formulation of the HMM state structure assumes that at most one electron is present within the detector bandwidth at a given instant, which is not necessarily true. In Project 8 Phase II T2 data, the observed event rate of approximately 1 per hour [15], with each event lasting several milliseconds, implied that the probability of multi-electron coincidences was negligible. Extensions to the above state space for multiple electrons are possible, though these are unlikely to provide meaningful improvements in the detection efficiency given the additional computational burden for most applications of CRES.

2.2. Application to raw spectrograms

Finally, one can reconfigure the Viterbi algorithm to reconstruct CRES signals directly in raw, continuous data, as opposed to the previously considered binary quantized sparse spectrograms. This is done by replacing the emission matrix, which stores the probability that each state emits each possible symbol, with a function containing the desired probability distribution function for an on-the-fly computation of this probability. Outside a constant penalty from removing a lookup table, such an infinite extension does not affect the overall computational time complexity.

We consider the typical case in which the noise spectrum is additive white Gaussian, as a result of thermal noise. Given that the initial phases of the signals are unknown a priori, the voltage Fourier magnitudes for the noise-only and signal-plus-noise hypotheses are given by Rayleigh and Rician distributions [5], respectively.

Equation (6)

where I0 is the modified Bessel function of the first kind, σ2 is the thermal noise variance, y is the voltage Fourier magnitude, and ν is the expected voltage amplitude under ${\mathcal{H}}_{1}$. The Rayleigh distribution, fRice(y|0, σ), corresponds to the voltage amplitude distribution under the ${\mathcal{H}}_{0}$ hypothesis.

Using these continuous distributions for the emission probabilities, we can likewise express the Viterbi maximal-likelihood detection criterion for CRES signal detection within raw spectrograms:

Equation (7)

where N is the number of time bins in the track candidate and Ntr is the total number of tracks in the event.

The versatility of this algorithm, namely that the same algorithm can be reapplied on different compression levels of the data series, is appealing from an experimental standpoint. For instance, one could use the sparse spectrogram, one bit version of the algorithm on a first pass through the data, in which data throughput is the limiting factor, followed by a higher resolution search on the track candidates, using the same algorithm applied directly to the raw spectrogram data.

In general, such as for non-thermal noise spectra, the Viterbi detection criterion is a sequential likelihood ratio test given by:

Equation (8)

Equations (4), (5) and (7) are specific instances of the likelihood for Johnson–Nyquist noise. The likelihood sum is evaluated on the time series in sequence, resetting to zero after the inequality is satisfied.

3. Results

The application of the Viterbi algorithm to CRES signals is validated using both Monte Carlo (MC) studies of HMM generated data and Project 8 Phase II 83mKr data [17, 18]. An example Viterbi track and event reconstruction is illustrated in figure 2. In this three-track Monte Carlo CRES event, the Viterbi algorithm correctly identifies and clusters the CRES tracks, despite the background noise. Without loss of generality, the inelastic scattering energy loss was uniformly distributed within three frequency bins, as opposed to the full effective energy loss spectrum resulting from the particular detector gas composition.

Figure 4 illustrates the reconstruction performance of the Viterbi algorithm on simulated HMM MC data via the track duration distribution of first tracks in reconstructed events. The underlying MC truth track duration distribution (black) is exponentially distributed given the time-independence of scattering in a constant density gaseous volume. The resulting Viterbi reconstruction of MC HMM data is highly efficient for tracks longer than ≈1 ms, at the given signal and noise parameters. Since temporally short tracks, with few high-amplitude bins, cannot be definitively distinguished from noise fluctuations, lower detection efficiency of short tracks is an inevitable consequence of false event rejection.

Figure 4.

Figure 4. Track duration distributions of first reconstructed tracks in events from application of the Viterbi algorithm. Data generated using MC HMM data with Phase II detector parameters (Pe = 0.35 fW, TN = 135 K, τ = 0.5 ms), where TN is the system noise temperature. Illustrated and reconstructed for both binary quantized (sparse) and continuous (raw) Fourier data with 40.96 μs histogram time bins. MC reconstruction efficiency ɛ inlaid for each event class.

Standard image High-resolution image

To decrease the computational burden of MC simulations, the event rate parameter T01 is decoupled between event generation and reconstruction, such that the Viterbi decoding is applied to only event-dense data and not lengthy noise-only segments. Similarly, in Project 8, finding the most probable hidden state sequence may not align exactly with the experimental goal of tritium endpoint neutrino mass inference. Even if a data segment is strictly more likely the result of a CRES electron than not (equations (7) and (8)), specific experimental goals may require greater rejection of false alarms. If desired, a CRES experiment may use an unphysical T01 in reconstruction for greater false alarm rejection. Unphysical values for T01 are equivalent to the general Neyman–Pearson decision rule $\mathrm{ln}\,\mathcal{L} > \eta $, which maximizes the detection probability for a given false alarm rate, which can be tuned with appropriate choice of η (T01). Dependence of the false alarm rate on these Viterbi parameters lacks closed analytic form outside specific emission spectra [19].

For the sparse spectrograms, the minimal detectable track length, ten bins and two bins for the single-track and multi-track events respectively, is consistent with the Viterbi detection criteria at the simulated detector parameters. However, given the infinite domain of noise fluctuations, it is possible, though highly unlikely, that the raw spectrogram decision rule given by equation (7) can be satisfied within a single sample, so the ultimate limit on the quickest detection of CRES signals is ill-defined for raw data. Instead, Viterbi detection of CRES signals in noise yields a probability distribution function for the detection latency. For non-negligible SNRs (ν2/σ2 ≳ 5), equation (7) implies a median time required for detection:

Equation (9)

where h is the final term (braces) of equations (7) and (8), and τSNR = kB TN/Pe, where kB is Boltzmann's constant, and Pe is the radiated Larmor power. The shortest detectable track duration in this MC study was found to agree with the chosen threshold for Viterbi decoding, td = 0.10 ms. It should be noted that the median detection time is distinct from the mode of the track duration distribution, which scales with the mean free time of the electron in the gas. In addition, we note the displacement between the modes arising from differences in the reconstruction performance of single-track and multi-track events. No false events were found within the MC sample.

Direct application of the Viterbi algorithm to Phase II CRES data is challenging due to detector design choices which introduce additional signal complexity. In particular, electron-waveguide-mode interactions in Phase II can cause significant variations in the electron radiation rate, resulting in additional energy loss that varies on an event-by-event basis. To handle this, we apply a linear dechirp $({\text{e}}^{-\text{i}\alpha {t}^{2}/2})$, α = 4π × 108 s−2, removing most of the cyclotron-radiation-induced frequency slope from the spectrogram. Additional energy losses of signals from cyclotron radiation are handled by the transition matrix in the same manner as scattering energy losses. Significant variations in the event-by-event electron radiation rate are not expected to be a prominent feature of future phased-antenna array tritium endpoint neutrino mass CRES experiments, since these are limited to these small-volume waveguide detectors [4].

Despite the simplified signal assumptions, such as fixed signal amplitudes, the Viterbi algorithm is applied to 10 h of Phase II 83mKr data. The Viterbi algorithm was implemented in C++ using the Katydid [20] repository for standardized Project 8 file pre-processing. Code optimizations, such as reducing the Viterbi search space by neglecting the possibility of electron energy gains from scattering reduced both time and memory costs by approximately 50%. Performance is assessed in figure 5 via diagnostics of reconstructed tracks and events.

Figure 5.

Figure 5. Performance of Viterbi algorithm on Phase II 83mKr sparse spectrogram data. Distribution of first track durations in reconstructed events with 40.96 μs histogram time bins (left), and start frequency distribution of reconstructed events (right). 17.825 keV 83mKr line highlighted in red.

Standard image High-resolution image

Figure 5 (left) illustrates the track duration distribution of the first track in events, which approximately mimics the idealized reconstruction performance for short tracks in MC simulations (figure 4). The right-hand plot of figure 5, illustrating the initial frequency bin of the first track in reconstructed events, adds further support to the reconstruction performance, reproducing the detector response to the monoenergetic 17.825 keV K-32 line of 83mKr. In addition, since the Viterbi algorithm is applied to the frequency bandwidth 5–95 MHz, the absence of reconstructed events below 40 MHz is consistent with the chosen 1/100 day false alarm rate. Use of the full information of the raw spectrogram yields an approximate ×3 improvement in the detection efficiency compared to the sparse spectrogram Viterbi reconstruction, increasing the number of detected events from 976 to 2713 without additional apparent false alarms.

Besides yielding analytical limits on the minimum track duration required for CRES signal detection, this work demonstrates that the Viterbi algorithm itself is a feasible and efficient reconstruction algorithm for CRES events.

4. Conclusion

CRES is a modern technique, well-suited for high-precision energy spectroscopy measurements of light charged particles, and of interest for next-generation electron spectroscopy experiments. Despite the potential experimental utility of the technique, a systematic examination of the general detection of CRES-like signals had not been performed directly. We used a novel application of the Viterbi algorithm, a dynamic programming algorithm which explicitly searches all possible hidden state sequences, in order to characterize maximal-likelihood CRES signal detection. Using the Viterbi decision rules, we can determine the expected minimal track duration that can be detected, as a function of the effective signal power and the system noise temperature, under arbitrary data quantization schemes and detector gas compositions.

Acknowledgments

This material is based upon work supported by the following sources: the U.S. Department of Energy Office of Science, Office of Nuclear Physics, under Award No. DE-SC0020433 to Case Western Reserve University (CWRU), under Award No. DE-SC0011091 to the Massachusetts Institute of Technology (MIT), under the Early Career Research Program to Pacific Northwest National Laboratory (PNNL), a multiprogram national laboratory operated by Battelle for the U.S. Department of Energy under Contract No. DE-AC05-76RL01830, under Early Career Award No. DE-SC0019088 to Pennsylvania State University, under Award No. DE-FG02-97ER41020 to the University of Washington, and under Award No. DE-SC0012654 to Yale University; the National Science Foundation under Award No. PHY-1806251 to MIT. This work has been supported by the Cluster of Excellence 'Precision Physics, Fundamental Interactions, and Structure of Matter' (PRISMA + EXC 2118/1) funded by the German Research Foundation (DFG) within the German Excellence Strategy (Project ID 39083149); Laboratory Directed Research and Development (LDRD) 18-ERD-028 at Lawrence Livermore National Laboratory (LLNL), prepared by LLNL under Contract DE-AC52-07NA27344, LLNL-JRNL-830010; the LDRD Program at PNNL; the University of Washington Royalty Research Foundation; and support from Yale University. A portion of the research was performed using Research Computing at PNNL. The isotope(s) used in this research were supplied by the United States Department of Energy Office of Science by the Isotope Program in the Office of Nuclear Physics. We further acknowledge support from Yale University, the PRISMA + Cluster of Excellence at the University of Mainz, and the Karlsruhe Institute of Technology (KIT) Center Elementary Particle and Astroparticle Physics (KCETA).

Data availability statement

The data that support the findings of this study are available upon reasonable request from the authors.

Please wait… references are loading.
10.1088/1367-2630/ac66f6