Skip to main content

quantwave_core/indicators/
butterworth.rs

1use crate::indicators::metadata::{IndicatorMetadata, ParamDef};
2use crate::traits::Next;
3use std::f64::consts::PI;
4
5/// 2-Pole Butterworth Filter
6///
7/// Based on John Ehlers' "Poles, Zeros, and Higher Order Filters".
8/// A second-order IIR filter with two poles and two zeros at Z=-1.
9#[derive(Debug, Clone)]
10pub struct Butterworth2 {
11    c1: f64,
12    b: f64,
13    aa: f64,
14    price_history: [f64; 2],
15    filt_history: [f64; 2],
16    count: usize,
17}
18
19impl Butterworth2 {
20    pub fn new(period: usize) -> Self {
21        let p = period as f64;
22        let a = (-1.414 * PI / p).exp();
23        let b = 2.0 * a * (1.414 * PI / p).cos();
24        let aa = a * a;
25        let c1 = (1.0 - b + aa) / 4.0;
26        Self {
27            c1,
28            b,
29            aa,
30            price_history: [0.0; 2],
31            filt_history: [0.0; 2],
32            count: 0,
33        }
34    }
35}
36
37impl Next<f64> for Butterworth2 {
38    type Output = f64;
39
40    fn next(&mut self, input: f64) -> Self::Output {
41        self.count += 1;
42        let res = if self.count < 3 {
43            input
44        } else {
45            self.b * self.filt_history[0] - self.aa * self.filt_history[1]
46                + self.c1 * (input + 2.0 * self.price_history[0] + self.price_history[1])
47        };
48
49        self.filt_history[1] = self.filt_history[0];
50        self.filt_history[0] = res;
51        self.price_history[1] = self.price_history[0];
52        self.price_history[0] = input;
53        res
54    }
55}
56
57/// 3-Pole Butterworth Filter
58///
59/// Based on John Ehlers' "Poles, Zeros, and Higher Order Filters".
60/// A third-order IIR filter with three poles and three zeros at Z=-1.
61#[derive(Debug, Clone)]
62pub struct Butterworth3 {
63    c1: f64,
64    b: f64,
65    c: f64,
66    bc: f64,
67    cc: f64,
68    price_history: [f64; 3],
69    filt_history: [f64; 3],
70    count: usize,
71}
72
73impl Butterworth3 {
74    pub fn new(period: usize) -> Self {
75        let p = period as f64;
76        let a = (-PI / p).exp();
77        let b = 2.0 * a * (1.738 * PI / p).cos();
78        let c = a * a;
79        let bc = b * c;
80        let cc = c * c;
81        let c1 = (1.0 - b + c) * (1.0 - c) / 8.0;
82        Self {
83            c1,
84            b,
85            c,
86            bc,
87            cc,
88            price_history: [0.0; 3],
89            filt_history: [0.0; 3],
90            count: 0,
91        }
92    }
93}
94
95impl Next<f64> for Butterworth3 {
96    type Output = f64;
97
98    fn next(&mut self, input: f64) -> Self::Output {
99        self.count += 1;
100        let res = if self.count < 4 {
101            input
102        } else {
103            (self.b + self.c) * self.filt_history[0] - (self.c + self.bc) * self.filt_history[1]
104                + self.cc * self.filt_history[2]
105                + self.c1
106                    * (input
107                        + 3.0 * self.price_history[0]
108                        + 3.0 * self.price_history[1]
109                        + self.price_history[2])
110        };
111
112        self.filt_history[2] = self.filt_history[1];
113        self.filt_history[1] = self.filt_history[0];
114        self.filt_history[0] = res;
115        self.price_history[2] = self.price_history[1];
116        self.price_history[1] = self.price_history[0];
117        self.price_history[0] = input;
118        res
119    }
120}
121
122pub const BUTTERWORTH2_METADATA: IndicatorMetadata = IndicatorMetadata {
123    name: "Butterworth2",
124    description: "2-pole Butterworth low-pass filter.",
125    usage: "Use to smooth price or intermediate indicator values with a flat passband and sharp rolloff. The 3-pole version provides steeper attenuation at the cost of marginally more lag.",
126    keywords: &["filter", "ehlers", "dsp", "smoothing", "low-pass"],
127    ehlers_summary: "Butterworth filters are maximally flat in the passband, introducing no ripple. Ehlers implements 2-pole and 3-pole Butterworth IIR designs in Cycle Analytics for Traders, noting that the SuperSmoother is actually a critically-damped 2-pole Butterworth variant.",
128    params: &[ParamDef {
129        name: "period",
130        default: "14",
131        description: "Critical period",
132    }],
133    formula_source: "https://github.com/lavs9/quantwave/blob/main/references/Ehlers%20Papers/Poles.pdf",
134    formula_latex: r#"
135\[
136a = \exp(-1.414\pi/P)
137\]
138\[
139b = 2a \cos(1.414\pi/P)
140\]
141\[
142f = bf_{t-1} - a^2f_{t-2} + \frac{1-b+a^2}{4}(g + 2g_{t-1} + g_{t-2})
143\]
144"#,
145    gold_standard_file: "butterworth2.json",
146    category: "Ehlers DSP",
147};
148
149pub const BUTTERWORTH3_METADATA: IndicatorMetadata = IndicatorMetadata {
150    name: "Butterworth3",
151    description: "3-pole Butterworth low-pass filter.",
152    usage: "Use to smooth price or intermediate indicator values with a flat passband and sharp rolloff. The 3-pole version provides steeper attenuation at the cost of marginally more lag.",
153    keywords: &["filter", "ehlers", "dsp", "smoothing", "low-pass"],
154    ehlers_summary: "Butterworth filters are maximally flat in the passband, introducing no ripple. Ehlers implements 2-pole and 3-pole Butterworth IIR designs in Cycle Analytics for Traders, noting that the SuperSmoother is actually a critically-damped 2-pole Butterworth variant.",
155    params: &[ParamDef {
156        name: "period",
157        default: "14",
158        description: "Critical period",
159    }],
160    formula_source: "https://github.com/lavs9/quantwave/blob/main/references/Ehlers%20Papers/Poles.pdf",
161    formula_latex: r#"
162\[
163a = \exp(-\pi/P)
164\]
165\[
166b = 2a \cos(1.738\pi/P)
167\]
168\[
169c = a^2
170\]
171\[
172f = (b+c)f_{t-1} - (c+bc)f_{t-2} + c^2f_{t-3} + \frac{(1-b+c)(1-c)}{8}(g + 3g_{t-1} + 3g_{t-2} + g_{t-3})
173\]
174"#,
175    gold_standard_file: "butterworth3.json",
176    category: "Ehlers DSP",
177};
178
179#[cfg(test)]
180mod tests {
181    use super::*;
182    use crate::traits::Next;
183    use proptest::prelude::*;
184
185    #[test]
186    fn test_butterworth_basic() {
187        let mut b2 = Butterworth2::new(14);
188        let mut b3 = Butterworth3::new(14);
189        for i in 0..20 {
190            let val = i as f64;
191            assert!(!b2.next(val).is_nan());
192            assert!(!b3.next(val).is_nan());
193        }
194    }
195
196    proptest! {
197        #[test]
198        fn test_butterworth2_parity(
199            inputs in prop::collection::vec(1.0..100.0, 10..100),
200        ) {
201            let p = 14;
202            let mut b2 = Butterworth2::new(p);
203            let streaming_results: Vec<f64> = inputs.iter().map(|&x| b2.next(x)).collect();
204
205            // Batch implementation
206            let mut batch_results = Vec::with_capacity(inputs.len());
207            let p_f = p as f64;
208            let a = (-1.414 * PI / p_f).exp();
209            let b = 2.0 * a * (1.414 * PI / p_f).cos();
210            let aa = a * a;
211            let c1 = (1.0 - b + aa) / 4.0;
212
213            let mut f_hist = [0.0; 2];
214            let mut g_hist = [0.0; 2];
215
216            for (i, &input) in inputs.iter().enumerate() {
217                let bar = i + 1;
218                let res = if bar < 3 {
219                    input
220                } else {
221                    b * f_hist[0] - aa * f_hist[1] + c1 * (input + 2.0 * g_hist[0] + g_hist[1])
222                };
223                f_hist[1] = f_hist[0];
224                f_hist[0] = res;
225                g_hist[1] = g_hist[0];
226                g_hist[0] = input;
227                batch_results.push(res);
228            }
229
230            for (s, b) in streaming_results.iter().zip(batch_results.iter()) {
231                approx::assert_relative_eq!(s, b, epsilon = 1e-10);
232            }
233        }
234
235        #[test]
236        fn test_butterworth3_parity(
237            inputs in prop::collection::vec(1.0..100.0, 10..100),
238        ) {
239            let p = 14;
240            let mut b3 = Butterworth3::new(p);
241            let streaming_results: Vec<f64> = inputs.iter().map(|&x| b3.next(x)).collect();
242
243            // Batch implementation
244            let mut batch_results = Vec::with_capacity(inputs.len());
245            let p_f = p as f64;
246            let a = (-PI / p_f).exp();
247            let b = 2.0 * a * (1.738 * PI / p_f).cos();
248            let c = a * a;
249            let bc = b * c;
250            let cc = c * c;
251            let c1 = (1.0 - b + c) * (1.0 - c) / 8.0;
252
253            let mut f_hist = [0.0; 3];
254            let mut g_hist = [0.0; 3];
255
256            for (i, &input) in inputs.iter().enumerate() {
257                let bar = i + 1;
258                let res = if bar < 4 {
259                    input
260                } else {
261                    (b + c) * f_hist[0] - (c + bc) * f_hist[1] + cc * f_hist[2]
262                        + c1 * (input + 3.0 * g_hist[0] + 3.0 * g_hist[1] + g_hist[2])
263                };
264                f_hist[2] = f_hist[1];
265                f_hist[1] = f_hist[0];
266                f_hist[0] = res;
267                g_hist[2] = g_hist[1];
268                g_hist[1] = g_hist[0];
269                g_hist[0] = input;
270                batch_results.push(res);
271            }
272
273            for (s, b) in streaming_results.iter().zip(batch_results.iter()) {
274                approx::assert_relative_eq!(s, b, epsilon = 1e-10);
275            }
276        }
277    }
278}