diff --git a/Cargo.lock b/Cargo.lock index 1380c7b..d1e00e0 100644 --- a/Cargo.lock +++ b/Cargo.lock @@ -348,8 +348,7 @@ dependencies = [ "criterion", "libm", "ndarray", - "ndarray-stats", - "rand 0.8.3", + "rand", "serde", ] @@ -426,17 +425,6 @@ dependencies = [ "version_check", ] -[[package]] -name = "getrandom" -version = "0.1.16" -source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "8fc3cb4d91f53b50155bdcfd23f6a4c39ae1969c2ae85982b135750cccaf5fce" -dependencies = [ - "cfg-if 1.0.0", - "libc", - "wasi 0.9.0+wasi-snapshot-preview1", -] - [[package]] name = "getrandom" version = "0.2.2" @@ -445,7 +433,7 @@ checksum = "c9495705279e7140bf035dde1f6e750c162df8b625267cd52cc44e0b156732c8" dependencies = [ "cfg-if 1.0.0", "libc", - "wasi 0.10.2+wasi-snapshot-preview1", + "wasi", ] [[package]] @@ -627,30 +615,6 @@ dependencies = [ "rawpointer", ] -[[package]] -name = "ndarray-stats" -version = "0.4.0" -source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "22c95a780960082c5746f6bf0ab22d4a3b8cee72bf580acfe9f1e10bc5ea8152" -dependencies = [ - "indexmap", - "itertools 0.9.0", - "ndarray", - "noisy_float", - "num-integer", - "num-traits", - "rand 0.7.3", -] - -[[package]] -name = "noisy_float" -version = "0.1.13" -source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "a14c16cde392a1cd18084ffd8348cb8937525130e62f0478d72dcc683698809d" -dependencies = [ - "num-traits", -] - [[package]] name = "num-complex" version = "0.3.1" @@ -775,19 +739,6 @@ version = "0.2.2" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "e2a38df5b15c8d5c7e8654189744d8e396bddc18ad48041a500ce52d6948941f" -[[package]] -name = "rand" -version = "0.7.3" -source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "6a6b1679d49b24bbfe0c803429aa1874472f50d9b363131f0e89fc356b544d03" -dependencies = [ - "getrandom 0.1.16", - "libc", - "rand_chacha 0.2.2", - "rand_core 0.5.1", - "rand_hc 0.2.0", -] - [[package]] name = "rand" version = "0.8.3" @@ -795,19 +746,9 @@ source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "0ef9e7e66b4468674bfcb0c81af8b7fa0bb154fa9f28eb840da5c447baeb8d7e" dependencies = [ "libc", - "rand_chacha 0.3.0", - "rand_core 0.6.1", - "rand_hc 0.3.0", -] - -[[package]] -name = "rand_chacha" -version = "0.2.2" -source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "f4c8ed856279c9737206bf725bf36935d8666ead7aa69b52be55af369d193402" -dependencies = [ - "ppv-lite86", - "rand_core 0.5.1", + "rand_chacha", + "rand_core", + "rand_hc", ] [[package]] @@ -817,16 +758,7 @@ source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "e12735cf05c9e10bf21534da50a147b924d555dc7a547c42e6bb2d5b6017ae0d" dependencies = [ "ppv-lite86", - "rand_core 0.6.1", -] - -[[package]] -name = "rand_core" -version = "0.5.1" -source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "90bde5296fc891b0cef12a6d03ddccc162ce7b2aff54160af9338f8d40df6d19" -dependencies = [ - "getrandom 0.1.16", + "rand_core", ] [[package]] @@ -835,16 +767,7 @@ version = "0.6.1" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "c026d7df8b298d90ccbbc5190bd04d85e159eaf5576caeacf8741da93ccbd2e5" dependencies = [ - "getrandom 0.2.2", -] - -[[package]] -name = "rand_hc" -version = "0.2.0" -source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "ca3129af7b92a17112d59ad498c6f81eaf463253766b90396d39ea7a39d6613c" -dependencies = [ - "rand_core 0.5.1", + "getrandom", ] [[package]] @@ -853,7 +776,7 @@ version = "0.3.0" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "3190ef7066a446f2e7f42e239d161e905420ccab01eb967c9eb27d21b2322a73" dependencies = [ - "rand_core 0.6.1", + "rand_core", ] [[package]] @@ -1183,12 +1106,6 @@ dependencies = [ "winapi-util", ] -[[package]] -name = "wasi" -version = "0.9.0+wasi-snapshot-preview1" -source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "cccddf32554fecc6acb585f82a32a72e28b48f8c4c1883ddfeeeaa96f7d8e519" - [[package]] name = "wasi" version = "0.10.2+wasi-snapshot-preview1" diff --git a/dsp/Cargo.toml b/dsp/Cargo.toml index 7e30348..de96d81 100644 --- a/dsp/Cargo.toml +++ b/dsp/Cargo.toml @@ -12,7 +12,6 @@ serde = { version = "1.0", features = ["derive"], default-features = false } criterion = "0.3" rand = "0.8" ndarray = "0.14" -ndarray-stats = "0.4" [[bench]] name = "micro" diff --git a/dsp/src/lockin.rs b/dsp/src/lockin.rs index b75038e..0f6b78f 100644 --- a/dsp/src/lockin.rs +++ b/dsp/src/lockin.rs @@ -8,6 +8,7 @@ pub struct Lockin { } impl Lockin { + /// Create a new Lockin with given IIR coefficients. pub fn new(ba: &iir_int::IIRState) -> Self { let mut iir = iir_int::IIR::default(); iir.ba.0.copy_from_slice(&ba.0); @@ -17,9 +18,10 @@ impl Lockin { } } - pub fn update(&mut self, signal: i32, phase: i32) -> Complex { + /// Update the lockin with a sample taken at a given phase. + pub fn update(&mut self, sample: i32, phase: i32) -> Complex { // Get the LO signal for demodulation. - let m = Complex::from_angle(phase); + let lo = Complex::from_angle(phase); // Mix with the LO signal, filter with the IIR lowpass, // return IQ (in-phase and quadrature) data. @@ -27,15 +29,18 @@ impl Lockin { Complex( self.iir.update( &mut self.iir_state[0], - ((signal as i64 * m.0 as i64) >> 32) as _, + ((sample as i64 * lo.0 as i64) >> 32) as _, ), self.iir.update( &mut self.iir_state[1], - ((signal as i64 * m.1 as i64) >> 32) as _, + ((sample as i64 * lo.1 as i64) >> 32) as _, ), ) } + /// Feed an iterator into the Lockin and return the latest I/Q data. + /// Initial stample phase and frequency (phase increment between samples) + /// are supplied. pub fn feed>( &mut self, signal: I, @@ -46,131 +51,10 @@ impl Lockin { signal .into_iter() - .map(|s| { + .map(|sample| { phase = phase.wrapping_add(frequency); - self.update(s, phase) + self.update(sample, phase) }) .last() } } - -#[cfg(test)] -mod test { - use crate::{ - iir_int::IIRState, - lockin::Lockin, - rpll::RPLL, - testing::{isclose, max_error}, - Complex, - }; - - use std::f64::consts::PI; - use std::vec::Vec; - - /// ADC full scale in machine units (16 bit signed). - const ADC_SCALE: f64 = ((1 << 15) - 1) as _; - - struct PllLockin { - harmonic: i32, - phase: i32, - lockin: Lockin, - } - - impl PllLockin { - pub fn new(harmonic: i32, phase: i32, iir: &IIRState) -> Self { - PllLockin { - harmonic, - phase, - lockin: Lockin::new(iir), - } - } - - pub fn update( - &mut self, - input: Vec, - phase: i32, - frequency: i32, - ) -> Complex { - let sample_frequency = frequency.wrapping_mul(self.harmonic); - let mut sample_phase = - self.phase.wrapping_add(phase.wrapping_mul(self.harmonic)); - input - .iter() - .map(|&s| { - let input = (s as i32) << 16; - let signal = - self.lockin.update(input, sample_phase.wrapping_neg()); - sample_phase = sample_phase.wrapping_add(sample_frequency); - signal - }) - .last() - .unwrap_or(Complex::default()) - } - } - - /// Single-frequency sinusoid. - #[derive(Copy, Clone)] - struct Tone { - // Frequency (in Hz). - frequency: f64, - // Phase offset (in radians). - phase: f64, - // Amplitude in dBFS (decibels relative to full-scale). - // A 16-bit ADC has a minimum dBFS for each sample of -90. - amplitude_dbfs: f64, - } - - /// Convert dBFS to a linear ratio. - fn linear(dbfs: f64) -> f64 { - 10f64.powf(dbfs / 20.) - } - - impl Tone { - fn eval(&self, time: f64) -> f64 { - linear(self.amplitude_dbfs) - * (self.phase + self.frequency * time).cos() - } - } - - /// Generate a full batch of samples with size `sample_buffer_size` starting at `time_offset`. - fn sample_tones( - tones: &Vec, - time_offset: f64, - sample_buffer_size: u32, - ) -> Vec { - (0..sample_buffer_size) - .map(|i| { - let time = 2. * PI * (time_offset + i as f64); - let x: f64 = tones.iter().map(|t| t.eval(time)).sum(); - assert!(-1. < x && x < 1.); - (x * ADC_SCALE) as i16 - }) - .collect() - } - - /// Total maximum noise amplitude of the input signal after 2nd order lowpass filter. - /// Constructive interference is assumed. - /// - /// # Args - /// * `tones` - Noise sources at the ADC input. - /// * `frequency` - Frequency of the signal of interest. - /// * `corner` - Low-pass filter 3dB corner cutoff frequency. - /// - /// # Returns - /// Upper bound of the total amplitude of all noise sources in linear units full scale. - fn sampled_noise_amplitude( - tones: &Vec, - frequency: f64, - corner: f64, - ) -> f64 { - tones - .iter() - .map(|t| { - let df = (t.frequency - frequency) / corner; - // Assuming a 2nd order lowpass filter: 40dB/decade. - linear(t.amplitude_dbfs - 40. * df.abs().max(1.).log10()) - }) - .sum::() - .max(1. / ADC_SCALE / 2.) // 1/2 LSB from quantization - } -} diff --git a/dsp/src/rpll.rs b/dsp/src/rpll.rs index 0006831..9759cf9 100644 --- a/dsp/src/rpll.rs +++ b/dsp/src/rpll.rs @@ -84,7 +84,6 @@ impl RPLL { mod test { use super::RPLL; use ndarray::prelude::*; - use ndarray_stats::QuantileExt; use rand::{prelude::*, rngs::StdRng}; use std::vec::Vec; @@ -254,7 +253,7 @@ mod test { h.shift_frequency = 10; h.shift_phase = 9; - h.measure(1 << 16, [5e-7, 3e-2, 3e-2, 2e-2]); + h.measure(1 << 16, [5e-7, 3e-2, 2e-5, 2e-2]); } #[test]