A QPSK modem with a measured BER curve
Build a QPSK transmitter and receiver in NumPy with root-raised-cosine pulse shaping, an AWGN channel and a matched filter, then measure the bit error rate from 0 to 10 dB Eb/N0 against theory.
| Runs on | Laptop CPU (macOS, Linux or Windows), Python and NumPy, under 15 seconds |
|---|---|
| Time | about 2 hours |
| Stack | python, numpy |
What you will build
This project builds a complete QPSK modem in one NumPy file and measures its bit error rate (BER). The transmitter maps random bits to Gray-coded QPSK symbols and shapes them with a root-raised-cosine (RRC) filter. The channel delays the signal by a few samples and adds white Gaussian noise. The receiver runs a matched filter, takes one sample per symbol at the correct sampling instant and decides each bit from the sign of the real or imaginary part.
The program then sweeps Eb/N0 from 0 dB to 10 dB, counts errors, and prints the measured BER next to the textbook result for Gray-coded QPSK in additive white Gaussian noise (AWGN). In our run the two agree to within the statistical spread of the error counts at every point. At 6 dB the modem makes 482 errors in 200,000 bits (BER 2.41 × 10⁻³ against 2.39 × 10⁻³ from theory).
A last check samples at the wrong sampling instant on purpose. Half a symbol of timing error turns a working link into a coin toss. That is why timing recovery exists, and it is the first thing you will need on a real radio.
The lessons on digital modulation, pulse shaping and matched filters and noise, SNR and capacity cover the pieces one at a time. The modulation and BER cheat sheet has the formulas in one place.
What you need
- Python 3.10 or later and NumPy. We ran it with Python 3.12 and NumPy 2.5.4 on an Apple Silicon Mac mini.
- Nothing else. There is no SciPy, no plotting library and no radio. All randomness comes from one seeded generator.
- About 2 hours to read the design, run the program and work through the variations.
Design
Symbols and energy
Each pair of bits picks one of four points, (±1 ± j)/√2. The first bit sets the sign of the real part and the second bit sets the sign of the imaginary part. That is Gray coding: neighbouring points differ in one bit, so the most likely symbol error costs one bit, not two. Every symbol has energy Es = 1, so each bit has Eb = 0.5.
The symbol rate is 1 Msym/s and there are 8 samples per symbol, so the sample rate is 8 MHz. The numbers only matter for labels. The BER depends on Eb/N0 and nothing else.
Pulse shaping and the matched filter
The transmitter puts each symbol on an impulse every 8 samples and convolves the impulse train with an RRC filter of roll-off 0.35, 10 symbols long (81 taps). The filter is scaled to unit energy. The occupied bandwidth is (1 + 0.35) × 1 MHz = 1.35 MHz.
The receiver convolves with the same filter reversed in time. That is the matched filter, the linear filter that maximises the SNR at the sampling instant in white noise. The two RRC filters in cascade form a raised-cosine pulse, which is zero at every other symbol instant. Sampling at the right moment therefore sees one symbol and no intersymbol interference. The program checks this with a noise-free run.
Noise scaling
Getting the noise level right is the step most first modems get wrong. With unit-energy symbols and a unit-energy filter, the transmitted energy per symbol, summed over its samples, is Es = 1. Complex white noise with variance N0 per sample then gives the correct Es/N0, and the matched filter does not change the ratio. The code converts Eb/N0 in dB to Es/N0 = 2 Eb/N0 for two bits per symbol, sets N0 = 1 / (Es/N0) and splits N0 equally between the real and imaginary parts.
Timing at a known offset
This modem does not recover timing. It knows the delay. The RRC filter has 81 taps, so each filter delays the signal by 40 samples, and the pair by 80 samples. The channel adds 3 more. The receiver takes its first sample at index 83 and then every 8th sample. A real receiver has to estimate that offset, and the last check shows how quickly the BER degrades when the offset is wrong.
The theoretical curve
For Gray-coded QPSK in AWGN with coherent detection, each rail is an independent BPSK decision. The bit error probability is Pb = Q(√(2 Eb/N0)), where Q is the tail probability of a standard normal distribution. The code computes Q from math.erfc as Q(x) = 0.5 erfc(x/√2).
Code
Save this as qpsk_modem.py. It is one block: filter design, modulator, channel, receiver and the three checks.
# qpsk_modem.py: QPSK modem with RRC pulse shaping, AWGN and a measured BER curve.
import math
import numpy as np
SPS = 8 # samples per symbol
ALPHA = 0.35 # RRC roll-off
SPAN = 10 # RRC length in symbols
CHANNEL_DELAY = 3 # extra delay in samples, known to the receiver
SYMBOL_RATE = 1e6 # 1 Msym/s, so the sample rate is 8 MHz
rng = np.random.default_rng(2026)
def rrc_taps(sps: int, alpha: float, span: int) -> np.ndarray:
"""Root-raised-cosine impulse response with unit energy."""
t = np.arange(-span * sps // 2, span * sps // 2 + 1) / sps # time in symbols
h = np.zeros_like(t)
for i, ti in enumerate(t):
if abs(ti) < 1e-12:
h[i] = 1.0 - alpha + 4 * alpha / math.pi
elif abs(abs(4 * alpha * ti) - 1.0) < 1e-12:
h[i] = (alpha / math.sqrt(2)) * ((1 + 2 / math.pi) * math.sin(math.pi / (4 * alpha))
+ (1 - 2 / math.pi) * math.cos(math.pi / (4 * alpha)))
else:
num = math.sin(math.pi * ti * (1 - alpha)) + 4 * alpha * ti * math.cos(math.pi * ti * (1 + alpha))
den = math.pi * ti * (1 - (4 * alpha * ti) ** 2)
h[i] = num / den
return h / np.sqrt(np.sum(h ** 2))
def qpsk_map(bits: np.ndarray) -> np.ndarray:
"""Gray-coded QPSK, unit symbol energy: bit 0 -> +1, bit 1 -> -1 on each rail."""
b = bits.reshape(-1, 2)
return ((1 - 2 * b[:, 0]) + 1j * (1 - 2 * b[:, 1])) / np.sqrt(2)
def qpsk_demap(z: np.ndarray) -> np.ndarray:
"""Hard decisions: the sign of each rail."""
return np.column_stack([(z.real < 0), (z.imag < 0)]).astype(np.int8).ravel()
def transmit(symbols: np.ndarray, h: np.ndarray) -> np.ndarray:
up = np.zeros(len(symbols) * SPS, dtype=np.complex128)
up[::SPS] = symbols # impulse train, one impulse per symbol
return np.convolve(up, h) # pulse shaping
def receive(rx: np.ndarray, h: np.ndarray, n_sym: int, offset: int) -> np.ndarray:
mf = np.convolve(rx, h[::-1].conj()) # matched filter
return mf[offset: offset + n_sym * SPS: SPS] # one sample per symbol
def run_link(ebn0_db: float, n_bits: int, h: np.ndarray, timing_error: int = 0) -> tuple[int, int]:
bits = rng.integers(0, 2, n_bits, dtype=np.int8)
sym = qpsk_map(bits)
tx = transmit(sym, h)
tx = np.concatenate([np.zeros(CHANNEL_DELAY, dtype=np.complex128), tx]) # channel delay
es_n0 = 10 ** (ebn0_db / 10) * 2 # 2 bits per symbol, Es = 1
n0 = 1.0 / es_n0 # noise power per complex sample
noise = np.sqrt(n0 / 2) * (rng.standard_normal(len(tx)) + 1j * rng.standard_normal(len(tx)))
rx = tx + noise
offset = (len(h) - 1) + CHANNEL_DELAY + timing_error # TX filter delay + MF delay + channel
z = receive(rx, h, len(sym), offset)
errors = int(np.count_nonzero(qpsk_demap(z) != bits))
return errors, n_bits
def q_func(x: np.ndarray) -> np.ndarray:
return 0.5 * np.array([math.erfc(v / math.sqrt(2)) for v in np.atleast_1d(x)])
h = rrc_taps(SPS, ALPHA, SPAN)
print(f"RRC: {len(h)} taps, alpha {ALPHA}, {SPS} samples/symbol, sample rate {SYMBOL_RATE * SPS / 1e6:.0f} MHz")
# Check 1: with no noise the RRC pair has no intersymbol interference.
bits = rng.integers(0, 2, 2000, dtype=np.int8)
sym = qpsk_map(bits)
z = receive(np.concatenate([np.zeros(CHANNEL_DELAY), transmit(sym, h)]), h, len(sym),
len(h) - 1 + CHANNEL_DELAY)
print(f"noise-free: max |z - s| = {np.max(np.abs(z - sym)):.2e}, bit errors = {np.count_nonzero(qpsk_demap(z) != bits)}")
# Check 2: BER against Eb/N0, stop at 100 errors or 20 million bits.
print("Eb/N0 (dB) bits errors BER sim BER theory")
results = []
for ebn0_db in range(0, 11):
errs, total = 0, 0
while errs < 100 and total < 20_000_000:
e, n = run_link(ebn0_db, 200_000, h)
errs, total = errs + e, total + n
ber = errs / total
theory = q_func(np.sqrt(2 * 10 ** (ebn0_db / 10)))[0]
results.append((ebn0_db, ber, theory))
print(f"{ebn0_db:>6} {total:>9} {errs:>6} {ber:.3e} {theory:.3e}")
# Check 3: sample at the wrong instant (Eb/N0 = 6 dB).
print("timing error (samples) BER at 6 dB")
for err in [0, 1, 2, 3, 4]:
e, n = run_link(6.0, 400_000, h, timing_error=err)
print(f"{err:>10} {e / n:.3e}")
worst = max(abs(b / t - 1) for _, b, t in results)
print(f"largest relative gap between simulated and theoretical BER: {100 * worst:.1f} %")
The sweep stops at each Eb/N0 once it has counted 100 errors or simulated 20 million bits. One hundred errors give a relative standard error of about 10 % on the BER estimate (one over the square root of the error count). The 20 million bit cap keeps the 10 dB point under a few seconds.
Run
python -m pip install numpy
python qpsk_modem.py
On our machine the whole run takes about 8 s.
Test
This is the output from our run. With the seed fixed at 2026 and the same NumPy version you should see the same numbers.
RRC: 81 taps, alpha 0.35, 8 samples/symbol, sample rate 8 MHz
noise-free: max |z - s| = 1.80e-02, bit errors = 0
Eb/N0 (dB) bits errors BER sim BER theory
0 200000 15857 7.928e-02 7.865e-02
1 200000 11082 5.541e-02 5.628e-02
2 200000 7277 3.639e-02 3.751e-02
3 200000 4472 2.236e-02 2.288e-02
4 200000 2520 1.260e-02 1.250e-02
5 200000 1185 5.925e-03 5.954e-03
6 200000 482 2.410e-03 2.388e-03
7 200000 174 8.700e-04 7.727e-04
8 600000 104 1.733e-04 1.909e-04
9 3000000 100 3.333e-05 3.363e-05
10 20000000 85 4.250e-06 3.872e-06
timing error (samples) BER at 6 dB
0 2.452e-03
1 6.370e-03
2 3.231e-02
3 1.130e-01
4 2.498e-01
largest relative gap between simulated and theoretical BER: 12.6 %
Read the output from the top.
The noise-free line. Zero bit errors, and the largest distance between a received sample and its transmitted symbol is 0.018. That small residue is intersymbol interference from cutting the RRC filter off at 10 symbols. A longer filter makes it smaller. If this number is near 0.5 or there are bit errors, the sampling offset is wrong or the filter is not normalised.
The BER table. The largest relative gap is 12.6 %, at 7 dB, where the estimate rests on 174 errors. With that many errors the standard error is about 7.6 %, so the gap is within two standard errors. At 10 dB, 85 errors give a standard error of about 11 %, and the measured BER sits 10 % above theory. If your curve sits a constant distance to the right of theory, check the noise scaling: a factor of 2 in N0 shifts the whole curve by 3 dB.
The timing sweep. At 6 dB the BER rises from 2.45 × 10⁻³ with correct timing to 6.37 × 10⁻³ with one sample (an eighth of a symbol) of error. Two samples of error raise it to 3.23 × 10⁻², and four samples, half a symbol, give 0.25. At half a symbol the sample falls between two symbols, and the two pulses contribute almost equally.
On a real radio
A real link breaks three assumptions this program makes, and each one needs a new block in the receiver.
- Unknown timing. The receiver does not know the delay. A preamble with good autocorrelation finds the frame start. A timing error detector, such as Gardner’s, then tracks the best sampling instant as the two clocks drift. Without it, the timing sweep above shows what happens.
- Carrier offset and phase. The transmitter and receiver oscillators differ by some parts per million. At 2.4 GHz, 10 ppm is 24 kHz of carrier frequency offset, which spins the constellation 0.15 rad per symbol at 1 Msym/s. A frequency estimate from the preamble and a phase-locked loop (a Costas loop for QPSK) remove it. QPSK also has a fourfold phase ambiguity, which a known preamble or differential coding resolves.
- Gain and the analogue front end. The received level is unknown and changes. Hard decisions on QPSK do not need the amplitude, but higher-order QAM does, so real receivers include automatic gain control. DC offset and IQ imbalance from the mixer also distort the constellation.
To run this modem over the air, write the transmitted samples to a file as complex64, play them through an SDR transmitter at 8 MHz, record the received samples and feed them to the receiver half of the program with the three blocks above added. Do this into a dummy load or over a cable with attenuators unless you hold a licence for the frequency you use. The lesson on recording and reading IQ shows the file side, and the first GNU Radio flowgraph lesson shows the radio side. GNU Radio’s digital module has ready-made symbol synchronisation and Costas loop blocks that implement the first two items.
Things to try
- Change the roll-off. Set
ALPHA = 0.2and thenALPHA = 0.8. The BER at correct timing does not change, because the matched filter output SNR does not depend on the pulse shape. The timing sweep does change. In our runs one sample of timing error gave a BER of 7.70 × 10⁻³ with a roll-off of 0.2 and 4.40 × 10⁻³ with 0.8, against 6.37 × 10⁻³ at 0.35. A smaller roll-off has larger tails, so the same timing error costs more. - Shorten the filter. Set
SPAN = 4(33 taps). The noise-free residue grew from 0.018 to 0.120 in our run because the truncated pulse is no longer close to Nyquist. The high end of the curve lifts away from theory: at 10 dB the measured BER was 9.18 × 10⁻⁶ against 3.87 × 10⁻⁶. - Break the matched filter. Replace the receive filter with a plain moving average over 8 samples. The curve shifts to the right. Measure by how many dB.
- Remove Gray coding. Map bit pairs 00, 01, 10 and 11 to the four points in order round the circle. At low Eb/N0 the BER rises, because adjacent-symbol errors now cost two bits about half the time.
- Add a carrier offset. Multiply the received samples by exp(j 2π f t) with f = 1 kHz at the 8 MHz sample rate before the matched filter. The constellation rotates by 0.36° per symbol and the BER collapses once it has turned 45°. Then add a simple decision-directed phase tracker.
- Measure an eye diagram. Reshape the matched filter output into rows of 16 samples (two symbols) and print the spread of values at each column. The narrowest spread marks the best sampling offset. The eye diagram is how you would find the offset by inspection.
Summary
- A QPSK modem is four pieces: Gray mapping, RRC pulse shaping, a matched filter and a sign slicer. Each fits in a few lines of NumPy.
- Scale the noise from Eb/N0 with care. With unit-energy symbols and filters, N0 per complex sample is 1 / (2 Eb/N0) for QPSK.
- The measured BER matched Pb = Q(√(2 Eb/N0)) from 7.9 × 10⁻² at 0 dB to about 4 × 10⁻⁶ at 10 dB, within the spread expected from the error counts.
- Count enough errors. One hundred errors give about 10 % relative uncertainty on a BER estimate.
- Timing matters. One eighth of a symbol of error multiplied the BER at 6 dB by about 2.6, and half a symbol made the link useless.
Background: the modulation and synchronisation and estimation topic hubs, and the OFDM link project, which reuses the same noise scaling on a multipath channel.