1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
use crate::float_trait::Float;
use crate::sorted_array::SortedArray;
use conv::{ConvAsUtil, ConvUtil, RoundToNearest};
use enum_dispatch::enum_dispatch;
use itertools::Itertools;
use schemars::JsonSchema;
use serde::{Deserialize, Serialize};
use std::fmt::Debug;
#[enum_dispatch]
trait NyquistFreqTrait: Send + Sync + Clone + Debug {
fn nyquist_freq<T: Float>(&self, t: &[T]) -> T;
}
#[enum_dispatch(NyquistFreqTrait)]
#[derive(Clone, Debug, Serialize, Deserialize, JsonSchema)]
#[non_exhaustive]
pub enum NyquistFreq {
Average(AverageNyquistFreq),
Median(MedianNyquistFreq),
Quantile(QuantileNyquistFreq),
}
#[derive(Clone, Debug, Serialize, Deserialize, JsonSchema)]
#[serde(rename = "Average")]
pub struct AverageNyquistFreq;
impl NyquistFreqTrait for AverageNyquistFreq {
fn nyquist_freq<T: Float>(&self, t: &[T]) -> T {
let n = t.len();
T::PI() * (n - 1).value_as().unwrap() / (t[n - 1] - t[0])
}
}
fn diff<T: Float>(x: &[T]) -> Vec<T> {
x.iter().tuple_windows().map(|(&a, &b)| b - a).collect()
}
#[derive(Clone, Debug, Serialize, Deserialize, JsonSchema)]
#[serde(rename = "Median")]
pub struct MedianNyquistFreq;
impl NyquistFreqTrait for MedianNyquistFreq {
fn nyquist_freq<T: Float>(&self, t: &[T]) -> T {
let sorted_dt: SortedArray<_> = diff(t).into();
let dt = sorted_dt.median();
T::PI() / dt
}
}
#[derive(Clone, Debug, Serialize, Deserialize, JsonSchema)]
#[serde(rename = "Quantile")]
pub struct QuantileNyquistFreq {
pub quantile: f32,
}
impl NyquistFreqTrait for QuantileNyquistFreq {
fn nyquist_freq<T: Float>(&self, t: &[T]) -> T {
let sorted_dt: SortedArray<_> = diff(t).into();
let dt = sorted_dt.ppf(self.quantile);
T::PI() / dt
}
}
#[derive(Clone, Debug)]
pub struct FreqGrid<T> {
pub step: T,
pub size: usize,
}
impl<T> FreqGrid<T>
where
T: Float,
{
pub fn from_t(t: &[T], resolution: f32, max_freq_factor: f32, nyquist: NyquistFreq) -> Self {
assert!(resolution.is_sign_positive() && resolution.is_finite());
let sizef = t.len().value_as::<T>().unwrap();
let duration = t[t.len() - 1] - t[0];
let step = T::two() * T::PI() * (sizef - T::one())
/ (sizef * resolution.value_as::<T>().unwrap() * duration);
let max_freq = nyquist.nyquist_freq(t) * max_freq_factor.value_as::<T>().unwrap();
let size = (max_freq / step).approx_by::<RoundToNearest>().unwrap();
Self { step, size }
}
}