Fourier Series & FFT Visualization | Epicycles, Butterflies

The idea: circles that add up to a wave

A periodic signal repeats itself every period. Fourier's claim is that such a signal is a sum of sine waves whose frequencies are whole multiples of the base frequency: one that turns once per period, one that turns twice, three times, and so on. This sum is the Fourier series.

A sine wave is the height of a point that goes round a circle. So a sum of sine waves is a chain of circles: each circle sits on the tip of the one before, circle n turns n times while the first turns once, and the height of the last tip is the sum. Circles riding on circles are called epicycles.

Three circles for n = 1, 3, 5: the height of the last tip, plotted against the angle θ, is the sum of three sine waves.

The coefficients are the radii

Write the series as f(θ) = b₁·sin θ + b₂·sin 2θ + b₃·sin 3θ + … The number b_n, the coefficient, is the radius of circle n: how much of the frequency n the signal contains. It is found by an average over one period, b_n = (1/π)·∫ f(θ)·sin(nθ) dθ.

signal	terms	radius of circle n	after 8 circles
square wave (±1)	odd n only: 1, 3, 5, …	4/(πn): 1.273, 0.424, 0.255, …	peak 1.18
sawtooth (π − θ)/π	every n: 1, 2, 3, …	2/(πn): 0.637, 0.318, 0.212, …	peak 1.07
heart (a 2D shape)	k = −4 … 4 except 0	12.5, 3, 2.5, 2.5, 1, 0.5, 0.5, 0.5	exact

The square wave is symmetric in a way that cancels every even term. The radii shrink like 1/n, so each new circle is smaller and the sum changes less. (Demo 1: square wave, Demo 3: sawtooth)

The Gibbs overshoot

Next to a jump the sum always overshoots. With 1 circle the square wave's sum peaks at 1.27, with 2 at 1.20, then 1.19, 1.18, … and it never gets below 1.179. The overshoot of 0.179 is about 9% of the jump from −1 to +1. More circles only make the bump narrower and push it closer to the jump. This is the Gibbs phenomenon. At every point away from a jump the sum still converges to the wave; it is the largest value that does not. (Demo 2: Gibbs)

The square wave near its jump at θ = 0, with the partial sums for 1, 3 and 8 terms. The 3- and 8-term sums wiggle and overshoot above 1; a dashed line at 1.179 marks the limit of the peak — More terms squeeze the bump towards the jump, but its height stays at about 1.18: 9% of the jump.

Drawing a shape: complex coefficients

A closed curve can be drawn the same way. Read each point (x, y) as one complex number z = x + iy. Then z(θ) = Σ c_k·e^ikθ: each term is a circle of radius |c_k| turning k times per turn, anticlockwise for k > 0 and clockwise for k < 0. The coefficient c_k is again an average, (1/2π)·∫ z(θ)·e^−ikθ dθ.

The heart curve x = 16 sin³θ, y = 13 cos θ − 5 cos 2θ − 2 cos 3θ − cos 4θ is made of only 8 terms, because sin³θ and cos kθ are themselves short sums of e^±ikθ. The page adds the largest first: c₋₁ = 12.5i, then c₋₃ = −3i, c₋₂ = c₂ = −2.5i, c₃ = i, and three of radius 0.5. After 8 circles the drawing is exact. A real drawing needs infinitely many terms; a few dozen look right. (Demo 4: heart)

From the series to the DFT

A computer has samples, not a formula. With N samples x[0] … x[N−1] the coefficient integral becomes a sum, the discrete Fourier transform (DFT):

X[k] = Σ x[n]·e^(−2πi·kn/N)      k = 0 … N−1,  sum over n = 0 … N−1

X[k] says how much of the frequency "k cycles per N samples" the signal contains. Done directly, each of the N outputs is a sum of N products: N² complex multiplications. That is 64 for N = 8, and about a million for N = 1024.

The FFT: split in half, then butterflies

The fast Fourier transform of Cooley and Tukey (1965) splits the samples into the even ones and the odd ones. Each half is a DFT of size N/2, call them E and O, and with the twiddle factor W = e^−2πi/N:

X[k]       = E[k] + W^k·O[k]
X[k + N/2] = E[k] − W^k·O[k]          k = 0 … N/2 − 1

One multiplication W^k·O[k] gives two outputs. This pair is drawn as a butterfly: a and b in, a + W·b and a − W·b out. Split the halves again until each piece is one sample. For N = 8 that gives 3 stages of 4 butterflies each, with twiddles W⁰ in stage 1, W⁰ and W² in stage 2, and W⁰ to W³ in stage 3.

The 8-point FFT flow graph: inputs in bit-reversed order x0, x4, x2, x6, x1, x5, x3, x7 on the left, three stages of crossing butterflies with a twiddle factor times W0 to W3 under the lower input of each butterfly, and outputs X0 to X7 in natural order on the right. Each butterfly takes a and b and gives a plus W times b and a minus W times b — The 8-point radix-2 FFT: bit-reversed input, 3 stages of 4 butterflies, outputs in order.

Splitting into even and odd again and again puts the inputs in bit-reversed order: write the index in binary and read it backwards.

row	binary	reversed	holds
0	000	000	x[0]
1	001	100	x[4]
2	010	010	x[2]
3	011	110	x[6]
4	100	001	x[1]
5	101	101	x[5]
6	110	011	x[3]
7	111	111	x[7]

Counting multiplications

Each of the log₂N stages has N/2 butterflies with one multiplication each: (N/2)·log₂N in total, against N² for the direct DFT.

N	direct DFT: N²	FFT: (N/2)·log₂N	FFT is faster by
8	64	12	5×
64	4,096	192	21×
1,024	1,048,576	5,120	205×
1,048,576	about 1.1 × 10¹²	10,485,760	about 100,000×

The page counts the multiplications by W⁰ = 1 too, as the formula does; real libraries skip them. The split is the same divide and conquer as merge sort: two halves plus linear work to combine them, so O(N log N) (see Quicksort vs Merge Sort).

Reading the spectrum

For real samples, X[N−k] is the complex conjugate of X[k], so the bars on the right half mirror the left half. A cosine is half e^+iωn and half e^−iωn, so one cycle per 8 samples puts 4 in bin 1 and 4 in bin 7, the frequency −1. (Demo 5) A single spike, the impulse, contains every frequency equally: all eight |X[k]| are 1. (Demo 6) A cosine of 1.5 cycles does not fit a whole number of times into 8 samples, so no bin matches it and the energy leaks into every bin: 1.00, 2.40, 2.80, 1.19, 1.00, 1.19, 2.80, 2.40. (Demo 7)

The FFT is also how convolution is done fast: transform both, multiply bin by bin, transform back (see Convolution in a CNN). Multiplying by W^k is a rotation by −45°·k, the complex form of a rotation matrix.

Taylor Series builds a function from powers of x that fit at one point, where a Fourier series uses sines that fit the whole period. The twiddle factors are roots of unity: see Complex Numbers for why multiplying by e^(−2πi/N) turns a point round the unit circle.

Each spinning circle is a point going round a circle; The Unit Circle: sin, cos, tan shows how one turn draws a sine and a cosine wave, and what a phase shift is.

Adding two waves of the same frequency that come from two different places, instead of many frequencies at one place, gives bright and dark fringes: see Wave Interference: Two Coherent Sources.

The Fourier transform weighs a signal by a spinning e^−iωt; the Laplace transform uses e^−st with s = σ + iω, which adds a decay so that even growing signals have a transform. See The Laplace Transform.

For where each coefficient comes from (the area of f·sin nx), cosine terms and a₀, and the Gibbs overshoot measured up to 99 harmonics, see Fourier Series, Harmonic by Harmonic.

JPEG uses a cousin of the DFT, the cosine-only DCT, on 8 × 8 blocks: see JPEG Compression: DCT and Quantization.

Feed the DFT 256 points around a heart, a star or the letter π and every coefficient becomes a spinning circle that helps draw it: see Fourier Epicycles: Drawing a Shape with Spinning Circles.

What the page leaves out

The series uses sine terms only; a general signal also has cosine terms and a constant a₀, and the page does not show when a series converges. The circles move between 16 positions per turn in straight lines, so the small fast circles jump. The FFT is the plain radix-2 version for N = 8 with separate columns; real libraries (FFTW, cuFFT) work in place, use real-input, radix-4 and split-radix variants, handle any N (Bluestein), and skip trivial twiddles. The inverse FFT, windowing against leakage and rounding errors are not shown.