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
use crate::{Echo, SuperSmoother, View};
#[derive(Debug, Clone)]
pub struct RoofingFilter<V> {
view: V,
super_smoother: SuperSmoother<Echo>,
window_len: usize,
i: usize,
alpha_1: f64,
val_1: f64,
val_2: f64,
hp_1: f64,
hp_2: f64,
}
impl<V> RoofingFilter<V>
where
V: View,
{
#[inline]
pub fn new(view: V, window_len: usize, super_smoother_len: usize) -> Self {
let alpha_1 = ((4.4422 / window_len as f64).cos() + (4.4422 / window_len as f64).sin()
- 1.0)
/ (4.4422 / window_len as f64).cos();
RoofingFilter {
view,
super_smoother: SuperSmoother::new(Echo::new(), super_smoother_len),
window_len,
i: 0,
alpha_1,
val_1: 0.0,
val_2: 0.0,
hp_1: 0.0,
hp_2: 0.0,
}
}
}
impl<V> View for RoofingFilter<V>
where
V: View,
{
fn update(&mut self, val: f64) {
self.view.update(val);
let val = self.view.last();
let hp = (1.0 - self.alpha_1 / 2.0).powi(2) * (val - 2.0 * self.val_1 + self.val_2)
+ 2.0 * (1.0 - self.alpha_1) * self.hp_1
- (1.0 - self.alpha_1).powi(2) * self.hp_2;
self.hp_2 = self.hp_1;
self.hp_1 = hp;
self.val_2 = self.val_1;
self.val_1 = val;
if self.i > self.window_len {
self.super_smoother.update(hp);
}
self.i += 1;
}
#[inline(always)]
fn last(&self) -> f64 {
if self.i < self.window_len * 2 {
0.0
} else {
self.super_smoother.last()
}
}
}
#[cfg(test)]
mod tests {
use super::*;
use crate::plot::plot_values;
use crate::test_data::TEST_DATA;
use crate::Echo;
#[test]
fn roofing_filter_plot() {
let mut rf = RoofingFilter::new(Echo::new(), 48, 10);
let mut out: Vec<f64> = Vec::new();
for v in &TEST_DATA {
rf.update(*v);
out.push(rf.last());
}
let filename = "img/roofing_filter.png";
plot_values(out, filename).unwrap();
}
}