Cyclic Autocorrelation Function (CAF): How BPSK and QPSK Differ in Cyclostationary Analysis

In the previous article, we established that a cyclostationary signal has an autocorrelation that is periodic in time $t$ at a fixed lag $\tau$. We also introduced the Spectral Correlation Function (SCF) as a way to visualise this periodic structure. The missing piece is the tool that sits between the raw time-varying autocorrelation and the SCF: the Cyclic Autocorrelation Function (CAF).

The CAF is the function that extracts the periodic component of $R_x(t, \tau)$ with respect to $t$, at each lag $\tau$. It is the foundation that the SCF is built on, and it is the tool that makes modulation classification possible — as we will see, BPSK and QPSK produce fundamentally different CAF signatures at a specific cycle frequency that makes them directly distinguishable.

What is the Cyclic Autocorrelation Function?

Recall from the previous article that for a cyclostationary signal, the time-varying autocorrelation:

$R_x(t, \tau) = \mathbb{E}[x(t + \tau/2)\, x^*(t - \tau/2)]$

is periodic in $t$ with period $T_0$. Since it is periodic in $t$, it has a Fourier series expansion. The CAF is defined as the Fourier series coefficient of $R_x(t, \tau)$ at cycle frequency $\alpha$:

$R_x^\alpha(\tau) = \lim_{T \to \infty} \frac{1}{T} \int_{-T/2}^{T/2} x\!\left(t + \frac{\tau}{2}\right) x^*\!\left(t - \frac{\tau}{2}\right) e^{-j2\pi\alpha t}\, dt$

where:

• $\tau$ is the lag — how far apart the two samples being compared are
• $\alpha$ is the cycle frequency — the specific periodic component being extracted from $R_x(t, \tau)$
• $e^{-j2\pi\alpha t}$ is a complex tone at frequency $\alpha$. Multiplying by this tone and averaging over $t$ extracts only the component of the autocorrelation that oscillates at exactly $\alpha$. All other components average to zero
• The limit over $T$ means we average over a long observation window to suppress noise

The result $R_x^\alpha(\tau)$ is a single complex-valued function of lag $\tau$, computed at one specific cycle frequency $\alpha$.

The key test for cyclostationarity: if $R_x^\alpha(\tau) \neq 0$ for some $\alpha \neq 0$, the signal is cyclostationary at that cycle frequency. Noise produces $R_x^\alpha(\tau) = 0$ for all $\alpha \neq 0$, because noise has no periodic structure to extract.

The discrete-time estimator

In practice the signal is sampled, so the continuous integral becomes a sum. The discrete-time estimator of the CAF is:

$\hat{R}_x^\alpha(\tau) = \frac{1}{N} \sum_{n=0}^{N-1} x(n+\tau)\, x^*(n)\, e^{-j2\pi\alpha n / f_s}$

where:

• $n$ is the sample index
• $N$ is the total number of samples — a larger $N$ gives a more accurate estimate
• $f_s$ is the sample rate in Hz
• $x(n+\tau)$ and $x^*(n)$ are the signal and its conjugate, separated by lag $\tau$
• $e^{-j2\pi\alpha n / f_s}$ is the complex tone at cycle frequency $\alpha$, evaluated at each sample

This estimator uses a one-sided lag convention, where one sample is shifted by $\tau$ and the other is not. The standard definition from Gardner uses a symmetric convention where one sample is shifted by $+\tau/2$ and the other by $-\tau/2$. The two are equivalent up to a phase factor that does not affect the magnitude, which is what we plot.

BPSK and QPSK signals

Before computing the CAF, it helps to see what the two signals look like.

Figure 1 — BPSK switches between two phases (0° and 180°). QPSK switches between four phases (45°, 135°, 225°, 315°). Both use the same carrier frequency and symbol rate.

The signal parameters used throughout this article are: sample rate $f_s = 1000$ Hz, carrier frequency $f_c = 50$ Hz, and symbol rate $R_s = 5$ symbols per second, giving a symbol period $T_s = 0.2$ seconds and 200 samples per symbol.

The two signals look visually similar in the time domain. A standard power spectrum would also show two similar peaks at $\pm f_c$. The CAF is the tool that reveals the difference between them.

CAF at $\alpha = 0$: ordinary autocorrelation

When $\alpha = 0$, the complex exponential $e^{-j2\pi \cdot 0 \cdot n} = 1$, so the CAF reduces to the ordinary autocorrelation:

$R_x^0(\tau) = \frac{1}{N} \sum_n x(n+\tau)\, x^*(n) = R_x(\tau)$

Figure 2 — CAF at $\alpha = 0$ for BPSK (left) and QPSK (right). Both show a peak at $\tau = 0$ — this is just the ordinary autocorrelation and does not distinguish between the two modulations.

At $\alpha = 0$, both BPSK and QPSK respond with a peak at $\tau = 0$. This tells us there is signal power, but nothing about the modulation type. To distinguish between them we need to evaluate the CAF at non-zero cycle frequencies.

CAF at $\alpha = R_s$: the baud-rate feature

The symbol rate $R_s = 1/T_s$ is a cycle frequency present in both BPSK and QPSK. Both signals have symbols that repeat with period $T_s$, so both exhibit cyclic structure at $\alpha = R_s$.

Figure 3 — CAF at $\alpha = R_s = 5$ Hz for BPSK (left) and QPSK (right). Both show a peak at $\tau = 0$, confirming that both signals have cyclic structure at the baud rate. This cycle frequency alone cannot distinguish between them.

The peak at $\tau = 0$ and at $\tau = \pm T_s$ (lag equal to one symbol period) confirms that both signals are cyclostationary at the baud rate. However since both respond, this cycle frequency by itself is not sufficient to classify the modulation type.

CAF at $\alpha = 2f_c$: the key discriminator

This is the result that the article has been building towards. BPSK has non-zero conjugate cyclostationarity — it produces a peak in the CAF at cycle frequency $\alpha = 2f_c$. QPSK, and all other complex-valued symmetric constellations such as 16-QAM and 8-PSK, have zero conjugate cyclostationarity — they produce no peak at $\alpha = 2f_c$.

Figure 4 — CAF at $\alpha = 2f_c = 100$ Hz for BPSK (left) and QPSK (right). BPSK shows a clear peak. QPSK is near zero. This single cycle frequency is the discriminator between the two modulations.

Why does BPSK produce a peak at $\alpha = 2f_c$?

BPSK symbols are real-valued: $a_k \in \{-1, +1\}$. Squaring a real symbol gives $a_k^2 = 1$ — a constant. This means when you compute the product $x(t) \cdot x(t+\tau)$ without taking the conjugate (the conjugate CAF), the random symbol component cancels and what remains is the carrier term $\cos(2\pi f_c t) \cdot \cos(2\pi f_c (t+\tau))$. Expanding this product gives a term at frequency $2f_c$, which produces the peak at $\alpha = 2f_c$.

Why does QPSK produce no peak at $\alpha = 2f_c$?

QPSK symbols are complex-valued, drawn from $\{e^{j\pi/4}, e^{j3\pi/4}, e^{j5\pi/4}, e^{j7\pi/4}\}$. For these symbols, $a_k^2$ is not constant — it takes values $\{e^{j\pi/2}, e^{j3\pi/2}, ...\}$ which are complex and average to zero over many symbols. The cancellation that happens for BPSK does not happen for QPSK. The conjugate spectral correlation function is zero for QPSK and all other higher-order PSK and QAM constellations.

The cyclic profile: a summary view

Rather than testing one $\alpha$ at a time, we can sweep $\alpha$ from 0 to $2f_c + R_s$ and plot $|R_x^\alpha(0)|$ at each value. This gives the cyclic profile — a compact view of which cycle frequencies are active for each signal.

Figure 5 — Cyclic profile for BPSK (top) and QPSK (bottom). BPSK shows peaks at harmonics of $R_s$ and at $2f_c$. QPSK shows peaks at harmonics of $R_s$ only. The presence or absence of the $2f_c$ peak is the modulation fingerprint.

The cyclic profile makes the difference between the two modulations immediately visible:

A received signal showing a peak at $\alpha = 2f_c$ is BPSK. A signal showing no peak there, only at baud-rate harmonics, is QPSK or a higher-order constellation. This is the basis of second-order modulation classification in SIGINT.

Python code

The complete code to reproduce all figures is provided in the Jupyter notebook below. Here is the CAF estimator function used throughout:

import numpy as np
import matplotlib.pyplot as plt

# ── Parameters ────────────────────────────────────────────────────────────────
fs        = 1000       # sample rate (Hz)
fc        = 50         # carrier frequency (Hz)
Rs        = 5          # symbol rate (symbols/sec)
sps       = int(fs/Rs) # samples per symbol = 200
n_symbols = 100        # number of symbols
N         = n_symbols * sps
np.random.seed(42)
t = np.arange(N) / fs

# ── Generate BPSK ─────────────────────────────────────────────────────────────
bits_bpsk = np.random.choice([-1, 1], size=n_symbols)
bpsk      = np.repeat(bits_bpsk, sps) * np.cos(2*np.pi*fc*t)

# ── Generate QPSK ─────────────────────────────────────────────────────────────
# Four phases: 45°, 135°, 225°, 315°
bits_qpsk = np.random.choice([0, 1, 2, 3], size=n_symbols)
phases    = np.array([np.pi/4, 3*np.pi/4, 5*np.pi/4, 7*np.pi/4])
qpsk_syms = np.exp(1j * phases[bits_qpsk])
qpsk_bb   = np.repeat(qpsk_syms, sps)
qpsk      = np.real(qpsk_bb) * np.cos(2*np.pi*fc*t) \
           - np.imag(qpsk_bb) * np.sin(2*np.pi*fc*t)

# ── CAF estimator ─────────────────────────────────────────────────────────────
def compute_caf(signal, alpha, fs, max_lag):
    """
    Estimate the CAF magnitude at a given cycle frequency alpha.

    signal  : input signal (1D array)
    alpha   : cycle frequency to test (Hz)
    fs      : sample rate (Hz)
    max_lag : maximum lag in samples — sweeps from -max_lag to +max_lag
    """
    N    = len(signal)
    n    = np.arange(N)
    # Complex tone at cycle frequency alpha
    # This extracts only the component of the autocorrelation
    # that oscillates at exactly alpha — all other frequencies average to zero
    tone = np.exp(-1j * 2*np.pi * alpha/fs * n)
    lags = np.arange(-max_lag, max_lag + 1)
    caf  = np.zeros(len(lags), dtype=complex)
    for i, tau in enumerate(lags):
        if tau >= 0:
            prod   = signal[tau:] * np.conj(signal[:N-tau]) * tone[tau:]
        else:
            prod   = signal[:N+tau] * np.conj(signal[-tau:]) * tone[:N+tau]
        caf[i] = np.mean(prod)
    return lags, np.abs(caf)

max_lag = 2 * sps

# Test at three key cycle frequencies
alpha_0   = 0       # ordinary autocorrelation
alpha_Rs  = Rs      # baud-rate feature — both BPSK and QPSK respond
alpha_2fc = 2*fc    # conjugate feature — BPSK only

lags, caf_bpsk_0   = compute_caf(bpsk, alpha_0,   fs, max_lag)
lags, caf_bpsk_Rs  = compute_caf(bpsk, alpha_Rs,  fs, max_lag)
lags, caf_bpsk_2fc = compute_caf(bpsk, alpha_2fc, fs, max_lag)

lags, caf_qpsk_0   = compute_caf(qpsk, alpha_0,   fs, max_lag)
lags, caf_qpsk_Rs  = compute_caf(qpsk, alpha_Rs,  fs, max_lag)
lags, caf_qpsk_2fc = compute_caf(qpsk, alpha_2fc, fs, max_lag)

Summary

The CAF is the Fourier series coefficient of the time-varying autocorrelation $R_x(t, \tau)$ with respect to time $t$. At each cycle frequency $\alpha$, it gives one function of lag $\tau$. A non-zero peak in $|R_x^\alpha(\tau)|$ at some $\alpha \neq 0$ confirms that the signal is cyclostationary at that frequency.

For BPSK and QPSK specifically:

• Both are cyclostationary at $\alpha = k \cdot R_s$ — the baud-rate harmonics
• Only BPSK is cyclostationary at $\alpha = 2f_c$ — the conjugate feature

This single difference at $\alpha = 2f_c$ is the second-order discriminator between BPSK and QPSK, and it is directly visible in the CAF without any demodulation or prior knowledge of the signal parameters.

Next in the SIGINT series: the Spectral Correlation Function (SCF) in depth — the full 2D surface, how to read the fingerprint across all $(\alpha, f)$ pairs, and what changes when you add noise.

Jupyter notebook: sigint_02_cyclic_autocorrelation_function.ipynb