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
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113

#[derive(Clone, Copy, PartialEq, Debug)]
pub enum Peak {
    Low,
    High,
}

pub struct PeaksDetector {
    threshold: f64,
    influence: f64,
    window: Vec<f64>,
}

impl PeaksDetector {
    pub fn new(lag: usize, threshold: f64, influence: f64) -> PeaksDetector {
        PeaksDetector {
            threshold,
            influence,
            window: Vec::with_capacity(lag),
        }
    }

    pub fn signal(&mut self, value: f64) -> Option<Peak> {
        if self.window.len() < self.window.capacity() {
            self.window.push(value);
            None
        } else if let (Some((mean, stddev)), Some(&window_last)) = (self.stats(), self.window.last()) {
            self.window.remove(0);
            if (value - mean).abs() > (self.threshold * stddev) {
                let next_value =
                    (value * self.influence) + ((1. - self.influence) * window_last);
                self.window.push(next_value);
                Some(if value > mean { Peak::High } else { Peak::Low })
            } else {
                self.window.push(value);
                None
            }
        } else {
            None
        }
    }

    pub fn stats(&self) -> Option<(f64, f64)> {
        if self.window.is_empty() {
            None
        } else {
            let window_len = self.window.len() as f64;
            let mean = self.window.iter().fold(0., |a, v| a + v) / window_len;
            let sq_sum = self.window.iter().fold(0., |a, v| a + ((v - mean) * (v - mean)));
            let stddev = (sq_sum / window_len).sqrt();
            Some((mean, stddev))
        }
    }
}

pub struct PeaksIter<I, F> {
    source: I,
    signal: F,
    detector: PeaksDetector,
}

pub trait PeaksFilter<I> where I: Iterator {
    fn peaks<F>(self, detector: PeaksDetector, signal: F) -> PeaksIter<I, F>
        where F: FnMut(&I::Item) -> f64;
}

impl<I> PeaksFilter<I> for I where I: Iterator {
    fn peaks<F>(self, detector: PeaksDetector, signal: F) -> PeaksIter<I, F>
        where F: FnMut(&I::Item) -> f64
    {
        PeaksIter { source: self, signal, detector, }
    }
}

impl<I, F> Iterator for PeaksIter<I, F> where I: Iterator, F: FnMut(&I::Item) -> f64 {
    type Item = (I::Item, Peak);

    fn next(&mut self) -> Option<Self::Item> {
        while let Some(item) = self.source.next() {
            let value = (self.signal)(&item);
            if let Some(peak) = self.detector.signal(value) {
                return Some((item, peak));
            }
        }
        None
    }
}

#[cfg(test)]
mod tests {
    use super::{Peak, PeaksDetector, PeaksFilter};

    #[test]
    fn sample_data() {
        let input = vec![
            1.0, 1.0, 1.1, 1.0, 0.9, 1.0, 1.0, 1.1, 1.0, 0.9, 1.0, 1.1, 1.0, 1.0, 0.9, 1.0, 1.0, 1.1, 1.0,
            1.0, 1.0, 1.0, 1.1, 0.9, 1.0, 1.1, 1.0, 1.0, 0.9, 1.0, 1.1, 1.0, 1.0, 1.1, 1.0, 0.8, 0.9, 1.0,
            1.2, 0.9, 1.0, 1.0, 1.1, 1.2, 1.0, 1.5, 1.0, 3.0, 2.0, 5.0, 3.0, 2.0, 1.0, 1.0, 1.0, 0.9, 1.0,
            1.0, 3.0, 2.6, 4.0, 3.0, 3.2, 2.0, 1.0, 1.0, 0.8, 4.0, 4.0, 2.0, 2.5, 1.0, 1.0, 1.0
        ];
        let output: Vec<_> = input
            .into_iter()
            .enumerate()
            .peaks(PeaksDetector::new(30, 5.0, 0.0), |e| e.1)
            .map(|((i, _), p)| (i, p))
            .collect();
        assert_eq!(output, vec![
            (45, Peak::High), (47, Peak::High), (48, Peak::High), (49, Peak::High),
            (50, Peak::High), (51, Peak::High), (58, Peak::High), (59, Peak::High),
            (60, Peak::High), (61, Peak::High), (62, Peak::High), (63, Peak::High),
            (67, Peak::High), (68, Peak::High), (69, Peak::High), (70, Peak::High),
        ]);
    }
}