[Feature]: Implement Inverse Short-Time Fourier Transform (iSTFT)

[gabyfle]

Initial Checklist

• I have searched the existing issues to ensure this feature has not already been requested.

Problem Statement & User Value

As a DSP Engineer, Machine Learning Engineer, I want to reconstruct a time-domain signal from its Short-Time Fourier Transform (STFT) representation, so that I can implement and apply spectral processing algorithms such as filtering, time-stretching, pitch-shifting, and audio synthesis.

Acceptance Criteria (Definition of Done)

API Definition

• A new public function SoundML.Transform.istft is available in the library's Transform module.
• The function signature is clearly defined in the interface (.mli) file, adhering to the design proposed in the "Proposed Solution" section.
• Comprehensive documentation (using odoc) is provided for the function. The documentation must detail each parameter, the return value, potential exceptions, and explicitly state the conditions required for perfect reconstruction.

Core Functionality & Correctness

• The function correctly implements a least-squares-based overlap-add reconstruction algorithm, consistent with established libraries like librosa and SciPy.
• Perfect Reconstruction Test: For any input signal $x$, the round-trip operation $istft(stft(x))$ must yield a reconstructed signal $x'$ that is numerically identical to $x$ (within a small tolerance for floating-point error, e.g., $1\times10^{-7}$). This test must pass under the following conditions:
  • The same n_fft, hop_length, and window parameters are used for both stft andistft.
  • The chosen window and hop_length satisfy the Constant Overlap-Add (COLA) constraint.
  • Boundary handling (e.g., center=true) is used consistently for both transforms.
• The function correctly handles both one-sided (for real-valued input signals) and two-sided STFT matrices.

Parameter Handling & Validation

• The function correctly infers n_fft and win_length from the input STFT matrix shape when these parameters are not explicitly provided, as documented by librosa.
• The function supports a variety of standard window types (e.g., Hann, Hamming, Rectangular) passed via a variant type.
• The function correctly handles a center flag to ensure proper alignment with the padding behavior of the forward stft function, which is critical for avoiding boundary artifacts.
• The function correctly implements an optional length parameter to trim or zero-pad the output signal to a user-specified length, enabling exact signal reconstruction.
• The function raises appropriate, clearly-defined exceptions for invalid inputs, such as a malformed STFT matrix or inconsistent parameter combinations (e.g., win_length > n_fft).

Utility Functions

• A public helper function SoundML.Utils.is_cola : window -> hop_length -> bool is provided. This utility allows users to programmatically check if their chosen window and hop length parameters satisfy the COLA constraint, mirroring the functionality of MATLAB's iscola and SciPy's check_COLA.

Performance

• The implementation must manage memory allocation efficiently within the overlap-add loop to prevent performance degradation when processing long audio signals.

Proposed Solution or Technical Approach (Optional)

1. Core Algorithm: Overlap-Add (OLA)

The fundamental method for iSTFT is the Overlap-Add (OLA) technique. The process involves these steps:

1. Each time frame (column) of the input STFT matrix is transformed back to a short time-domain signal segment using an Inverse Fast Fourier Transform (IFFT).
2. Each resulting time-domain segment is multiplied by a synthesis window.
3. These windowed segments are then overlapped and summed (added) together at intervals determined by the hop_length to reconstruct the full time-domain signal.
   The reconstruction formula can be expressed as:
   $x[n] = \frac{\sum_{t} y_t[n] \cdot w_s[n - tH]}{\sum_{t} w_a[n - tH] \cdot w_s[n - tH]}$
   where $y_t$ is the IFFT of the $t$-th frame, $w_a$ is the analysis window used in the forward STFT, $w_s$ is the synthesis window, and $H$ is the hop size. For perfect reconstruction, the denominator must be a non-zero constant.

2. COLA and NOLA Constraints

For the iSTFT to be a true inverse, the windowing and overlap parameters must satisfy certain constraints.

• NOLA (Nonzero Overlap-Add): This is the minimum requirement for invertibility. It ensures the denominator in the OLA equation is never zero, preventing division by zero but not necessarily preserving the signal's original amplitude.
• COLA (Constant Overlap-Add): This is a stricter condition where the denominator is a non-zero constant. Satisfying COLA ensures that all windowing artifacts are perfectly cancelled out, leading to perfect reconstruction of an unmodified STFT.
• Weak COLA is sufficient for perfect reconstruction of unmodified spectra but relies on alias cancellation that can be disrupted by spectral modifications.
• Strong COLA is preferred when the spectrum is modified (e.g., for audio effects), as it prevents aliasing altogether.

3. High-Level Algorithm Steps:

1. Initialization: Infer parameters (n_fft, win_length, etc.) from the input matrix and optional arguments. Generate the synthesis window array win. Pre-calculate the squared window win_sq = win * win. Initialize an output signal buffer y and a window normalization buffer win_sum_sq to zeros, with a length determined by the STFT dimensions and hop length.
2. Frame Processing Loop: Iterate through each time frame (column t) of the stft matrix.
3. IFFT: For each frame, perform an n_fft-point IFFT. If the input is one-sided (from a real signal), first reconstruct the full conjugate-symmetric spectrum.
4. Overlap-Add: Multiply the resulting time-domain segment by the synthesis window. Determine the start sample index n = t * hop_length. Add the windowed segment to the output buffer y and the squared window win_sq to the normalization buffer win_sum_sq at the appropriate offsets.
5. Normalization: After the loop, normalize the output signal: y = y / win_sum_sq. To prevent division by zero at the signal edges where the window sum might be zero, replace zeros in win_sum_sq with a small epsilon before division.
6. Finalization: If the length parameter is provided, trim or pad the final signal y to the specified length.

Additional Context (Mockups, Links, etc.)

Reference table of common window functions and their COLA-compliant overlap percentages should be included in the documentation.

Window Function	COLA-Compliant Overlap
Rectangular	0%, 50%, 75%
Hann	50%, 66.7% (2/3), 75%
Hamming	50%
Bartlett	50%

This information is derived from the documentation of SciPy and librosa.

References

• Smith, J.O. (2011). Spectral Audio Signal Processing. W3K Publishing. The definitive textbook on STFT, overlap-add methods, and windowing.
• Griffin, D., & Lim, J. (1984). "Signal estimation from modified short-time Fourier transform." IEEE Transactions on Acoustics, Speech, and Signal Processing, 32(2), 236-243. The seminal paper describing the least-squares approach to STFT inversion and the Griffin-Lim algorithm.
• Librosa: librosa.istft documentation.
• SciPy: scipy.signal.istft documentation. A key reference for the NOLA constraint and the least-squares implementation method.
• MATLAB: istft documentation. Strict, COLA-focused implementation and the utility of providing a COLA-checking function.

Code of Conduct

• I agree to follow this project's Code of Conduct.