Harmonic analysis

Averaging the partial sums cancels exactly the oscillation that makes $S_N f$ misbehave: while $D_N$ has wild negative lobes, its running average $F_N$ is a smooth, nonnegative bump that only ever adds up pieces of $f$ with a positive weight. So even though $S_N f$ can fail to converge, its running average always settles down — the same trick that tames a jittery sequence of partial sums in calculus by looking at its Cesàro average instead.

Step 1 (Fejér kernel is a positive approximate identity). Averaging $D_0,\ldots,D_N$ and using a telescoping trigonometric identity gives the closed form $F_N(x) = \frac{1}{N+1}\left(\dfrac{\sin\left(\tfrac{N+1}{2}x\right)}{\sin(x/2)}\right)^2 \ge 0$. Since $F_N = \frac{1}{N+1}\sum_{k=0}^N D_k$ and $\frac{1}{2\pi}\int D_k = 1$ for every $k$, also $\frac{1}{2\pi}\int F_N = 1$; being nonnegative, $\frac{1}{2\pi}\int |F_N| = 1$ too, so the $L^1$ norm never blows up.

Step 2 (mass concentrates at $0$). For any fixed $\delta>0$, on $\delta \le |x| \le \pi$ the denominator $\sin(x/2)^2$ is bounded below by a positive constant, so $F_N(x) = O(1/N)$ uniformly there; the mass outside $[-\delta,\delta]$ tends to $0$ as $N\to\infty$.

Step 3 (approximate-identity estimate). Write $\sigma_N f(x) - f(x) = \frac{1}{2\pi}\int_{-\pi}^{\pi} F_N(y)\big(f(x-y)-f(x)\big)\,dy$ using $\frac{1}{2\pi}\int F_N = 1$. Split the integral into $|y|<\delta$ and $\delta\le|y|\le\pi$. On the first piece, uniform continuity of $f$ makes $|f(x-y)-f(x)|<\varepsilon$ for $\delta$ small, and $F_N\ge0$ with total mass $1$ bounds that piece by $\varepsilon$. On the second piece, $F_N=O(1/N)$ and $f$ is bounded, so that piece $\to 0$ as $N\to\infty$. Both bounds are uniform in $x$, so $\sigma_N f\to f$ uniformly.

Given Kolmogorov's 1926 example of an $L^1$ function that diverges everywhere, it seemed entirely plausible that $L^2$ — barely a stronger condition — would fail the same way, or at least almost everywhere. Carleson's theorem is a genuine surprise: the tiny extra assumption of square-integrability rules out divergence everywhere except a set of measure zero. It closed a problem that had been open since Luzin conjectured it in 1913, and is regarded as one of the deepest theorems of twentieth-century analysis.

A full proof is one of the hardest arguments in twentieth-century analysis; Fefferman gave a celebrated simplification in 1973 whose strategy can be sketched in three moves.

Step 1 (control by a maximal operator). It suffices to bound the Carleson maximal operator $S^*f(x) = \sup_N |S_N f(x)|$ on $L^2$, since a weak-type bound on $S^*$ plus density of nice functions (for which convergence is easy) implies the full a.e. convergence statement.

Step 2 (time–frequency decomposition into tiles). Decompose $f$ using wave packets adapted to dyadic tiles in the time–frequency plane — rectangles of area $\sim 1$, exactly the atoms discussed in the 'Time and frequency' section below. Each tile carries a piece of $f$ localized to both a time interval and a frequency band.

Step 3 (organize tiles into trees and sum). Because the cutoff frequency $N$ can vary with $x$, tiles relevant to $S^*f(x)$ form 'trees' ordered by time-frequency containment. Fefferman's key combinatorial estimate bounds the total energy carried by all trees, controlling $S^*f$ in $L^2$ and completing the proof.

Physically, $\int |f|^2$ is the energy of a signal (think: power dissipated by a voltage waveform). Plancherel says you can compute that energy either by scanning the signal in time or by scanning its spectrum in frequency — a spectrum analyzer and an oscilloscope must agree on total power.

Step 1 (a special case). First check the identity for a Gaussian $f(x) = e^{-\pi x^2}$. A direct computation (completing the square in the exponent) gives $\hat f(\xi) = e^{-\pi \xi^2}$, so both sides of $\int_{-\infty}^{\infty} |f(x)|^2\,dx = \int_{-\infty}^{\infty} |\hat f(\xi)|^2\,d\xi$ equal the same Gaussian integral $\int e^{-2\pi x^2}\,dx$.

Step 2 (extend by linearity and the convolution theorem). For nice (Schwartz) functions $f$, write $g = f * \tilde f$ where $\tilde f(x) = \overline{f(-x)}$; then $g(0) = \int |f(x)|^2\,dx$, and by the convolution theorem $\widehat{f * g} = \hat f \cdot \hat g$ together with Fourier inversion, $g(0) = \int \hat g(\xi)\,d\xi = \int |\hat f(\xi)|^2\,d\xi$. This proves the identity for all Schwartz functions, which include Gaussians and all smooth rapidly-decaying functions.

Step 3 (density argument). Schwartz functions are dense in $L^2(\mathbb{R})$: every $f \in L^2$ is a limit $f_n \to f$ of Schwartz functions in the $L^2$ norm. Since the Fourier transform is an isometry on this dense subspace by Step 2, it extends uniquely to a bounded operator on all of $L^2(\mathbb{R})$ that still satisfies $\int_{-\infty}^{\infty} |f(x)|^2\,dx = \int_{-\infty}^{\infty} |\hat f(\xi)|^2\,d\xi$ — this extension is what 'the Fourier transform of an $L^2$ function' means when the defining integral $\hat f(\xi) = \int_{-\infty}^{\infty} f(x)\, e^{-2\pi i x \xi}\,dx$ may not converge absolutely.

It is the bridge between sampling a signal at integer points and periodizing its spectrum: the left side is what you get by adding up samples $f(n)$, the right side is what you get from the spectrum. This single identity underlies the sampling theorem, lattice sums in crystallography, and the functional equation of theta functions and the Riemann zeta function.

Step 1 (periodize). For a Schwartz function $f$, define $F(x) = \sum_{n=-\infty}^{\infty} f(x+n)$. Rapid decay of $f$ makes this sum converge absolutely and uniformly, and $F$ is a smooth function with period 1, so it has its own Fourier series on the circle.

Step 2 (compute the Fourier coefficients of the periodization). The $k$-th Fourier coefficient of $F$ is $c_k = \int_0^1 F(x) e^{-2\pi i k x}\,dx = \sum_{n} \int_0^1 f(x+n) e^{-2\pi i k x}\,dx$. Substituting $y = x+n$ in each term and using periodicity of $e^{-2\pi i k x}$ in $n$ reassembles the pieces into a single integral over all of $\mathbb{R}$: $c_k = \int_{-\infty}^{\infty} f(y) e^{-2\pi i k y}\,dy = \hat f(k)$.

Step 3 (evaluate at x=0). Since $F$ is smooth, its Fourier series converges to it pointwise, in particular at $x=0$: $F(0) = \sum_k c_k e^{0}$, i.e. $\sum_{n=-\infty}^{\infty} f(n) = \sum_{k=-\infty}^{\infty} \hat f(k)$. The left side is literally $F(0) = \sum_n f(n)$ by definition, which finishes the proof.

The $L^2$ case is essentially free: the multiplier $-i\,\mathrm{sgn}(\xi)$ has modulus exactly $1$ almost everywhere, so $H$ is literally an isometry on $L^2$ by Plancherel — a $90^\circ$ phase rotation of every frequency changes nothing about total energy. But extending boundedness to other $p$ is a genuinely different and much harder problem: there is no analogue of Parseval's identity outside $L^2$, so no algebraic shortcut is available, and the proof has to fall back on real-variable estimates about the size and geometry of the set where $Hf$ is large.

Step 1 ($L^2$, directly). By Plancherel and the multiplier formula, $\|Hf\|_{L^2}^2 = \int |\widehat{Hf}(\xi)|^2\,d\xi = \int |\mathrm{sgn}(\xi)|^2|\hat f(\xi)|^2\,d\xi = \|f\|_{L^2}^2$, using $|\mathrm{sgn}(\xi)|=1$ for $\xi\ne0$ (a single point does not affect the integral). So $H$ is an isometry on $L^2(\mathbb{R})$, in particular bounded with $C_2=1$.

Step 2 (general $1<p<\infty$, forward reference). The full range of exponents does not follow from Step 1 by any soft argument. It is a genuine theorem of real-variable Calderón–Zygmund theory: the Hilbert transform's kernel $1/(\pi x)$ satisfies exactly the smoothness and cancellation conditions of a Calderón–Zygmund kernel, so the general singular-integral machinery developed in the next section applies to it. That machinery produces a weak-type $(1,1)$ bound from the $L^2$ bound proved in Step 1, and Marcinkiewicz interpolation between the weak-$(1,1)$ bound and the $L^2$ bound (then a duality argument for $p>2$) completes the proof for every $1<p<\infty$ — this half of the argument is carried out in full in the Calderón–Zygmund section immediately below, rather than repeated here.

$Mf(x)$ bounds, in one stroke, every average of $f$ that could ever be taken over an interval around $x$ — the running mean at every possible scale. Controlling that single quantity turns out to be exactly what is needed to prove the Lebesgue differentiation theorem (the average of $f$ over shrinking intervals around $x$ converges to $f(x)$ for almost every $x$): once $Mf$ is known to be finite almost everywhere, a short soft argument upgrades that to the full differentiation statement. This is the real-variable engine behind the Calderón–Zygmund theory below, playing the role that Plancherel's identity played for the easy $L^2$ estimate.

Step 1 (a cover by good balls). Fix $\lambda>0$. For every $x$ with $Mf(x)>\lambda$, by definition of the supremum there is some interval $B_x$ centered at $x$ with $\frac{1}{|B_x|}\int_{B_x}|f| > \lambda$. These balls $\{B_x\}$ cover the set $\{Mf>\lambda\}$.

Step 2 (Vitali $5r$-covering lemma). From any collection of balls of bounded radius, one can always extract a countable disjoint subcollection $\{B_i\}$ such that the $5$-times dilated balls $\{5B_i\}$ still cover the union of the whole original collection. Apply this to $\{B_x\}$ to get a disjoint subfamily $\{B_i\}$ with $\{Mf>\lambda\}\subseteq\bigcup_i 5B_i$.

Step 3 (sum the disjoint pieces). Each $B_i$ satisfies $\lambda|B_i| < \int_{B_i}|f|$ by construction. Since the $B_i$ are pairwise disjoint, summing over $i$ gives $\lambda\sum_i|B_i| < \sum_i\int_{B_i}|f| \le \int_{\mathbb{R}}|f| = \|f\|_1$, so $\sum_i|B_i| < \|f\|_1/\lambda$. Finally $|\{Mf>\lambda\}| \le \sum_i|5B_i| = 5\sum_i|B_i| < \dfrac{5}{\lambda}\|f\|_1$, which is exactly $|\{x : Mf(x) > \lambda\}| \le \frac{5}{\lambda}\|f\|_1$.