Digital signal processing (DSP) is a subfield of electrical engineering and applied mathematics that focuses on the manipulation, analysis, and transformation of discrete-time signals—typically represented as sequences of numbers—using digital computation. Unlike analog signal processing, DSP operates on sampled and quantized data, enabling high precision, flexibility, and reproducibility. Its applications range from audio compression and image enhancement to radar systems and biomedical instrumentation.
1 Fundamentals of digital signals and systems
1.1 Discrete-time signals and sequences
A discrete-time signal is a function defined only at integer time indices, often obtained by sampling a continuous-time signal at uniform intervals. Such signals are represented as sequences \(x[n]\), where \(n \in \mathbb{Z}\). Common examples include unit impulse \(\delta[n]\), unit step \(u[n]\), and sinusoidal sequences. Sequences may be finite or infinite in length and are classified as energy signals or power signals based on their energy or average power.
1.2 Sampling and quantization
Converting a continuous-time signal to a discrete-time representation involves two steps: sampling (taking values at discrete instants) and quantization (approximating each sample by a finite-precision number). The fidelity of the digital representation depends critically on both the sampling rate and the number of quantization levels.
1.2.1 The sampling theorem (Nyquist–Shannon)
The Nyquist–Shannon sampling theorem states that a bandlimited continuous-time signal with maximum frequency \(f_{\text{max}}\) can be perfectly reconstructed from its samples if the sampling frequency \(f_s\) satisfies \(f_s > 2f_{\text{max}}\). The critical frequency \(f_s/2\) is called the Nyquist frequency. Violating this condition leads to aliasing, an irreversible distortion.
1.2.2 Quantization error and noise
Quantization maps each sample to the nearest level in a finite set. The difference between the original sample value and its quantized version is quantization error, often modeled as additive white noise uniformly distributed between \(\pm \Delta/2\), where \(\Delta\) is the step size. The signal-to-quantization-noise ratio (SQNR) improves by approximately 6 dB per additional bit.
1.2.3 Aliasing and anti-aliasing filters
When a signal contains frequencies above the Nyquist frequency, those components are folded into lower frequencies, creating aliased artifacts. To prevent this, an analog low-pass anti-aliasing filter is applied before sampling. In digital-to-analog conversion, a reconstruction (smoothing) filter removes spectral replicas.
1.3 Linear time-invariant (LTI) systems
LTI systems are the cornerstone of classical DSP. A system is linear if it obeys superposition and scaling; it is time-invariant if its behavior does not change over time. LTI systems are fully characterized by their impulse response.
1.3.1 Impulse response and convolution
The impulse response \(h[n]\) is the output of an LTI system when the input is a unit impulse \(\delta[n]\). For any input \(x[n]\), the output \(y[n]\) is given by the convolution sum: \(y[n] = x[n] * h[n] = \sum_{k=-\infty}^{\infty} x[k] h[n-k]\). Convolution is commutative, associative, and distributive.
1.3.2 Difference equations and system function
Many LTI systems are described by linear constant-coefficient difference equations: \(\sum_{k=0}^{N} a_k y[n-k] = \sum_{k=0}^{M} b_k x[n-k]\). The system function (transfer function) in the z-domain is \(H(z) = \frac{Y(z)}{X(z)} = \frac{\sum_{k=0}^{M} b_k z^{-k}}{\sum_{k=0}^{N} a_k z^{-k}}\), which directly reveals stability and frequency response.
2 Transform-domain analysis
2.1 The discrete-time Fourier transform (DTFT)
The DTFT maps a discrete-time sequence \(x[n]\) into a continuous, periodic function of frequency: \(X(e^{j\omega}) = \sum_{n=-\infty}^{\infty} x[n] e^{-j\omega n}\). It is defined for absolutely summable sequences and provides the spectral representation of a signal.
2.1.1 Properties and symmetry
Key properties include linearity, time shifting (multiplication by \(e^{-j\omega n_0}\)), convolution (product in frequency), and Parseval’s theorem (energy preservation). For real-valued sequences, the DTFT exhibits conjugate symmetry: \(X(e^{j\omega}) = X^*(e^{-j\omega})\). This implies that the magnitude spectrum is even and the phase spectrum is odd.
2.1.2 Frequency response of LTI systems
| The frequency response \(H(e^{j\omega})\) of an LTI system is the DTFT of its impulse response. For a stable system, \(H(e^{j\omega}) = H(z)\) evaluated on the unit circle \(z = e^{j\omega}\). The magnitude \( | H(e^{j\omega}) | \) and phase \(\angle H(e^{j\omega})\) describe how the system attenuates or amplifies different frequencies. |
|---|
2.2 The discrete Fourier transform (DFT)
The DFT is a finite, discrete version of the DTFT, defined for a sequence of length \(N\): \(X[k] = \sum_{n=0}^{N-1} x[n] e^{-j2\pi kn/N}\), for \(k=0,\dots,N-1\). The DFT is widely used in practice because it is computable and invertible. The inverse DFT (IDFT) reconstructs the original sequence.
2.2.1 Fast Fourier transform (FFT) algorithms
The FFT is an efficient algorithm to compute the DFT, reducing complexity from \(O(N^2)\) to \(O(N \log N)\). The most common FFT algorithms are radix-2, which require \(N\) to be a power of two.
###### 2.2.1.1 Radix-2 decimation-in-time
Decimation-in-time (DIT) recursively splits the input sequence into even- and odd-indexed samples, computes DFTs of half the length, and recombines them using twiddle factors. The Cooley–Tukey algorithm is the classic implementation.
###### 2.2.1.2 Radix-2 decimation-in-frequency
Decimation-in-frequency (DIF) instead splits the output frequency points into even and odd groups. It uses a similar butterfly structure but with twiddle factors applied before the recursion. Both DIT and DIF have identical computational cost.
2.2.2 Spectral leakage and windowing
Because the DFT assumes periodic extension of the finite-length sequence, discontinuities at the boundaries cause spectral leakage—energy spreading from a true frequency into adjacent bins. Windowing (e.g., Hamming, Hann, Blackman) multiplies the signal by a smoothly tapering function to reduce leakage at the cost of broader mainlobes.
2.3 The z-transform
The z-transform generalizes the DTFT and is defined as \(X(z) = \sum_{n=-\infty}^{\infty} x[n] z^{-n}\), where \(z\) is a complex variable. It is especially useful for analyzing LTI systems and solving difference equations.
2.3.1 Region of convergence (ROC)
The ROC is the set of \(z\) for which the series converges. For finite-length sequences the ROC is the entire \(z\)-plane except possibly \(z=0\) or \(z=\infty\). For right-sided sequences the ROC is outside a circle; for left-sided sequences it is inside. The ROC cannot contain any poles and must be a connected region.
2.3.2 Poles, zeros, and stability
The zeros of a system function are the roots of the numerator polynomial; the poles are the roots of the denominator. An LTI system is stable (BIBO) if and only if its ROC includes the unit circle. For causal systems, stability requires all poles to lie inside the unit circle.
3 Digital filter design
3.1 Finite impulse response (FIR) filters
FIR filters have an impulse response that is finite in length (i.e., \(h[n]=0\) for \(n<0\) and \(n>M\)). They are always stable and can be designed with exact linear phase, making them ideal for applications requiring phase distortionlessness.
3.1.1 Linear-phase design
A linear-phase FIR filter has a phase response that is a linear function of frequency, corresponding to a constant group delay. This is achieved by making the impulse response symmetric or antisymmetric about its center. The design problem reduces to specifying the magnitude response under the linear-phase constraint.
3.1.2 Window-based design method
The window method begins with an ideal (infinite-length) impulse response \(h_{\text{id}}[n]\) derived from the desired frequency response. It is then truncated and multiplied by a window function \(w[n]\) to produce a finite-length, causal filter. The choice of window trades off between mainlobe width (transition band) and sidelobe suppression (stopband attenuation).
3.1.3 Optimal (Parks–McClellan) design
The Parks–McClellan algorithm (also known as the Remez exchange algorithm) produces equiripple FIR filters that minimize the maximum deviation from the desired response in the passband and stopband. It yields the lowest possible order for a given set of specifications and is widely used for high-performance filters.
3.2 Infinite impulse response (IIR) filters
IIR filters have an impulse response that theoretically continues indefinitely, due to feedback. They can achieve sharp frequency transitions with much lower order than FIR filters, but they may be unstable and typically exhibit nonlinear phase.
3.2.1 Analog prototype transformations
Classical analog lowpass filters serve as prototypes for digital IIR designs. The most common types are Butterworth (maximally flat passband), Chebyshev (equiripple in either passband or stopband), and elliptic (equiripple in both). Their analog transfer functions have known pole-zero locations.
###### 3.2.1.1 Butterworth
Butterworth filters have a monotonically decreasing magnitude response with no ripple. The –3 dB cutoff frequency defines the passband edge. The frequency response is maximally flat at dc.
###### 3.2.1.2 Chebyshev (Type I and II)
Chebyshev Type I filters have equiripple in the passband and monotonic stopband; Type II (inverse Chebyshev) has a monotonic passband and equiripple stopband. They achieve a sharper rolloff than Butterworth for the same order.
###### 3.2.1.3 Elliptic
Elliptic (Cauer) filters have equiripple in both passband and stopband, yielding the steepest transition for a given order. Their design involves elliptic functions and is more complex, but they are very efficient.
3.2.2 Bilinear transform and impulse invariance
The bilinear transform maps the analog frequency variable \(s\) to the digital variable \(z\) via \(s = \frac{2}{T} \frac{1-z^{-1}}{1+z^{-1}}\). It warps the frequency axis but preserves stability and maps the entire analog spectrum to the unit circle. Impulse invariance directly samples the analog impulse response, but it can cause aliasing unless the analog filter is bandlimited.
3.3 Filter implementation structures
A given transfer function can be realized in various structures, each with different sensitivity to coefficient quantization and computational efficiency.
3.3.1 Direct forms (I, II)
Direct Form I implements the difference equation directly, requiring separate delay lines for input and output. Direct Form II uses a single delay line by factoring the denominator first, reducing memory requirements. Both are sensitive to quantization errors in high-order filters.
3.3.2 Cascade and parallel forms
Higher-order filters are often decomposed into second-order sections (biquads) in cascade or parallel. Cascade form minimizes coefficient sensitivity and is numerically robust. Parallel form combines second-order sections whose outputs are summed, useful for certain applications.
3.3.3 Quantization effects in fixed-point arithmetic
When filters are implemented with fixed-point arithmetic, coefficient quantization alters the pole-zero locations, potentially causing instability or passband distortion. Limit cycles (oscillations due to rounding) and overflow nonlinearities must be managed through scaling and careful word-length choice.
4 Advanced topics and applications
4.1 Multirate digital signal processing
Multirate DSP involves changing the sampling rate of a signal, enabling efficient processing of wideband or narrowband components.
4.1.1 Decimation and interpolation
Decimation reduces the sampling rate by removing samples after low-pass filtering to avoid aliasing. Interpolation increases the sampling rate by inserting zeros and then low-pass filtering to remove spectral images. Both operations can be combined in rational sampling rate changes.
4.1.2 Polyphase filter banks
Polyphase decomposition splits a filter into multiple phases that can be processed at a lower rate, greatly reducing computational cost. Polyphase filter banks are used in subband coding, transmultiplexers, and wavelet transforms.
4.1.3 Sampling rate conversion
To convert between arbitrary sampling rates (e.g., 44.1 kHz to 48 kHz), a combination of interpolation and decimation is used, often implemented efficiently with a polyphase structure or a fractional delay filter.
4.2 Adaptive filtering
Adaptive filters adjust their coefficients in real time to minimize an error signal, enabling them to track changing environments.
4.2.1 Least mean squares (LMS) algorithm
LMS is a stochastic gradient algorithm that updates filter coefficients using the instantaneous error: \(w_{n+1} = w_n + \mu e[n] x[n]\). It is simple, robust, and widely used, but its convergence speed depends on the step size \(\mu\) and the eigenvalue spread of the input correlation matrix.
4.2.2 Recursive least squares (RLS) algorithm
RLS minimizes the sum of squared errors from the present back to the initial time, using a recursive update that involves the inverse correlation matrix. It converges faster than LMS but has higher computational cost (\(O(N^2)\) per sample) and may be prone to numerical instability.
4.2.3 Applications in noise cancellation and echo suppression
Adaptive filters are key in acoustic echo cancellation (e.g., in hands-free telephony), active noise control (canceling ambient noise with an inverted anti-noise signal), and adaptive beamforming for microphone arrays.
4.3 Spectral estimation
Spectral estimation aims to identify the frequency content of a signal from a finite number of samples.
4.3.1 Nonparametric methods (periodogram, Welch)
The periodogram computes the squared magnitude of the DFT, but it has high variance and poor resolution for short data lengths. Welch’s method averages overlapping, windowed periodograms to reduce variance at the cost of lower frequency resolution.
4.3.2 Parametric methods (AR, MA, ARMA models)
Parametric methods assume the signal follows a known model (autoregressive, moving average, or autoregressive moving average). The model parameters are estimated from the data, and the spectrum is computed from the model’s transfer function. AR models are especially popular for high-resolution peak detection.
4.3.3 Subspace methods (MUSIC, ESPRIT)
Subspace methods decompose the data covariance matrix into signal and noise subspaces. MUSIC (Multiple Signal Classification) uses the orthogonality between signal and noise subspaces to estimate frequencies of sinusoids in noise. ESPRIT (Estimation of Signal Parameters via Rotational Invariance Techniques) exploits shift invariance and offers lower computational load.
4.4 Digital signal processors and real-time implementation
Real-time DSP requires hardware and software optimized for repetitive arithmetic operations, especially multiply-accumulate (MAC) instructions.
4.4.1 Fixed-point vs. floating-point architectures
Fixed-point processors use integer arithmetic with implicit scaling, offering low power and cost but requiring careful scaling to avoid overflow. Floating-point processors simplify programming by providing a wider dynamic range, at the expense of higher power and cost. Many modern DSPs and FPGAs support both modes.
4.4.2 Pipelining and parallel processing
Pipelining stages memory accesses and arithmetic operations to increase throughput. Parallel processing (multiple MAC units, SIMD instructions, or multicore DSPs) further accelerates computations such as FFTs and filter banks.
4.4.3 Code optimization for embedded systems
Optimizing DSP code involves minimizing memory access, using circular buffers for delay lines, unrolling loops, and exploiting instruction-level parallelism. Compilers for embedded DSPs often provide intrinsic functions for multiply-accumulate and saturating arithmetic.