Computational Audio and Image Analysis with the Spectrograms Library

Different implementations use different scaling conventions for the DFT and inverse DFT. The three most common conventions are:

1. No forward scaling (this library): Forward DFT uses coefficient 1, inverse DFT uses $1/N$
2. Symmetric scaling: Both forward and inverse use $1/\sqrt{N}$ (unitary transform)
3. No inverse scaling: Forward uses $1/N$, inverse uses 1

This library uses Convention 1: The forward DFT has no normalization factor (coefficient of 1), and the inverse DFT includes the $1/N$ factor. This convention is standard in most signal processing literature and libraries (FFTW, NumPy, SciPy).

Energy implications: With this scaling, Parseval's theorem takes the form:

$$\sum_{n=0}^{N-1} |x[n]|^2 = \frac{1}{N}\sum_{k=0}^{N-1} |X[k]|^2$$

When comparing spectral amplitudes between implementations, always verify the scaling convention to ensure energy-consistent analysis.

Let $W_N^{k} = e^{-j2\pi k/N}$ be the twiddle factor. These complex exponentials satisfy:

$$\begin{aligned} W_N^{k+N/2} &= -W_N^{k} \quad \text{(symmetry)} \ W_N^{k+N} &= W_N^{k} \quad \text{(periodicity)} \end{aligned}$$

Periodicity means the complex exponential repeats every $N$ samples. Symmetry means that values half an FFT bin apart differ only in sign. These properties enable the FFT's recursive decomposition.

For a real-valued signal $x[n]$, the DFT satisfies Hermitian symmetry:

$$X[k] = X^*[N-k]$$

where $X^*$ denotes complex conjugate. This means:

• $\operatorname{Re}(X[k]) = \operatorname{Re}(X[N-k])$ (real parts are symmetric)
• $\operatorname{Im}(X[k]) = -\operatorname{Im}(X[N-k])$ (imaginary parts are antisymmetric)
• $X[0]$ and $X[N/2]$ (for even $N$) are purely real

This property follows directly from the DFT definition and the fact that $e^{j\theta} + e^{-j\theta} = 2\cos\theta$ is real for real $\theta$.

Critical constraint for R2C transforms: For even $N$, the Nyquist frequency bin $X[N/2]$ must be purely real (i.e., $\operatorname{Im}(X[N/2]) = 0$). This follows from Hermitian symmetry since $X[N/2] = X^$[N - N/2] = X^[N/2], which requires the imaginary part to be zero. Violating this constraint results in non-real values after inverse FFT, causing numerical errors and incorrect signal reconstruction.

When the DFT of a finite signal segment is computed, the signal is implicitly multiplied by a rectangular window (ones inside the segment, zeros outside). This abrupt truncation causes spectral leakage: energy from one frequency "leaks" into adjacent frequency bins.

Consider a pure sinusoid at frequency $f_0$. If sampling captures exactly an integer number of periods, the DFT produces a perfect spike at $f_0$. But with non-integer periods, the energy spreads across many bins.

This occurs because multiplying by a rectangular window in the time domain is equivalent to convolving with a sinc function in the frequency domain, which has slowly-decaying sidelobes.

Window functions can be defined with two distinct conventions:

Symmetric windows use $N-1$ in the denominator:

$$w_{\text{sym}}[n] = f\left(\frac{n}{N-1}\right), \quad n \in [0, N-1]$$

Periodic windows use $N$ in the denominator:

$$w_{\text{per}}[n] = f\left(\frac{n}{N}\right), \quad n \in [0, N-1]$$

The difference is subtle but important:

• Symmetric: Values at endpoints match (e.g., both zero for Hanning). Standard for spectral analysis.
• Periodic: The window is periodic with period $N$, so $w[N] = w[0]$. Required for perfect reconstruction in overlap-add systems.

For the Constant Overlap-Add (COLA) property:

$$\sum_{m=-\infty}^{\infty} w[n - mH] = C \quad \text{(constant)}$$

where $H$ is the hop size. Periodic windows satisfy COLA for specific overlap ratios, enabling perfect signal reconstruction.

Library convention: All window implementations in this library use the symmetric formulation ($N-1$), which is standard for analysis applications. For synthesis requiring COLA, use appropriate overlap ratios (e.g., 50% for Hanning).

The modified Bessel function of the first kind, order zero, $I_0(x)$ is defined by:

$$I_0(x) = \sum_{k=0}^{\infty} \frac{1}{(k!)^2}{\left(\frac{x}{2}\right)^{2k}}$$

For computational efficiency, a truncated power series is used. The series converges rapidly; typically 20–30 terms achieve machine-precision accuracy for $x < 10$.

When using the STFT for signal modification or synthesis (e.g., time-stretching, pitch-shifting, noise reduction), the time-domain signal must be reconstructed from modified spectral frames. The Constant Overlap-Add (COLA) constraint ensures perfect reconstruction is possible.

For a window function $w[n]$ and hop size $H$, the COLA property states that the sum of overlapping windows must be constant:

$$\sum_{m=-\infty}^{\infty} w[n - mH] = C \quad \text{(constant for all } n\text{)}$$

where $C$ is a non-zero constant (typically normalized to 1).

Physical interpretation: When frames are overlapped by hop size $H$ and summed (overlap-add), every sample in the reconstructed signal receives equal weighting from the windows. This prevents amplitude modulation artifacts.

Why pad $N/2$ zeros? Without centering, the first window is positioned with its left edge at sample 0. This means the window's center is at sample $N/2$, causing a temporal shift in the resulting spectrogram.

By padding $N/2$ zeros before the signal, the first window's center aligns with sample 0. This ensures:

• Frame $m$ is centered at time $t = mH/f_s$ (no offset)
• Edge effects are symmetric at signal boundaries
• Improved temporal localization for transient events

The same $N/2$ padding is applied at the signal's end, ensuring symmetric treatment of boundaries. This is standard practice in audio analysis libraries (e.g., librosa's center=True).

Padding type specification: This library uses zero-padding (also called zero-extension) for both the initial and final $N/2$ samples:

$$\tilde{x}[n] = \begin{cases} 0 & n < 0 \ x[n] & 0 \leq n < L \ 0 & n \geq L \end{cases}$$