shape_runtime/intrinsics/
fft.rs

1//! FFT (Fast Fourier Transform) intrinsics — full migration to typed marshal layer.
2//!
3//! Per the intrinsics-typed-CC migration's per-file table, all 5 surviving fft
4//! intrinsics (`fft`, `psd`, `dominant_frequency`, `bandpass`, `harmonics`)
5//! migrate to `register_typed_fn_N` typed entries via
6//! [`create_fft_intrinsics_module`]. The 6th legacy `intrinsic_ifft` was
7//! deleted as orphan-cleanup per supervisor sign-off relayed 2026-05-07
8//! (N3 → DELETE-NOW): zero stdlib/package consumers verified pre-deletion;
9//! N3 architectural decision (polymorphic-input dispatch via TypedObject
10//! FFT-result vs (real_arr, imag_arr) two-array form) deferred pending
11//! future consumer with similar polymorphic-input shape. Same precedent
12//! as scan.rs deletion at `663b63a`.
13//!
14//! Migrated entries take `Arc<Vec<f64>>` (series + kernel)
15//! and scalars (frequencies, sample_rate, num_harmonics); outputs project
16//! through `ConcreteReturn::ArrayF64` (psd, bandpass) or
17//! `TypedReturn::TypedObject(...)` (fft, dominant_frequency, harmonics).
18//!
19//! Provides FFT and related spectral analysis functions for:
20//! - Medical signal processing (ECG, EEG)
21//! - Power electronics (harmonic analysis)
22//! - Audio/vibration analysis
23//! - General frequency domain analysis
24
25use crate::marshal::{
26    register_typed_fn_1, register_typed_fn_2_full, register_typed_fn_4_full,
27};
28use crate::module_exports::{ModuleExports, ModuleParam};
29use crate::typed_module_exports::{ConcreteReturn, ConcreteType, TypedReturn};
30use rustfft::{FftPlanner, num_complex::Complex};
31use std::sync::Arc;
32
33// ───────────────────── Module factory (4 typed entries) ─────────────────────
34
35/// Create the fft intrinsics module with 4 typed-marshal entry points
36/// (`fft`, `psd`, `dominant_frequency`, `bandpass`, `harmonics`).
37/// `ifft` remains as legacy `IntrinsicFn` body in this module until the N3
38/// sub-decision (polymorphic input: TypedObject vs (real, imag) arrays)
39/// resolves.
40pub fn create_fft_intrinsics_module() -> ModuleExports {
41    let mut module = ModuleExports::new("std::core::intrinsics::fft");
42    module.description =
43        "FFT and frequency-domain analysis intrinsics (typed entries; ifft stays legacy pending N3 sub-decision)"
44            .to_string();
45
46    register_typed_fn_1::<_, Arc<Vec<f64>>>(
47        &mut module,
48        "__intrinsic_fft",
49        "Forward FFT: real series → { real, imag, magnitude, phase, frequencies, n }",
50        "series",
51        "Array<number>",
52        ConcreteType::TypedObject,
53        |series, _ctx| {
54            let data = series.as_slice();
55            let n = data.len();
56            if n == 0 {
57                return Ok(TypedReturn::TypedObject(empty_fft_pairs()));
58            }
59            let mut buffer: Vec<Complex<f64>> =
60                data.iter().map(|&x| Complex::new(x, 0.0)).collect();
61            let mut planner = FftPlanner::new();
62            let fft = planner.plan_fft_forward(n);
63            fft.process(&mut buffer);
64            let real: Vec<f64> = buffer.iter().map(|c| c.re).collect();
65            let imag: Vec<f64> = buffer.iter().map(|c| c.im).collect();
66            let magnitude: Vec<f64> = buffer.iter().map(|c| c.norm()).collect();
67            let phase: Vec<f64> = buffer.iter().map(|c| c.im.atan2(c.re)).collect();
68            let frequencies: Vec<f64> = (0..n).map(|i| i as f64 / n as f64).collect();
69            Ok(TypedReturn::TypedObject(vec![
70                ("real".to_string(), ConcreteReturn::ArrayF64(real)),
71                ("imag".to_string(), ConcreteReturn::ArrayF64(imag)),
72                ("magnitude".to_string(), ConcreteReturn::ArrayF64(magnitude)),
73                ("phase".to_string(), ConcreteReturn::ArrayF64(phase)),
74                (
75                    "frequencies".to_string(),
76                    ConcreteReturn::ArrayF64(frequencies),
77                ),
78                ("n".to_string(), ConcreteReturn::F64(n as f64)),
79            ]))
80        },
81    );
82
83    register_typed_fn_1::<_, Arc<Vec<f64>>>(
84        &mut module,
85        "__intrinsic_psd",
86        "Power spectral density: scaled |FFT(x)|^2",
87        "series",
88        "Array<number>",
89        ConcreteType::ArrayNumber,
90        |series, _ctx| {
91            let data = series.as_slice();
92            let n = data.len();
93            if n == 0 {
94                return Ok(TypedReturn::Concrete(ConcreteReturn::ArrayF64(vec![])));
95            }
96            let mut buffer: Vec<Complex<f64>> =
97                data.iter().map(|&x| Complex::new(x, 0.0)).collect();
98            let mut planner = FftPlanner::new();
99            let fft = planner.plan_fft_forward(n);
100            fft.process(&mut buffer);
101            let scale = 1.0 / n as f64;
102            let psd: Vec<f64> = buffer.iter().map(|c| c.norm_sqr() * scale).collect();
103            Ok(TypedReturn::Concrete(ConcreteReturn::ArrayF64(psd)))
104        },
105    );
106
107    register_typed_fn_2_full::<_, Arc<Vec<f64>>, f64>(
108        &mut module,
109        "__intrinsic_dominant_frequency",
110        "Dominant frequency: { frequency, magnitude, bin } of strongest spectrum bin",
111        [
112            ModuleParam {
113                name: "series".to_string(),
114                type_name: "Array<number>".to_string(),
115                required: true,
116                description: "Input series".to_string(),
117                ..Default::default()
118            },
119            ModuleParam {
120                name: "sample_rate".to_string(),
121                type_name: "number".to_string(),
122                required: false,
123                description: "Sample rate (default 1.0)".to_string(),
124                default_snippet: Some("1.0".to_string()),
125                ..Default::default()
126            },
127        ],
128        ConcreteType::TypedObject,
129        |series, sample_rate, _ctx| {
130            let data = series.as_slice();
131            let n = data.len();
132            if n == 0 {
133                return Err("dominant_frequency(): empty series".to_string());
134            }
135            let mut buffer: Vec<Complex<f64>> =
136                data.iter().map(|&x| Complex::new(x, 0.0)).collect();
137            let mut planner = FftPlanner::new();
138            let fft = planner.plan_fft_forward(n);
139            fft.process(&mut buffer);
140            let half_n = n / 2;
141            let (max_bin, max_mag) = buffer[1..=half_n]
142                .iter()
143                .enumerate()
144                .map(|(i, c)| (i + 1, c.norm()))
145                .max_by(|a, b| a.1.partial_cmp(&b.1).unwrap())
146                .unwrap_or((0, 0.0));
147            let frequency = max_bin as f64 * sample_rate / n as f64;
148            Ok(TypedReturn::TypedObject(vec![
149                ("frequency".to_string(), ConcreteReturn::F64(frequency)),
150                ("magnitude".to_string(), ConcreteReturn::F64(max_mag)),
151                ("bin".to_string(), ConcreteReturn::F64(max_bin as f64)),
152            ]))
153        },
154    );
155
156    register_typed_fn_4_full::<_, Arc<Vec<f64>>, f64, f64, f64>(
157        &mut module,
158        "__intrinsic_bandpass",
159        "Bandpass filter via FFT: zero frequencies outside [low_freq, high_freq]",
160        [
161            ModuleParam {
162                name: "series".to_string(),
163                type_name: "Array<number>".to_string(),
164                required: true,
165                description: "Input series".to_string(),
166                ..Default::default()
167            },
168            ModuleParam {
169                name: "low_freq".to_string(),
170                type_name: "number".to_string(),
171                required: true,
172                description: "Low cutoff frequency".to_string(),
173                ..Default::default()
174            },
175            ModuleParam {
176                name: "high_freq".to_string(),
177                type_name: "number".to_string(),
178                required: true,
179                description: "High cutoff frequency".to_string(),
180                ..Default::default()
181            },
182            ModuleParam {
183                name: "sample_rate".to_string(),
184                type_name: "number".to_string(),
185                required: false,
186                description: "Sample rate (default 1.0)".to_string(),
187                default_snippet: Some("1.0".to_string()),
188                ..Default::default()
189            },
190        ],
191        ConcreteType::ArrayNumber,
192        |series, low_freq, high_freq, sample_rate, _ctx| {
193            let data = series.as_slice();
194            let n = data.len();
195            if n == 0 {
196                return Ok(TypedReturn::Concrete(ConcreteReturn::ArrayF64(vec![])));
197            }
198            let mut buffer: Vec<Complex<f64>> =
199                data.iter().map(|&x| Complex::new(x, 0.0)).collect();
200            let mut planner = FftPlanner::new();
201            let fft = planner.plan_fft_forward(n);
202            fft.process(&mut buffer);
203            let freq_resolution = sample_rate / n as f64;
204            for (i, c) in buffer.iter_mut().enumerate() {
205                let freq = if i <= n / 2 {
206                    i as f64 * freq_resolution
207                } else {
208                    (n - i) as f64 * freq_resolution
209                };
210                if freq < low_freq || freq > high_freq {
211                    *c = Complex::new(0.0, 0.0);
212                }
213            }
214            let ifft = planner.plan_fft_inverse(n);
215            ifft.process(&mut buffer);
216            let scale = 1.0 / n as f64;
217            let result: Vec<f64> = buffer.iter().map(|c| c.re * scale).collect();
218            Ok(TypedReturn::Concrete(ConcreteReturn::ArrayF64(result)))
219        },
220    );
221
222    register_typed_fn_4_full::<_, Arc<Vec<f64>>, f64, i64, f64>(
223        &mut module,
224        "__intrinsic_harmonics",
225        "Harmonic analysis: extract harmonics of a fundamental frequency",
226        [
227            ModuleParam {
228                name: "series".to_string(),
229                type_name: "Array<number>".to_string(),
230                required: true,
231                description: "Input series".to_string(),
232                ..Default::default()
233            },
234            ModuleParam {
235                name: "fundamental_freq".to_string(),
236                type_name: "number".to_string(),
237                required: true,
238                description: "Fundamental frequency".to_string(),
239                ..Default::default()
240            },
241            ModuleParam {
242                name: "num_harmonics".to_string(),
243                type_name: "int".to_string(),
244                required: true,
245                description: "Number of harmonics to extract".to_string(),
246                ..Default::default()
247            },
248            ModuleParam {
249                name: "sample_rate".to_string(),
250                type_name: "number".to_string(),
251                required: false,
252                description: "Sample rate (default 1.0)".to_string(),
253                default_snippet: Some("1.0".to_string()),
254                ..Default::default()
255            },
256        ],
257        ConcreteType::TypedObject,
258        |series, fundamental, num_harmonics, sample_rate, _ctx| {
259            if num_harmonics < 0 {
260                return Err("harmonics(): num_harmonics must be non-negative".to_string());
261            }
262            let data = series.as_slice();
263            let n = data.len();
264            if n == 0 {
265                return Err("harmonics(): empty series".to_string());
266            }
267            let num_harmonics = num_harmonics as usize;
268            let mut buffer: Vec<Complex<f64>> =
269                data.iter().map(|&x| Complex::new(x, 0.0)).collect();
270            let mut planner = FftPlanner::new();
271            let fft = planner.plan_fft_forward(n);
272            fft.process(&mut buffer);
273            let freq_resolution = sample_rate / n as f64;
274            let mut harmonic_freqs = Vec::with_capacity(num_harmonics);
275            let mut magnitudes = Vec::with_capacity(num_harmonics);
276            let mut phases = Vec::with_capacity(num_harmonics);
277            for h in 1..=num_harmonics {
278                let target_freq = fundamental * h as f64;
279                let bin = (target_freq / freq_resolution).round() as usize;
280                if bin < n / 2 {
281                    let c = buffer[bin];
282                    harmonic_freqs.push(target_freq);
283                    magnitudes.push(c.norm());
284                    phases.push(c.im.atan2(c.re));
285                }
286            }
287            Ok(TypedReturn::TypedObject(vec![
288                (
289                    "harmonics".to_string(),
290                    ConcreteReturn::ArrayF64(harmonic_freqs),
291                ),
292                (
293                    "magnitudes".to_string(),
294                    ConcreteReturn::ArrayF64(magnitudes),
295                ),
296                ("phases".to_string(), ConcreteReturn::ArrayF64(phases)),
297            ]))
298        },
299    );
300
301    module
302}
303
304// ───────────────────── Empty-FFT helper ─────────────────────
305
306fn empty_fft_pairs() -> Vec<(String, ConcreteReturn)> {
307    vec![
308        ("real".to_string(), ConcreteReturn::ArrayF64(vec![])),
309        ("imag".to_string(), ConcreteReturn::ArrayF64(vec![])),
310        ("magnitude".to_string(), ConcreteReturn::ArrayF64(vec![])),
311        ("phase".to_string(), ConcreteReturn::ArrayF64(vec![])),
312        ("frequencies".to_string(), ConcreteReturn::ArrayF64(vec![])),
313        ("n".to_string(), ConcreteReturn::F64(0.0)),
314    ]
315}
316
shape_runtime/intrinsics/fft.rs

shape_runtime/intrinsics/
fft.rs
