Project Grimoire · Inner Sanctum

The Fourier Transform

Discrete Spectral Analysis · C & Python

I — The Working

The Discrete Fourier Transform decomposes a finite sequence of time-domain samples into their constituent sinusoidal frequencies — transforming a signal from the time domain into the frequency domain.

This implementation takes a composite signal built from two frequencies — 4 Hz and 12 Hz — and reveals them precisely within the magnitude spectrum, demonstrating the DFT's power to unmask hidden periodicity.

Sample Count

N = 64

Sample Rate

Fs = 64 Hz

Components

f₁=4, f₂=12

II — The Signal

Time-Domain Waveform

x[n] = sin(2π·4·n/N) + 0.5·cos(2π·12·n/N)

Live Render

Composite x[n]

sin(2π·4·n/N)

0.5·cos(2π·12·n/N)

III — The Formula

X[k] = ∑_n=0^N−1 x[n] · e^−j2πknN

Discrete Fourier Transform · k = 0, 1, …, N−1

IV — The Spectrum

Magnitude Spectrum |X[k]|

One-sided, bins k = 0 … N/2

DFT Output

V — The Algorithm

01
A composite test signal is synthesised: sin(2π·4·n/N) + 0.5·cos(2π·12·n/N) — two sinusoids at frequencies 4 and 12 cycles per window, with amplitudes 1 and 0.5 respectively.
02
The DFT kernel iterates over all N output bins k. For each bin it accumulates the inner product of the input with a complex exponential e^(−j2πkn/N) — the cosine and sine components computed separately and combined.
03
The magnitude at each bin is computed as |X[k]| = √(Re² + Im²). The algorithm is O(N²) — each of the N output bins requires N multiplications. The FFT reduces this to O(N log N) via the Cooley–Tukey butterfly factorisation.
04
Because the input is real-valued, the spectrum is conjugate-symmetric: the peak at k = 4 mirrors at k = N−4 = 60, and k = 12 mirrors at k = 52. Only bins 0…N/2 carry unique information.
05
The full spectrum is written to dft_output.txt and then read by the Python companion script, which renders three plots: the time-domain signal, the two-sided spectrum, and the one-sided spectrum with annotated peaks.

VI — The Inscription

dft.c

Core transform — C99, compiled with gcc -lm

/* DFT — O(N²) naive implementation */
void dft(const double *x_real, const double *x_imag,
         double *X_real, double *X_imag, int n_samples)
{
    for (int k = 0; k < n_samples; k++) {
        double sum_real = 0.0;
        double sum_imag = 0.0;
        for (int n = 0; n < n_samples; n++) {
            double angle = -2.0 * PI * k * n / n_samples;
            sum_real += x_real[n] * cos(angle)
                      - x_imag[n] * sin(angle);
            sum_imag += x_real[n] * sin(angle)
                      + x_imag[n] * cos(angle);
        }
        X_real[k] = sum_real;
        X_imag[k] = sum_imag;
    }
}

/* Test signal: 4 Hz + 0.5 × 12 Hz, N=64 samples */
for (int n = 0; n < N; n++) {
    double t = (double)n / Fs;
    x_real[n] = 1.0 * sin(2.0 * PI * f1 * t)
              + 0.5 * cos(2.0 * PI * f2 * t);
}

VII — The Findings

Peak at k = 4

|X[4]| ≈ 32

Peak at k = 12

|X[12]| ≈ 16

Complexity

O(N²)