Skip to main content

harmos_signal/spectrum/
mod.rs

1//! The signal in frequency: the transform that gets there, the readings taken
2//! off the bins, and the reconstruction back.
3//!
4//! [`Spectrum`] is the one value this anchor hands out, and it is a transparent
5//! record rather than a container: three public fields, no identity, no cache,
6//! no lock. Everything that *reads* a spectrum takes plain slices instead, so a
7//! caller holding its own series type calls every kernel here without ever
8//! constructing one.
9
10mod decibel;
11mod fourier;
12mod harmonic;
13mod indicator;
14mod taper;
15
16pub use decibel::{Decibel, Scale};
17pub use fourier::{Estimate, Welch, reconstruct, transform};
18pub use harmonic::{fundamental, harmonics, peaks, thd};
19pub use indicator::{band_energy, band_rms, centroid, occupied_bandwidth, papr, rolloff};
20pub use taper::{Band, Taper};
21
22/// One single-window transform: a bin step, and a complex amplitude per bin.
23///
24/// The bins run from DC upwards in steps of [`step`](Self::step), and the last
25/// one is Nyquist. Amplitude normalization is already applied, so
26/// [`magnitude`] of a bin is the amplitude of the component that landed there
27/// and not a quantity awaiting a scale factor.
28///
29/// The fields are public because this is a record and not an owner: a caller
30/// with a spectrum of its own builds one to call [`reconstruct`], and a caller
31/// that received one destructures it into whatever it keeps its data in.
32#[derive(Clone, Debug, PartialEq)]
33pub struct Spectrum {
34    /// The frequency distance between neighbouring bins, in Hz.
35    pub step: f64,
36    /// The real part of each bin.
37    pub real: Vec<f64>,
38    /// The imaginary part of each bin.
39    pub imag: Vec<f64>,
40}
41
42impl Spectrum {
43    /// The sampling rate the samples behind this spectrum were taken at, in Hz.
44    ///
45    /// Derived, not remembered: a single-sided spectrum of `bins` bins came
46    /// from `(bins - 1) * 2` samples, which is also the length
47    /// [`reconstruct`] answers with.
48    #[must_use]
49    pub fn rate(&self) -> f64 {
50        match self.real.len() {
51            0 | 1 => 0.0,
52            bins => self.step * ((bins - 1) * 2) as f64,
53        }
54    }
55}
56
57/// The magnitude of each bin: `hypot(real, imag)`.
58///
59/// Takes slices rather than a [`Spectrum`] on purpose — a caller whose spectrum
60/// lives in its own type passes its own two columns.
61#[must_use]
62pub fn magnitude(real: &[f64], imag: &[f64]) -> Vec<f64> {
63    real.iter()
64        .zip(imag)
65        .map(|(real, imag)| real.hypot(*imag))
66        .collect()
67}
68
69/// The phase of each bin in radians, `-π..=π`.
70#[must_use]
71pub fn phase(real: &[f64], imag: &[f64]) -> Vec<f64> {
72    real.iter()
73        .zip(imag)
74        .map(|(real, imag)| imag.atan2(*real))
75        .collect()
76}
77
78/// The frequency axis of `count` bins spaced `step` Hz apart, in Hz.
79#[must_use]
80pub fn bins(step: f64, count: usize) -> Vec<f64> {
81    (0..count).map(|index| index as f64 * step).collect()
82}
83
84#[cfg(test)]
85mod tests {
86    use std::f64::consts::FRAC_PI_2;
87
88    use super::*;
89
90    #[test]
91    fn magnitude_is_the_hypotenuse() {
92        assert_eq!(magnitude(&[3.0, 0.0], &[4.0, 1.0]), vec![5.0, 1.0]);
93    }
94
95    #[test]
96    fn phase_is_the_angle_off_the_real_axis() {
97        let phase = phase(&[1.0, 0.0, 0.0], &[0.0, 1.0, -1.0]);
98
99        assert!((phase[0] - 0.0).abs() < 1e-12);
100        assert!((phase[1] - FRAC_PI_2).abs() < 1e-12);
101        assert!((phase[2] + FRAC_PI_2).abs() < 1e-12);
102    }
103
104    #[test]
105    fn the_axis_starts_at_dc() {
106        assert_eq!(bins(2.5, 4), vec![0.0, 2.5, 5.0, 7.5]);
107        assert!(bins(2.5, 0).is_empty());
108    }
109
110    #[test]
111    fn the_rate_is_twice_the_last_bin() {
112        let spectrum = Spectrum {
113            step: 1.0,
114            real: vec![0.0; 501],
115            imag: vec![0.0; 501],
116        };
117
118        assert!((spectrum.rate() - 1000.0).abs() < 1e-9);
119    }
120
121    #[test]
122    fn a_spectrum_too_short_to_have_come_from_samples_has_no_rate() {
123        assert_eq!(
124            Spectrum {
125                step: 1.0,
126                real: Vec::new(),
127                imag: Vec::new()
128            }
129            .rate(),
130            0.0
131        );
132        assert_eq!(
133            Spectrum {
134                step: 1.0,
135                real: vec![1.0],
136                imag: vec![0.0]
137            }
138            .rate(),
139            0.0
140        );
141    }
142}