Reconstruction and Multirate Processing

Sinc interpolation, the zero-order hold and its droop, linear interpolation, decimation, interpolation and polyphase resampling

1982 — the first CD players

The first compact disc player on sale, Sony's CDP-101, appears in Japan in October 1982 (Compact Disc Digital Audio). Philips has a problem. During the negotiations it argued for 14-bit samples and had already built a 14-bit digital-to-analog converter; Sony insisted on 16 bits, and won. Philips's answer is arithmetic, not a new chip: it runs its 14-bit converter at four times the CD's sampling rate, and gets 16-bit quality (Compact Disc Digital Audio). Oversampling — computing new samples between the old ones, then converting faster — also eases a problem every converter has. A digital-to-analog converter does not produce the smooth interpolation of the sampling theorem. It holds each value until the next one arrives: a staircase. This lesson is about getting back from numbers to signals, and about changing a signal's sampling rate without damaging it.

Lesson 5-1 recovered a bandlimited signal from its samples with an ideal low-pass filter. Applied to the impulse train

x_p(t)

, that filter turns each sample into a pulse; the pulse's shape is the filter's impulse response, and the sum of the pulses is the reconstruction. Ideal and practical reconstruction differ only in the pulse.

Derivation

Sinc interpolation

The ideal reconstruction filter has gain

T

for

|\omega| < \pi/T

and zero elsewhere. Its impulse response is

h_r(t) = \frac{1}{2\pi}\int_{-\pi/T}^{\pi/T} T\,e^{j\omega t}\,d\omega = \frac{\sin(\pi t/T)}{\pi t/T} = \operatorname{sinc}\!\left(\frac{t}{T}\right).

Filtering

x_p(t) = \sum_n x(nT)\,\delta(t - nT)

places one sinc at every sample:

x(t) = \sum_{n=-\infty}^{\infty} x(nT)\,\operatorname{sinc}\!\left(\frac{t - nT}{T}\right).

At a sampling instant

t = mT

\operatorname{sinc}(m - n)

is 1 for

n = m

and 0 for every other

n

, so the curve passes through every sample; between the instants, every sample contributes.

x(t) = \sum_{n=-\infty}^{\infty} x[n]\,\operatorname{sinc}\!\left(\frac{t - nT}{T}\right)

The Whittaker–Shannon interpolation formula: a bandlimited signal rebuilt exactly from its samples, one sinc per sample.

💡The zero-order hold

The sinc reaches infinitely far in both directions, so exact reconstruction needs every sample, future ones included — the non-causal ideal low-pass of Lesson 4-4 again. A converter instead holds each sample for one period: its pulse is a rectangle of width

T

. Normalised to the ideal gain

T

, its frequency response is

\frac{H_0(f)}{T} = e^{-j\pi fT}\,\operatorname{sinc}(fT):

a delay of half a sample, and a droop that reaches

\operatorname{sinc}(\tfrac12) = 2/\pi

, or −3.92 dB, at the Nyquist frequency. Worse, the sinc only dents the spectral copies. On a CD, a 20 kHz tone comes out 3.17 dB low, and its image at 24.1 kHz only 1.6 dB below that, so an analog smoothing filter after the converter must remove the images.

Straight lines and two cures. Joining the samples with straight lines — linear interpolation, a triangular pulse of width

2T

— has the response

\operatorname{sinc}^2(fT)

: twice the droop in decibels, −7.84 dB at the Nyquist frequency, but images that fall away twice as fast. Droop itself has two cures. An inverse-sinc filter can boost the samples before conversion — a 31-tap FIR filter flattens a CD's audio band to within 0.02 dB — or oversampling can move the band onto the flat top of the sinc: at four times 44.1 kHz, the droop at 20 kHz is only 0.18 dB.

💡Changing the sampling rate

Decimation by $M$ . Low-pass the sequence to $\Omega = \pi/M$ , then keep every $M$ th sample. Skipping the filter aliases (Lesson 5-2).
Interpolation by $L$ . Insert $L - 1$ zeros between samples. This squeezes the spectrum by $L$ and leaves $L - 1$ unwanted copies in the new band; a low-pass filter with gain $L$ and cut-off $\pi/L$ removes them.
A rational factor $L/M$ . Interpolate by $L$ , filter once with the narrower of the two cut-offs, decimate by $M$ . A polyphase implementation never multiplies the inserted zeros and never computes the samples it is about to throw away, so all the work happens at the lower rate.

📐From 44.1 kHz to 48 kHz

Studio and video audio runs at 48 kHz, CDs at 44.1 kHz. The ratio is

48\,000/44\,100 = 160/147

, so conversion interpolates by 160 and decimates by 147. Done literally, the intermediate rate would be

44\,100 \times 160 = 7.056

MHz. A polyphase resampler such as SciPy's resample_poly computes only the 48 000 output samples a second, each from a short slice of the filter.

The art of being wise is the art of knowing what to overlook.

– William James, The Principles of Psychology (1890), chapter 22, quoted on Wikiquote

🕰Older than the theorem

The interpolation formula predates the sampling theorem. It goes back to Émile Borel in 1898 and to E. T. Whittaker, who published it in 1915 and called it the cardinal series; Shannon used it in his 1949 formulation of the theorem (Whittaker–Shannon interpolation formula).

⚠Upsampling makes images, downsampling makes aliases

Inserting zeros adds no information; it exposes spectral copies that the interpolation filter must remove. Dropping samples does the opposite: the copies move closer and overlap unless a filter clears the way first. Rate changes in both directions need the low-pass filter — naive slicing and zero-stuffing are never enough on their own.

💡Where you will use this

The DTFT of Lesson 6-1 is the periodic spectrum these filters trim, and zero padding (Lesson 6-3) interpolates spectra the way upsampling interpolates signals. The low-pass filters inside every resampler are the FIR designs of Lesson 9-1, and Lesson 9-3 implements them efficiently.

📚Further reading

Oppenheim and Willsky, Signals and Systems, sections 7.2 and 7.5 — reconstruction by interpolation, and decimation and interpolation of sequences.
Ronald E. Crochiere and Lawrence R. Rabiner, Multirate Digital Signal Processing (Prentice-Hall, 1983) — the classic account of polyphase resampling.
Zero-order hold on Wikipedia.

↓But there is a problem

Inside the computer a signal is now a list of numbers, and its spectrum must come from those numbers alone. What does a frequency even mean for a sequence? In 2020, dashboards around the world quietly replaced daily case counts with seven-day averages. The reason is a frequency response — and it has zeros in exactly the right places.