# Discrete Hilbert Transform

*The Hilbert transform is the workhorse of numerical Kramers–Kronig and
of every analytic-signal computation in DSP. Three standard
implementations cover the field: FFT-based, FIR, and IIR allpass.*

This chapter uses the **engineering convention** `e^{+jωt}`, so the
analytic signal of a real x[n] is x[n] + j ℋ[x][n], with ℋ the discrete
Hilbert transform. The Fourier-domain multiplier is −j sgn(ω) in this
convention (the physics convention gives +i sgn(ω); see Chapter 16).

## The three standard implementations

### (a) FFT multiplier (Marple, 1999)

The most common numerical route. Given a real sequence x[n] of length
N, its DFT X[k] satisfies Hermitian symmetry. The analytic signal is
obtained by zeroing the negative-frequency half and doubling the
positive half:

```
X_a[k] =  X[0]                     k = 0
X_a[k] =  2 X[k]                   k = 1, 2, …, N/2 − 1
X_a[k] =  X[N/2]                   k = N/2        (if N even)
X_a[k] =  0                        k = N/2 + 1, …, N − 1
```

An inverse FFT of X_a gives the complex analytic signal x_a[n], and
its imaginary part is ℋ[x][n].

Cost: **O(N log N)**. Exact under DFT periodicity; no filter design,
no transient. The only implicit assumption is that x[n] is a single
period of a length-N periodic signal — which is exactly what the DFT
already assumes. MATLAB's `hilbert()` and NumPy's
`scipy.signal.hilbert()` are this algorithm.

Edge effects arise when x[n] is a finite window of a longer signal:
DFT periodisation wraps the end to the start, distorting the Hilbert
transform near n = 0 and n = N − 1. Mitigations are overlap-add on
blocks, zero-padding before the DFT, or windowing.

### (b) FIR Hilbert filter

An odd-length Type-III FIR filter with **anti-symmetric** impulse
response

```
h[n] = −h[N − 1 − n]
```

can realise a Hilbert transformer. The impulse response of the ideal
Hilbert filter is

```
h_ideal[n] =  2/(π n)    for odd n ≠ 0
h_ideal[n] =  0          for even n  (including n = 0)
```

which is anti-symmetric about n = 0 and decays as 1/n. Truncating to
N taps with N typically 64–256 and applying a window (Blackman,
Kaiser) or the Parks–McClellan equiripple method gives a stable FIR
Hilbert transformer with passband error of roughly 0.01 dB over
[ω_1, π − ω_1] for N = 128 and modest stopband transitions near DC
and Nyquist.

Group delay is (N − 1)/2 samples, a fixed latency. Implementation
cost is **O(N)** per sample, or O(N log N) per block via overlap-save.

### (c) IIR allpass pair (Ansari–Harris, 1987)

For low-latency real-time audio the preferred implementation is a pair
of allpass filters A_0(z) and A_1(z) designed so that their phase
responses differ by 90° over a specified band. The output of A_1 is
the Hilbert transform of the output of A_0 up to a small band-edge
error.

Each allpass is a cascade of second-order sections of the form

```
A(z) = (a + z^{−2}) / (1 + a z^{−2})
```

with coefficients a chosen to place the phase-quadrature frequencies
where required. Regalia–Mitra gives a design procedure based on
half-band filter factorisation; 4+4 allpass sections typically give
90° ± 0.5° from 0.02 f_s to 0.48 f_s.

Latency is a few samples — dramatically less than the FIR route — at
the price of phase-only (not magnitude-preserving) quadrature.

## Numerical issues

- **Gibbs phenomenon near band edges.** The ideal multiplier −j sgn(ω)
  is discontinuous at ω = 0 and ω = π, so any finite approximation
  rings. Windowing (FIR) or cos² tapers (FFT) trade ripple for
  transition width.
- **Latency.** FFT: one block length, typically hundreds of samples.
  FIR: (N − 1)/2 samples. IIR allpass: a handful of samples. Picking
  the right implementation is almost always a latency decision.
- **DC and Nyquist handling.** The multiplier is zero at DC and
  Nyquist, so any DC offset in x[n] is simply dropped by the Hilbert
  transform. If the K-K application needs a separate DC subtraction
  (see Chapter 15), handle it before the Hilbert step.
- **Aliasing after non-linear post-processing.** Squaring the analytic
  signal (for envelope power) doubles the bandwidth; upsample before
  non-linearities if the baseband x[n] is close to Nyquist.

## Applications

### Single-sideband (SSB) modulation

The Weaver and Hartley SSB modulators use the analytic signal to
generate a carrier-modulated waveform with only one sideband, for
bandwidth-efficient voice transmission. The Hartley modulator computes

```
s_USB(t) = x(t) cos ω_c t − ℋ[x](t) sin ω_c t
```

with ℋ the Hilbert transform. Sideband rejection is set by the
accuracy of the 90° quadrature — around 40 dB rejection for a 64-tap
FIR, 60+ dB for a well-designed IIR allpass pair.

### Envelope detection

In audio compressors, radar pulse detection, and biomedical signal
processing, the instantaneous amplitude

```
A(t) = |x(t) + j ℋ[x](t)|
```

gives a clean time-varying envelope without rectification artefacts.

### Instantaneous frequency and phase

The time derivative of arg(x + jℋx) is the instantaneous angular
frequency, used in FM demodulation, pitch tracking, and empirical mode
decomposition (EMD) / Hilbert–Huang transform.

### Numerical K-K

The direct-method K-K integral

```
χ″(ω) = −(1/π) P ∫_{−∞}^{∞}  χ′(ω′) / (ω′ − ω) dω′
```

is — up to a sign — the Hilbert transform of χ′ evaluated at ω.
Implementing K-K numerically is therefore implementing the Hilbert
transform of the measured (and extrapolated, see Chapter 29) data,
with the multiplier-in-FFT method being the workhorse. A typical code
path is: interpolate measured data to a uniform ω grid → FFT →
multiply by −j sgn(ω) → inverse FFT → imaginary part.

## Pseudocode

```python
import numpy as np

def hilbert_fft(x):
    """Analytic signal of a real array x, engineering convention."""
    N = len(x)
    X = np.fft.fft(x)
    H = np.zeros(N, dtype=complex)
    H[0] = 1
    if N % 2 == 0:
        H[1:N//2] = 2
        H[N//2]   = 1
    else:
        H[1:(N+1)//2] = 2
    return np.fft.ifft(X * H)

def hilbert_transform(x):
    """Discrete Hilbert transform, engineering convention."""
    return np.imag(hilbert_fft(x))
```

The Ansari IIR variant and the Parks–McClellan FIR design are one-line
calls in `scipy.signal` (`firwin`, `remez`, `iirfilter` with a suitable
design routine).

## See also

- [19_minimum_phase_DSP.md](19_minimum_phase_DSP.md)
- [27_audio_room_acoustics.md](27_audio_room_acoustics.md)
- [29_numerical_implementation.md](29_numerical_implementation.md)

## References

- A. V. Oppenheim and R. W. Schafer, *Discrete-Time Signal Processing*,
  3rd ed., Pearson, 2010 — Ch. 12.
- S. L. Marple, "Computing the Discrete-Time Analytic Signal via FFT",
  *IEEE Transactions on Signal Processing*, vol. 47, no. 9,
  pp. 2600–2603, 1999.
- R. Ansari, "IIR discrete-time Hilbert transformers", *IEEE Trans.
  ASSP*, vol. 35, no. 8, pp. 1116–1119, 1987.
- P. A. Regalia and S. K. Mitra, "Tunable digital frequency response
  equalization filters", *IEEE Trans. ASSP*, vol. 35, pp. 118–120,
  1987.