The cyclic autocorrelation for rectangular-pulse BPSK can be derived as a relatively simple closed-form expression (see My Papers [6] for example or The Literature [R1]). It can be estimated in a variety of ways, which we will discuss in future posts. The non-conjugate cycle frequencies for the signal are harmonics of the bit rate, , and the conjugate cycle frequencies are the non-conjugate cycle frequencies offset by the doubled carrier, or .

Recall that the simulated rectangular-pulse BPSK signal has samples per bit, or a bit rate of , and a carrier offset of , all in normalized units (meaning the sampling rate is unity). We’ve previously selected a sampling rate of MHz to provide a little physical realism. This means the bit rate is kHz and the carrier offset frequency is kHz. From these numbers, we see that the non-conjugate cycle frequencies are kHz, and that the conjugate cycle frequencies are kHz, or $100 + k 100$ kHz.

The blindly estimated CAF for our noisy rectangular-pulse BPSK signal is shown here:

The non-conjugate CAF for is the conventional autocorrelation function. Here we can see it is triangular in shape, with an additional inflation of the peak for , due to the presence of noise, which has non-zero power. For comparison, here is the numerically evaluated theoretical formula for the CAF for rectangular-pulse BPSK:

The match is excellent.

The theoretical formulas for the cyclic autocorrelation and spectral correlation function for BPSK signals (and other digital QAM/PSK signals) can be found in several places. One of the first places it was published is in the book Statistical Spectral Analysis by W. A. Gardner (The Literature [R1], Chapter 12). The formulas for both the th-order reduced-dimension cyclic cumulant and the th-order cyclic polyspectrum for PAM/PSK/QAM can be found in My Papers [6] and my dissertation. The -th order cyclic cumulant reduces to a lag-shifted version of the cyclic autocorrelation for .

### Interpreting the CAF

The non-conjugate CAF can be interpreted as the correlation between the lag product and the complex-valued sine wave . From that point of view, the CAF is the complex amplitude of the additive sine wave that is present in the lag product. That is, that sine wave has amplitude and phase . Note that the amplitude and phase depend on the delay . When the CAF is zero, there is no finite-strength additive sine-wave component in the lag product for the value of chosen; otherwise there is.

The conjugate CAF can be interpreted as the correlation between the lag product and the complex-valued sine wave .

Some signals may have only non-zero non-conjugate CAF values (such as large-alphabet digital QAM), some only non-zero conjugate CAF values (such as analog amplitude modulation, discounting the ever-present non-conjugate cycle frequency of zero), and some have both (such as our favorite signal, the rectangular-pulse BPSK signal).

### Comparison with a Bandwidth-Efficient BPSK Signal

The rectangular-pulse signal has infinite bandwidth, and therefore it possesses an infinite number of cycle frequencies. However, the strength of the cyclic autocorrelation for most of those cycle frequencies is very small, so that in practical terms, the rectangular-pulse signal possesses ten or so significant features. That’s actually quite a lot compared to practical real-world signals like BPSK with square-root raised-cosine (SRRC) pulses. Here is the estimated CAF for a BPSK with SRRC pulses and a pulse roll-off factor of 0.3 (the excess bandwidth is 30%):

The non-conjugate CAF possesses exactly three cycle frequencies, (only two of those are shown in the plot), and the conjugate CAF also possesses three cycle frequencies . Note also that the width of the CAF is much larger than for the rectangular-pulse signal (the axes are different between the plots for the two signals). This is a direct consequence of the fact that the rectangular-pulse signal is strictly time-limited whereas the SRRC-pulse signal is strictly bandlimited.

hello, could you please explain why cyclic frequencies in your example are equal to k* (0.1) MHz?

LikeLike

I’ve updated the post to answer your question. Leave a comment if unclear!

LikeLike

please I need an example code for cyclostastionary features detection for bpsk modulation

LikeLike

Hi Ben. The CSP Blog is a self-help blog. I don’t give out much code. Let me know if you get stuck, though, trying to write your own. That’s the only way you’ll really understand CSP.

LikeLike

this is my own code but i don’t know if it’s good to calculate the cyclic autoucorrelation function or not;

please help me where is the probléme in this code

thanks dear in advance.

function [caf2, taus, alphas] = caf2(s)

N1=length(s);

N=floor(N1/3);

caf2 = zeros(2*N+1,N);

prods = zeros(2*N+1,N);

taus=-N:N;

alphas=-pi+2*pi/(N):2*pi/(N):pi;

alphas=alphas/pi;

for tau = taus

for alpha=alphas

stat = (s(N+1:2*N)).*exp(-j*pi*alpha*t(N+1:2*N));

shifted = conj(s(N+1-tau:2*N-tau)).*exp(j*pi*alpha*t(N+1-tau:2*N-tau));

prod=stat.*shifted.*exp(-2*j*pi*alpha*t(N+1:2*N));

prods(tau+N+1,:)=prod;

end

caf2(tau+N+1,:)=fftshift(fft(prod));

end

end

LikeLike

Your code does not run to completion:

>> caf2(s)

Undefined function or variable ‘t’.

Error in caf2 (line 11)

stat = (s(N+1:2*N)).*exp(-j*pi*alpha*t(N+1:2*N));

So you need to fix that first.

I see that you are trying to compute a statistic for a large number of cycle frequencies and delays. I would like to suggest the following: Divide and conquer. Do the simplest things first, then build up to the most complicated. For example, try to estimate the cyclic autocorrelation for just a single cycle frequency (say, zero), and a few delays (tau) near zero. Then try to extend to all delays. Then extend so that you can specify any cycle frequency. Finally, when that is all correct, put it in a loop and call it for as many cycle frequencies as you want.

LikeLike

hey

i would like to know about the theorical formula for bpsk cyclic autocorrelation ?

is there any listed reference?

LikeLike

I’ve updated the post with references that should help you find the answer to your question. Let me know if you still have problems or questions.

LikeLike

thanks ser

please comment :

in witch frequency we have peak

Fs = 5000;

t = 0:1/Fs:1-1/Fs;

s = (cos(2*pi*600*t)+randn(size(t)));

N1=length(s);

N=floor(N1/3);

caf2 = zeros(2*N+1,N);

prods = zeros(2*N+1,N);

taus=-0.5:0.5;

alpha=0.2;

stat = (s(N+1:2*N)).*exp(-j*pi*alpha*t(N+1:2*N));

shifted = conj(s(N+1-taus:2*N-taus)).*exp(j*pi*alpha*t(N+1-taus:2*N-taus));

prod=stat.*shifted.*exp(-2*j*pi*alpha*t(N+1:2*N));

caf2=(fft(prod));

thanks in advance

LikeLike

Ben, that code also does not run correctly:

>> caf2_2

Subscript indices must either be real positive integers or logicals.

Error in caf2_2 (line 11)

shifted = conj(s(N+1-taus:2*N-taus)).*exp(j*pi*alpha*t(N+1-taus:2*N-taus));

Is it running without errors on your system? I think you should at least make sure the code runs before sending it to someone else.

LikeLike

thanks mister but my code run without errors

caf2 not caf_2

Fs = 5000;

t = 0:1/Fs:1-1/Fs;

s = (cos(2*pi*600*t)+randn(size(t)));

N1=length(s);

N=floor(N1/3);

caf2 = zeros(2*N+1,N);

prods = zeros(2*N+1,N);

taus=-0.5:0.5;

alpha=0.2;

stat = (s(N+1:2*N)).*exp(-j*pi*alpha*t(N+1:2*N));

shifted = conj(s(N+1-taus:2*N-taus)).*exp(j*pi*alpha*t(N+1-taus:2*N-taus));

prod=stat.*shifted.*exp(-2*j*pi*alpha*t(N+1:2*N));

caf2=fftshift(fft(prod));

plot(abs(caf2))

LikeLike

The code you posted on April 4 produces the following output using MATLAB R2015b:

>> bm_code_040417

Subscript indices must either be real positive integers or logicals.

Error in bm_code_040417 (line 11)

shifted = conj(s(N+1-taus:2*N-taus)).*exp(j*pi*alpha*t(N+1-taus:2*N-taus));

The code you posted on April 7 produces the following output:

>> bm_code_040717

Subscript indices must either be real positive integers or logicals.

Error in bm_code_040717 (line 11)

shifted = conj(s(N+1-taus:2*N-taus)).*exp(j*pi*alpha*t(N+1-taus:2*N-taus));

I’m pretty sure the problem is with

taus=-0.5:0.5;

which yields

>> taus

taus =

-0.500000000000000 0.500000000000000

So you are mixing fractions with integers in your indexing.

This isn’t the kind of debugging I normally help with as it doesn’t have anything to do with understanding CSP.

LikeLike

I use rectangle pulse for BPSK. But the CAF shows followings: 1)there are only 3 lines (represent alpha=-2fc, 0, 2fc)when delay tau=0, and it do nothing with symbol rate; 2)CAF quickly decreases and tends to constant as delay increases when alpha is 2fc. Are those right? Why can’t find cycle frequencies related with symbol rate? Thanks in advance.

LikeLike

Did you use your own rectangular-pulse BPSK signal? What happens if you use the one I provided at the end of the post on rectangular-pulse BPSK? See the end of that post for details.

LikeLike

Thanks for your reply. After using your BPSK signal code, the CAF shows more peaks. But it’s still a big difference from your diagram. I have two questions here:

1) What is the CA formula of BPSK signal? I only find the cyclic spctrum formula in the book Statistical Spectral Analysis; 2) Should s[t] multiply by e^-j*pi*alpha*t before calculating delay product, like code written by Ben Mohamed aboved?

LikeLike

You can look at the references I cite for the formula for the cyclic autocorrelation for PSK/QAM (one of the references

is my doctoral dissertation, which can be downloaded from the CSP Blog here). Or, look at formula (4) in the post on digital QAM and PSK. Specialize that formula to the case of . That formula is slightly more general than the formula for the “symmetrical” cyclic autocorrelation. The symmetrical autocorrelation uses the lags and , whereas the formula (4) uses generic lags .

You don’t have to multiply by before averaging the lag product. You can either multiply the lag product by or multiply by and by and then multiply them together.

I suggest that you first try to get the cyclic autocorrelation correct for a single cycle frequency before you try to write code that computes the entire cyclic autocorrelation surface. I have seen many attempts at the latter that fail due to not understanding how to do the former. By the way, the surface I show in the post includes cyclic autocorrelation estimates for only those cycle frequencies that are known to exist for the rectangular-pulse BPSK signal I used.

LikeLike