quantwave_core/indicators/
griffiths_dominant_cycle.rs1use crate::indicators::metadata::{IndicatorMetadata, ParamDef};
2use crate::traits::Next;
3use crate::indicators::high_pass::HighPass;
4use crate::indicators::super_smoother::SuperSmoother;
5use std::collections::VecDeque;
6use std::f64::consts::PI;
7
8#[derive(Debug, Clone)]
13pub struct GriffithsDominantCycle {
14 lb: usize,
15 ub: usize,
16 length: usize,
17 mu: f64,
18 hp: HighPass,
19 ss: SuperSmoother,
20 peak: f64,
21 signal_window: VecDeque<f64>,
22 coef: Vec<f64>,
23 prev_cycle: f64,
24}
25
26impl GriffithsDominantCycle {
27 pub fn new(lower_bound: usize, upper_bound: usize, length: usize) -> Self {
28 Self {
29 lb: lower_bound,
30 ub: upper_bound,
31 length,
32 mu: 1.0 / (length as f64),
33 hp: HighPass::new(upper_bound),
34 ss: SuperSmoother::new(lower_bound),
35 peak: 0.1,
36 signal_window: VecDeque::with_capacity(length + 1),
37 coef: vec![0.0; length + 1],
38 prev_cycle: (lower_bound + upper_bound) as f64 / 2.0,
39 }
40 }
41}
42
43impl Default for GriffithsDominantCycle {
44 fn default() -> Self {
45 Self::new(18, 40, 40)
46 }
47}
48
49impl Next<f64> for GriffithsDominantCycle {
50 type Output = f64;
51
52 fn next(&mut self, input: f64) -> Self::Output {
53 let hp_val = self.hp.next(input);
54 let lp_val = self.ss.next(hp_val);
55
56 self.peak *= 0.991;
57 if lp_val.abs() > self.peak {
58 self.peak = lp_val.abs();
59 }
60
61 let signal = if self.peak != 0.0 {
62 lp_val / self.peak
63 } else {
64 0.0
65 };
66
67 self.signal_window.push_front(signal);
68 if self.signal_window.len() > self.length {
69 self.signal_window.pop_back();
70 }
71
72 if self.signal_window.len() < self.length {
73 return self.prev_cycle;
74 }
75
76 let mut xx = vec![0.0; self.length + 1];
77 for (i, val) in xx.iter_mut().enumerate().skip(1).take(self.length) {
78 *val = self.signal_window[self.length - i];
79 }
80
81 let mut x_bar = 0.0;
82 for count in 1..=self.length {
83 x_bar += xx[self.length - count] * self.coef[count];
84 }
85
86 for count in 1..=self.length {
87 self.coef[count] += self.mu * (xx[self.length] - x_bar) * xx[self.length - count];
88 }
89
90 let mut max_pwr = 0.0;
92 let mut cycle = self.prev_cycle;
93
94 for period_idx in self.lb..=self.ub {
95 let period = period_idx as f64;
96 let mut real = 0.0;
97 let mut imag = 0.0;
98
99 for count in 1..=self.length {
100 let angle = 2.0 * PI * (count as f64) / period;
101 real += self.coef[count] * angle.cos();
102 imag += self.coef[count] * angle.sin();
103 }
104
105 let denom = (1.0 - real).powi(2) + imag.powi(2);
106 let pwr = 0.1 / denom;
107
108 if pwr > max_pwr {
109 max_pwr = pwr;
110 cycle = period;
111 }
112 }
113
114 if cycle > self.prev_cycle + 2.0 {
116 cycle = self.prev_cycle + 2.0;
117 } else if cycle < self.prev_cycle - 2.0 {
118 cycle = self.prev_cycle - 2.0;
119 }
120
121 self.prev_cycle = cycle;
122 cycle
123 }
124}
125
126pub const GRIFFITHS_DOMINANT_CYCLE_METADATA: IndicatorMetadata = IndicatorMetadata {
127 name: "GriffithsDominantCycle",
128 description: "Dominant cycle estimation using Griffiths adaptive spectral analysis.",
129 usage: "Use as a robust dominant cycle estimator less sensitive to amplitude changes than DFT-based methods, making it reliable across different market volatility regimes.",
130 keywords: &["cycle", "dominant-cycle", "ehlers", "dsp", "spectral"],
131 ehlers_summary: "The Griffiths method computes the dominant cycle by solving the real-roots of an autocorrelation polynomial. Adapted by Ehlers in Cycle Analytics for Traders, it remains stable even when market amplitude changes rapidly, unlike power-spectrum methods that can shift with volatility.",
132 params: &[
133 ParamDef {
134 name: "lower_bound",
135 default: "18",
136 description: "Lower period bound",
137 },
138 ParamDef {
139 name: "upper_bound",
140 default: "40",
141 description: "Upper period bound",
142 },
143 ParamDef {
144 name: "length",
145 default: "40",
146 description: "LMS filter length",
147 },
148 ],
149 formula_source: "https://github.com/lavs9/quantwave/blob/main/references/traderstipsreference/TRADERS’%20TIPS%20-%20JANUARY%202025.html",
150 formula_latex: r#"
151\[
152Pwr(Period) = \frac{0.1}{(1-Real)^2 + Imag^2}
153\]
154\[
155Real = \sum coef_i \cos(2\pi i / Period)
156\]
157\[
158Imag = \sum coef_i \sin(2\pi i / Period)
159\]
160"#,
161 gold_standard_file: "griffiths_dominant_cycle.json",
162 category: "Ehlers DSP",
163};
164
165#[cfg(test)]
166mod tests {
167 use super::*;
168 use crate::traits::Next;
169 use proptest::prelude::*;
170
171 #[test]
172 fn test_griffiths_dc_basic() {
173 let mut gdc = GriffithsDominantCycle::new(18, 40, 40);
174 for i in 0..200 {
175 let val = gdc.next((2.0 * PI * i as f64 / 30.0).sin());
177 if i > 150 {
178 assert!(val > 25.0 && val < 35.0);
180 }
181 }
182 }
183
184 proptest! {
185 #[test]
186 fn test_griffiths_dc_parity(
187 inputs in prop::collection::vec(1.0..100.0, 100..200),
188 ) {
189 let lb = 18;
190 let ub = 40;
191 let length = 40;
192 let mut gdc = GriffithsDominantCycle::new(lb, ub, length);
193 let streaming_results: Vec<f64> = inputs.iter().map(|&x| gdc.next(x)).collect();
194
195 let mut batch_results = Vec::with_capacity(inputs.len());
197 let mut hp = HighPass::new(ub);
198 let mut ss = SuperSmoother::new(lb);
199 let lp_vals: Vec<f64> = inputs.iter().map(|&x| ss.next(hp.next(x))).collect();
200
201 let mut peak = 0.1;
202 let mut signals = Vec::new();
203 let mut coef = vec![0.0; length + 1];
204 let mu = 1.0 / length as f64;
205 let mut prev_cycle = (lb + ub) as f64 / 2.0;
206
207 for (i, &lp_val) in lp_vals.iter().enumerate() {
208 peak *= 0.991;
209 if lp_val.abs() > peak {
210 peak = lp_val.abs();
211 }
212 let signal = if peak != 0.0 { lp_val / peak } else { 0.0 };
213 signals.push(signal);
214
215 if signals.len() < length {
216 batch_results.push(prev_cycle);
217 continue;
218 }
219
220 let mut xx = vec![0.0; length + 1];
221 for j in 1..=length {
222 xx[j] = signals[i - (length - j)];
223 }
224
225 let mut x_bar = 0.0;
226 for count in 1..=length {
227 x_bar += xx[length - count] * coef[count];
228 }
229
230 for count in 1..=length {
231 coef[count] += mu * (xx[length] - x_bar) * xx[length - count];
232 }
233
234 let mut max_pwr = 0.0;
235 let mut cycle = prev_cycle;
236
237 for period_idx in lb..=ub {
238 let period = period_idx as f64;
239 let mut real = 0.0;
240 let mut imag = 0.0;
241 for count in 1..=length {
242 let angle = 2.0 * PI * (count as f64) / period;
243 real += coef[count] * angle.cos();
244 imag += coef[count] * angle.sin();
245 }
246 let denom = (1.0 - real).powi(2) + imag.powi(2);
247 let pwr = 0.1 / denom;
248 if pwr > max_pwr {
249 max_pwr = pwr;
250 cycle = period;
251 }
252 }
253
254 if cycle > prev_cycle + 2.0 {
255 cycle = prev_cycle + 2.0;
256 } else if cycle < prev_cycle - 2.0 {
257 cycle = prev_cycle - 2.0;
258 }
259
260 prev_cycle = cycle;
261 batch_results.push(cycle);
262 }
263
264 for (s, b) in streaming_results.iter().zip(batch_results.iter()) {
265 approx::assert_relative_eq!(s, b, epsilon = 1e-10);
266 }
267 }
268 }
269}