microdsp 0.1.3

DSP algorithms and utilities. Performant, embedded friendly and no_std compatible.
Documentation 
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

#[derive(Copy, Clone)]
/// A key maximum, i.e an NSDF maximum that may or may not correspond
/// to the pitch period.
pub struct KeyMax {
    /// The index into the NSDF array corresponding to this maximum
    pub lag_index: usize,
    /// The NSDF value at `lag_index` for this maximum
    pub value_at_lag_index: f32,
    /// The NSDF value for this maximum, approximated using parabolic interpolation.
    pub value: f32,
    /// The lag, in samples, for this maximum, approximated using parabolic interpolation.
    /// Not necessarily an integer.
    pub lag: f32,
}

impl KeyMax {
    pub(crate) fn new() -> Self {
        KeyMax {
            lag_index: 0,
            value_at_lag_index: 0.0,
            value: 0.0,
            lag: 0.0,
        }
    }

    pub(crate) fn set(&mut self, nsdf: &[f32], lag_index: usize) {
        self.lag_index = lag_index;
        let value_at_lag_index = nsdf[lag_index];
        self.value_at_lag_index = value_at_lag_index;

        // Use parabolic interpolation to approximate
        // the true maximum using the left and right neighbors
        let left_index = 0.max(lag_index - 1);
        let right_index = (nsdf.len() - 1).min(lag_index + 1);
        let left = nsdf[left_index];
        let right = nsdf[right_index];

        // Compute coefficients of a parabola ax^2 + bx + c passing through
        // (-1, left), (0, max), (1, right)
        let a = 0.5 * (right - 2.0 * value_at_lag_index + left);
        let b = 0.5 * (right - left);
        let c = value_at_lag_index;
        // Find the x value where the derivative is zero, i.e where the parabola has its maximum
        let x_max = if a != 0.0 { -b / (2.0 * a) } else { 0.0 };
        let value = a * x_max * x_max + b * x_max + c;
        let lag = (lag_index as f32) + x_max;

        self.value = value;
        self.lag = lag;
    }
}

#[cfg(test)]
mod tests {
    use super::*;

    #[test]
    fn key_maximum_interpolation() {
        {
            let nsdf: [f32; 4] = [0.0, 0.0, 3.0, 0.0];
            let mut key_max = KeyMax::new();
            key_max.set(&nsdf, 2);
            assert!((key_max.lag - 2.0).abs() <= f32::EPSILON);
            assert!((key_max.value - 3.0).abs() <= f32::EPSILON);
        }

        {
            let nsdf: [f32; 3] = [-2.0, 0.0, -1.0];
            let mut key_max = KeyMax::new();
            key_max.set(&nsdf, 1);
            assert!((key_max.lag - 1.1666666_f32).abs() <= f32::EPSILON);
        }
    }
}