Skip to main content

qualia_core_db/audio/
istft.rs

1//! Real inverse STFT — resynthesise a time-domain signal from the two-sided
2//! complex frames produced by [`crate::audio::stft::forward_stft`].
3//!
4//! This is the cold-path inverse companion to the forward STFT: it reuses the
5//! already-built, caller-buffered WOLA resynthesiser in
6//! [`qualia_audio::features::framing::istft`] (weighted overlap-add normalised by
7//! the running Σ analysis·synthesis window product), so reconstruction is exact
8//! over the fully-overlapped interior regardless of whether the window+hop
9//! strictly satisfy COLA.
10//!
11//! Heap allocation is fine here: ISTFT runs at ingest / edit (cold path), never
12//! in the zero-heap U3 hot worklet. Native only — `qualia-audio` is a non-wasm
13//! dependency of `qualia-core-db`.
14
15use crate::audio::stft_bake::StftBakeError;
16
17/// Hann window coefficient `w[i] = 0.5*(1 - cos(2π·i/(N-1)))` for a length-`n`
18/// frame — identical to [`crate::audio::stft`]'s analysis window, so the same
19/// window drives both the analysis and synthesis args of the WOLA resynthesis.
20#[inline]
21fn hann(i: usize, n: usize) -> f32 {
22    if n <= 1 {
23        return 1.0;
24    }
25    0.5 * (1.0 - (core::f32::consts::TAU * i as f32 / (n - 1) as f32).cos())
26}
27
28/// Inverse STFT of the two-sided complex frames from
29/// [`crate::audio::stft::forward_stft`] back into a single time-domain signal.
30///
31/// - `spec`: `num_frames` frames, each a length-`frame_size` two-sided spectrum
32///   of `[re, im]` bins (exactly the [`crate::audio::stft::forward_stft`] output).
33/// - `frame_size`: `N`, a power of two; the analysis/synthesis frame length.
34/// - `hop`: sample advance between frames (`> 0`). A COLA Hann at `hop = N/2`
35///   gives exact interior reconstruction.
36///
37/// The frames are flattened into a `num_frames × 2·N` interleaved-complex buffer
38/// (`[re, im, …]`), a length-`N` Hann window is built for both the analysis and
39/// synthesis arguments, and the WOLA resynthesiser
40/// [`qualia_audio::features::framing::istft`] is invoked with caller-owned
41/// scratch/norm/out buffers (`out_len = (num_frames − 1)·hop + N`).
42///
43/// The core-db forward sign (`exp(-2πi·k·j/N)`, [`crate::audio::stft`]) matches
44/// `qualia-audio`'s inverse sign (`fft_radix2(_, false)` uses `e^{-j2πk/L}`), so
45/// no conjugation is required — the forward output feeds straight in.
46///
47/// Returns the reconstructed samples (length `out_len`), or an empty `Vec` when
48/// `spec` is empty. Errors map [`qualia_audio::types::AudioError`] →
49/// [`StftBakeError`].
50#[cfg(not(target_arch = "wasm32"))]
51pub fn inverse_stft(
52    spec: &[Vec<[f32; 2]>],
53    frame_size: usize,
54    hop: usize,
55) -> Result<Vec<f32>, StftBakeError> {
56    if !frame_size.is_power_of_two() || frame_size == 0 || hop == 0 {
57        return Err(StftBakeError::InvalidFrameCount);
58    }
59    let num_frames = spec.len();
60    if num_frames == 0 {
61        return Ok(Vec::new());
62    }
63    // Every frame must be a full two-sided spectrum of `frame_size` bins.
64    for frame in spec {
65        if frame.len() != frame_size {
66            return Err(StftBakeError::InvalidFrameCount);
67        }
68    }
69
70    let two_n = 2 * frame_size;
71    // Flatten frames into `num_frames × 2N` interleaved complex `[re, im, …]`.
72    // No conjugation: forward and inverse share the `e^{-j2πk/L}` sign convention.
73    let mut spectra = vec![0.0_f32; num_frames * two_n];
74    for (f, frame) in spec.iter().enumerate() {
75        let seg = &mut spectra[f * two_n..(f + 1) * two_n];
76        for (k, &[re, im]) in frame.iter().enumerate() {
77            seg[2 * k] = re;
78            seg[2 * k + 1] = im;
79        }
80    }
81
82    // Length-N Hann for both analysis and synthesis window args.
83    let window: Vec<f32> = (0..frame_size).map(|i| hann(i, frame_size)).collect();
84
85    let out_len = (num_frames - 1)
86        .checked_mul(hop)
87        .and_then(|v| v.checked_add(frame_size))
88        .ok_or(StftBakeError::InvalidFrameCount)?;
89
90    let mut ifft_scratch = vec![0.0_f32; frame_size];
91    let mut norm = vec![0.0_f32; out_len];
92    let mut out = vec![0.0_f32; out_len];
93
94    let written = qualia_audio::features::framing::istft(
95        &mut spectra,
96        frame_size,
97        hop,
98        &window,
99        &window,
100        &mut ifft_scratch,
101        &mut norm,
102        &mut out,
103    )
104    .map_err(map_audio_err)?;
105
106    out.truncate(written);
107    Ok(out)
108}
109
110/// Map a `qualia-audio` [`qualia_audio::types::AudioError`] onto the core-db
111/// sidecar [`StftBakeError`] surface.
112#[cfg(not(target_arch = "wasm32"))]
113#[inline]
114fn map_audio_err(e: qualia_audio::types::AudioError) -> StftBakeError {
115    use qualia_audio::types::AudioError;
116    match e {
117        AudioError::OutputBufferTooSmall | AudioError::WorkspaceTooSmall => {
118            StftBakeError::OutputTooSmall
119        }
120        _ => StftBakeError::InvalidFrameCount,
121    }
122}
123
124#[cfg(all(test, not(target_arch = "wasm32")))]
125mod tests {
126    use super::*;
127    use crate::audio::stft::forward_stft;
128    use core::f32::consts::TAU;
129
130    /// GOLDEN round trip: `forward_stft` → `inverse_stft` reconstructs the original
131    /// signal over the steady-state (fully-overlapped) interior with a COLA Hann at
132    /// hop = N/2 (mean-abs-error < 1e-2). Also proves no conjugate fix is needed.
133    #[test]
134    fn istft_round_trip_reconstructs_interior() {
135        const N: usize = 256;
136        const HOP: usize = N / 2; // COLA Hann
137        let fs = 16_000.0f32;
138        let len = 4096usize;
139        let signal: Vec<f32> = (0..len)
140            .map(|i| {
141                let t = i as f32 / fs;
142                0.6 * (TAU * 440.0 * t).sin() + 0.35 * (TAU * 1234.0 * t).cos()
143            })
144            .collect();
145
146        let spec = forward_stft(&signal, N, HOP).expect("forward stft");
147        let recon = inverse_stft(&spec, N, HOP).expect("inverse stft");
148
149        // Steady-state interior [N, out_len - N): every sample seen by full overlap.
150        let out_len = recon.len();
151        assert!(out_len >= 2 * N + 1, "not enough overlap for interior test");
152        let (lo, hi) = (N, out_len - N);
153        let mut mae = 0.0f64;
154        for k in lo..hi {
155            mae += (recon[k] - signal[k]).abs() as f64;
156        }
157        mae /= (hi - lo) as f64;
158        assert!(mae < 1e-2, "round-trip steady-state MAE {mae} exceeds 1e-2");
159    }
160
161    #[test]
162    fn istft_empty_spec_yields_empty() {
163        let empty: Vec<Vec<[f32; 2]>> = Vec::new();
164        let out = inverse_stft(&empty, 256, 128).expect("empty ok");
165        assert!(out.is_empty());
166    }
167
168    #[test]
169    fn istft_output_length_matches_formula() {
170        const N: usize = 64;
171        const HOP: usize = 32;
172        let signal = vec![0.2f32; N * 5];
173        let spec = forward_stft(&signal, N, HOP).expect("forward");
174        let num_frames = spec.len();
175        let recon = inverse_stft(&spec, N, HOP).expect("inverse");
176        assert_eq!(recon.len(), (num_frames - 1) * HOP + N);
177    }
178
179    #[test]
180    fn istft_rejects_non_power_of_two() {
181        let spec = vec![vec![[0.0f32; 2]; 48]];
182        assert_eq!(
183            inverse_stft(&spec, 48, 16),
184            Err(StftBakeError::InvalidFrameCount)
185        );
186    }
187
188    #[test]
189    fn istft_rejects_wrong_frame_width() {
190        // Frame length (10) disagrees with frame_size (256).
191        let spec = vec![vec![[0.0f32; 2]; 10]];
192        assert_eq!(
193            inverse_stft(&spec, 256, 128),
194            Err(StftBakeError::InvalidFrameCount)
195        );
196    }
197}