Cyclic Prefix in OFDM: hands-on demo in Matlab

Synopsis: Cyclic prefix in OFDM, tricks a natural channel to perform circular convolution. This simplifies equalizer design at the receiver. Hands-on demo in Matlab.

Cyclic Prefix-ed OFDM

A cyclic-prefixed OFDM (CP-OFDM) transceiver architecture is typically implemented using inverse discrete Fourier transform (IDFT) and discrete Fourier transform (DFT) blocks (refer Figure 13.3). In an OFDM transmitter, the modulated symbols are assigned to individual subcarriers and sent to an IDFT block. The output of the IDFT block, generally viewed as time domain samples, results in an OFDM symbol. Such OFDM symbols are then transmitted across a channel with certain channel impulse response (CIR). On the other hand, the receiver applies DFT over the received OFDM symbols for further demodulation of the information in the individual subcarriers.

This article is part of the book
Wireless Communication Systems in Matlab (second edition), ISBN: 979-8648350779 available in ebook (PDF) format and Paperback (hardcopy) format.

In reality, a multipath channel or any channel in nature acts as a linear filter on the transmitted OFDM symbols. Mathematically, the transmitted OFDM symbol (denoted as gets linearly convolved with the CIR and gets corrupted with additive white gaussian noise – designated as . Denoting linear convolution as ‘‘, the received signal in discrete-time can be represented as

The idea behind using OFDM is to combat frequency selective fading, where different frequency components of a transmitted signal can undergo different levels of fading. The OFDM divides the available spectrum in to small chunks called sub-channels. Within the subchannels, the fading experienced by the individual modulated symbol can be considered flat. This opens-up the possibility of using a simple frequency domain equalizer to neutralize the channel effects in the individual subchannels.

Circular convolution and designing a simple frequency domain equalizer

From one of the DFT properties, we know that the circular convolution of two sequences in time domain is equivalent to the multiplication of their individual responses in frequency domain.

Let and are two sequences of length with their DFTs denoted as and respectively. Denoting circular convolution as ,

If we ignore channel noise in the OFDM transmission, the received signal is written as

We can note that the channel performs a linear convolution operation on the transmitted signal. Instead, if the channel performs circular convolution (which is not the case in nature) then the equation would have been

By applying the DFT property given in equation 2,

As a consequence, the channel effects can be neutralized at the receiver using a simple frequency domain equalizer (actually this is a zero-forcing equalizer when viewed in time domain) that just inverts the estimated channel response and multiplies it with the frequency response of the received signal to obtain the estimates of the OFDM symbols in the subcarriers as

Demonstrating the role of Cyclic Prefix

The simple frequency domain equalizer shown in equation 6 is possible only if the channel performs circular convolution. But in nature, all channels perform linear convolution. The linear convolution can be converted into circular convolution by adding Cyclic Prefix (CP) in the OFDM architecture. The addition of CP makes the linear convolution imparted by the channel appear as circular convolution to the DFT process at the receiver (Reference [1]).

Let’s understand this by demonstration.

To simply stuffs, we will create two randomized vectors for and . is of length , the channel impulse response is of length and we intend to use -point DFT/IDFT wherever applicable.

N=8; %period of DFT
s = randn(1,8);
h = randn(1,3);
>> s =   -0.0155    2.5770    1.9238   -0.0629   -0.8105    0.6727   -1.5924   -0.8007
>> h =   -0.4878   -1.5351    0.2355

Now, convolve the vectors and linearly and circularly. The outputs are plotted in Figure 1. We note that the linear convolution and the circular convolution does not yield identical results.

lin_s_h = conv(h,s) %linear convolution of h and s
cir_s_h = cconv(h,s,N) % circular convolution of h and s with period N
lin_s_h = 0.0076 -1.2332 -4.8981 -2.3158  0.9449  0.9013 -0.4468  2.9934  0.8542 -0.1885 
cir_s_h = 0.8618 -1.4217 -4.8981 -2.3158  0.9449  0.9013 -0.4468  2.9934
Difference between linear convolution and circular convolution
Figure 1: Difference between linear convolution and circular convolution

Let’s append a cyclic prefix to the signal by copying last symbols from and pasting it to its front. Since the channel delay is 3 symbols (CIR of length 3), we need to add atleast 2 CP symbols.

Ncp = 2; %number of symbols to copy and paste for CP
s_cp = [s(end-Ncp+1:end) s]; %copy last Ncp symbols from s and prefix it.
s_cp = -1.5924 -0.8007 -0.0155  2.5770  1.9238 -0.0629 -0.8105  0.6727 -1.5924 -0.8007

Lets assume that we send the cyclic-prefixed OFDM symbol through the channel which perform linear filtering

lin_scp_h = conv(h,s_cp) %linear convolution of CP-OFDM symbol s_cp and the CIR h
lin_scp_h = 0.7767  2.8350  0.8618 -1.4217 -4.8981 -2.3158  0.9449  0.9013 -0.4468  2.9934  0.8542 -0.1885

Compare the outputs due to cir_s_h and lin_scp_h . We can immediately recognize that that the middle section of lin_scp_h is exactly same as the cir_s_h vector. This is shown in the figure 2.

cir_s_h   = 0.8618  -1.4217 -4.8981 -2.3158  0.9449  0.9013 -0.4468  2.9934
lin_scp_h = 0.7767   2.8350  0.8618 -1.4217 -4.8981 -2.3158  0.9449  0.9013 -0.4468 2.9934 0.8542 -0.1885

We have just tricked the channel to perform circular convolution by adding a cyclic extension to the OFDM symbol. At the receiver, we are only interested in the middle section of the lin_scp_h which is the channel output. The first block in the OFDM receiver removes the additional symbols from the front and back of the channel output. The resultant vector is the received symbol r after the removal of cyclic prefix in the receiver.

r = lin_scp_h(Ncp+1:N+Ncp)%cut from index Ncp+1 to N+Ncp

Verifying DFT property

The DFT property given in equation 5 can be re-written as

To verify this in Matlab, take N-point DFTs of the received signal and CIR. Then, we can see that the IDFT of product of DFTs of and will be equal to the N-point circular convolution of and

R = fft(r,N); %frequency response of received signal
H = fft(h,N); %frequency response of CIR
S = fft(s,N); %frequency response of OFDM signal (non CP)

r1 = ifft(S.*H); %IFFT of product of individual DFTs

display(['IFFT(DFT(H)*DFT(S)) : ',num2str(r1)])
display([cconv(s,h): ', numstr(r)])
IFFT(DFT(H)*DFT(S)) : 0.86175  -1.4217  -4.8981  -2.3158  0.94491  0.90128 -0.44682  2.9934
cconv(s,h):           0.86175  -1.4217  -4.8981  -2.3158  0.94491  0.90128 -0.44682  2.9934

Rate this article: Note: There is a rating embedded within this post, please visit this post to rate it.

Reference:

[1] H. Sari, G. Karam, and I. Jeanclaude, Transmission Techniques for Digital Terrestrial TV Broadcasting, IEEE Commun. Magazine, Vol. 33, pp. 100-109, Feb. 1995.↗

Topics in this chapter

Orthogonal Frequency Division Multiplexing (OFDM)
● Introduction
Understanding the role of cyclic prefix in a CP-OFDM system
 □ Circular convolution and designing a simple frequency domain equalizer
 □ Demonstrating the role of cyclic prefix
 □ Verifying DFT property
Discrete-time implementation of baseband CP-OFDM
Performance of MPSK-CP-OFDM and MQAM-CP-OFDM on AWGN channel
● Performance of MPSK-CP-OFDM and MQAM-CP-OFDM on frequency selective Rayleigh channel

Books by the author


Wireless Communication Systems in Matlab
Second Edition(PDF)

Note: There is a rating embedded within this post, please visit this post to rate it.
Checkout Added to cart

Digital Modulations using Python
(PDF ebook)

Note: There is a rating embedded within this post, please visit this post to rate it.
Checkout Added to cart

Digital Modulations using Matlab
(PDF ebook)

Note: There is a rating embedded within this post, please visit this post to rate it.
Checkout Added to cart
Hand-picked Best books on Communication Engineering
Best books on Signal Processing

Methods to compute linear convolution

Mathematical details of convolution, its relationship to polynomial multiplication and the application of Toeplitz matrices in computing linear convolution are discussed in the previous article. A short survey of different techniques to compute discrete linear convolution (with Matlab code) is given here.

Definition

Given an LTI (Linear Time Invariant) system with impulse response \(h[n]\) and an input sequence \(x[n]\), the output of the system \(y[n]\) is obtained by convolving the input sequence and impulse response.

\[ y[k]=h[n] \ast x[n] = \sum_{i= -\infty}^{\infty} x[i] h[k-i] \]

where, the sequence \(x[n]\) is of length \(N\) and \(h[n]\) is of length \(M\).

This article is part of the following books
Digital Modulations using Matlab : Build Simulation Models from Scratch, ISBN: 978-1521493885
Digital Modulations using Python ISBN: 978-1712321638
Wireless communication systems in Matlab ISBN: 979-8648350779
All books available in ebook (PDF) and Paperback formats

Method 1 – Brute-Force Method:

The above equation can simply be implemented using nested for-loops, which consume the most computational time of all the methods given here. Typically, the computational complexity is \(O(n^2)\) time.

Python

def conv_brute_force(x,h):
	"""
	Brute force method to compute convolution
	Parameters:
	x, h : numpy vectors
	Returns:
		y : convolution of x and h
	"""
	N=len(x)
	M=len(h)
	y = np.zeros(N+M-1) #array filled with zeros
	
	for i in np.arange(0,N):
		for j in np.arange(0,M):
			y[i+j] = y[i+j] + x[i] * h[j]
return y

Matlab

for i = 0:n+m,
   yi = 0;
for i = 0 :n,
   for k = 0 to m,
      y[i+k] = y[i+k] + x[i] · h[k];

Method 2 – Using Toeplitz Matrix:

When the sequences \(h[n]\) and \(x[n]\) are represented as matrices, the linear convolution operation can be equivalently represented as

\[y=h \ast x=x \ast h=Toeplitz(h) X\]

Assume that the sequence \(h[n]\) is of length 4 given by \(h[n]=[h_0,h_1,h_2,h_3]\) and the sequence \(x[n]\) is of length 3 given by \(x[n]=[x_0,x_1,x_2]\). The convolution \(h[n] \ast x[n]\) is given by

\[y[k]=h[n] \ast x[n] = \displaystyle{\sum_{i=-\infty}^{\infty}x[i]h[k-i]}, \quad\quad k=0,1,\cdots,5\]

Equivalent representation of the above convolution can be written as

\[\begin{bmatrix}y[0]\\y[1]\\y[2]\\y[3]\\y[4]\\y[5]\end{bmatrix}=\begin{bmatrix}h[0] &0 & 0 \\h[1] &h[0] & 0 \\h[2] & h[1] & h[0] \\h[3] & h[2] & h[1] \\0 & h[3] & h[2] \\0 & 0 & h[3]\end{bmatrix}\begin{bmatrix}x[0]\\x[1]\\x[2] \end{bmatrix}\]

The matrix representing the incremental delays of \(h[n]\) used in the above equation is a special form of matrix called Toeplitz matrix. Toeplitz matrix have constant entries along their diagonals. Toeplitz matrices are used to model systems that posses shift invariant properties. The property of shift invariance is evident from the matrix structure itself. Since we are modelling a Linear Time Invariant system[1], Toeplitz matrices are our natural choice. On a side note, a special form of Toeplitz matrix called “circulant matrix” is used in applications involving circular convolution and Discrete Fourier Transform (DFT)[2].

For python code: refer the book – Digital modulations using Python

Matlab has inbuilt function to compute Toeplitz matrix from given vector. Toeplitz matrix of sequence \(h\) is given as

\[Toeplitz (h) = \begin{bmatrix}h[0] &0 & 0 \\h[1] &h[0] & 0 \\h[2] & h[1] & h[0] \\h[3] & h[2] & h[1] \\0 & h[3] & h[2] \\0 & 0 & h[3]\end{bmatrix}\]

In Matlab, it is constructed by calling Toeplitz function[3] with two arguments. The first argument is the first column of the above matrix -\( [h_0,h_1,h_2,h_3,0,0]\) and the second argument is the first row of the above matrix – \([h_0,0,0]\).

y=toeplitz([h0 h1 h2 h3 0 0],[h0 0 0])*x.';

In the generalized case, where \(x\) and \(h\) are of any arbitrary lengths (say \(N\) and \(M\)), the code can be modified as follows [here it is assumed that \( length (x) \geq length(h)\) ]

k=toeplitz([h zeros(1,length(x)-1)],[h(1) zeros(1,length(x)-1])*x.'

Method 3: Using FFT:

Computation of convolution using FFT (Fast Fourier Transform) has the advantage of reduced computational complexity when the length of inputs are large. To compute convolution, take FFT of the two sequences \(x\) and \(h\) with FFT length set to convolution output length \(length (x)+length(h)-1\), multiply the results and convert back to time-domain using IFFT (Inverse Fast Fourier Transform). Note that FFT is a direct implementation of circular convolution in time domain. Here we are attempting to compute linear convolution using circular convolution (or FFT) with zero-padding either one of the input sequence. This causes inefficiency when compared to circular convolution. Nevertheless, this method still provides \(O\left(\frac{N}{log_2N}\right)\) savings over brute-force method.

\[y(n) = IFFT[FFT_L (X) \ast FFT_L (H) ], \quad \; L=2^{nextpower2(N+M-1)}\]

Usually, the following algorithm is suffice which ignores additional zeros in the output terms.

\[y(n) = IFFT [ FFT_L (X) \ast FFT_L (H)], \quad \; L=N+M-1\]

Python

from scipy.fftpack import fft,ifft
y=ifft(fft(x,L)*(fft(h,L))) #Convolution using FFT and IFFT

Matlab

y=ifft(fft(x,L).*(fft(h,L)))

Other Methods:

If the input sequence is of infinite length or very very large as in many real time applications, block processing methods like Overlap-Add and Overlap-Save methods can be used to compute convolution in a faster and efficient way.

Standard Algorithms for Fast Convolution:

Standard algorithms for computing convolution that greatly reduce the computational complexity are

1) Cook-Toom Algorithm
2) Modified Cook-Toom Algorithm
3) Winograd Algorithm
4) Modified Winograd Algorithm
5) Iterated Convolution

Refer [4] -Fast Algorithms for Signal Processing by Richard E. Blahut for more details on the above algorithms (see reference section below).

Matlab Code

Implementing method 2 (convolution using Toeplitz matrix transformation) and method 3 (convolution using FFT) and comparing against Matlab’s standard conv function:

%Create random vectors for test
x=rand(1,7)+1i*rand(1,7);
h=rand(1,3)+1i*rand(1,3);
L=length(x)+length(h)-1;

%Convolution Using Toeplitz matrix
y1=toeplitz([h zeros(1,length(x)-1)],[h(1) zeros(1,length(x)-1)])*x.'
%Convolution FFT
y2=ifft(fft(x,L).*(fft(h,L))).'
%Matlab's standard function
y3=conv(h,x).'

Simulation output:

Comparing method 2 and method 3 against Matlab’s standard convolution function returns identical results

\[y1 =\begin{bmatrix} 0.2428 + 0.4536i\\ 0.1512 + 0.6668i\\ 0.8663 + 1.1188i\\ 0.5332 + 1.0162i\\ 0.8589 + 0.9822i\\ 0.4633 + 1.1006i\\ 0.4024 + 1.4706i\\ 0.3225 + 0.9790i\\ 0.2324 + 0.8227i\end{bmatrix}; y2 = \begin{bmatrix}0.2428 + 0.4536i\\ 0.1512 + 0.6668i\\ 0.8663 + 1.1188i\\ 0.5332 + 1.0162i\\ 0.8589 + 0.9822i\\ 0.4633 + 1.1006i\\ 0.4024 + 1.4706i\\ 0.3225 + 0.9790i\\ 0.2324 + 0.8227i\\ \end{bmatrix};y3 = \begin{bmatrix}0.2428 + 0.4536i\\ 0.1512 + 0.6668i\\ 0.8663 + 1.1188i\\ 0.5332 + 1.0162i\\ 0.8589 + 0.9822i\\ 0.4633 + 1.1006i\\ 0.4024 + 1.4706i\\ 0.3225 + 0.9790i\\ 0.2324 + 0.8227i\end{bmatrix}\]

References:

[1] Reddi.S.S,”Eigen Vector properties of Toeplitz matrices and their application to spectral analysis of time series, Signal Processing, Vol 7, North-Holland, 1984, pp 46-56.↗
[2] Robert M. Gray,”Toeplitz and circulant matrices – an overview”,Deptartment of Electrical Engineering,Stanford University,Stanford 94305,USA.↗
[3] Matlab documentation help on Toeplitz command.↗
[4] Richard E. Blahut,”Fast Algorithms for Signal Processing”,ISBN-13: 978-0521190497,Cambridge University Press,August 16, 2010.↗

Topics in this chapter

Essentials of Signal Processing
● Generating standard test signals
 □ Sinusoidal signals
 □ Square wave
 □ Rectangular pulse
 □ Gaussian pulse
 □ Chirp signal
Interpreting FFT results - complex DFT, frequency bins and FFTShift
 □ Real and complex DFT
 □ Fast Fourier Transform (FFT)
 □ Interpreting the FFT results
 □ FFTShift
 □ IFFTShift
Obtaining magnitude and phase information from FFT
 □ Discrete-time domain representation
 □ Representing the signal in frequency domain using FFT
 □ Reconstructing the time domain signal from the frequency domain samples
● Power spectral density
Power and energy of a signal
 □ Energy of a signal
 □ Power of a signal
 □ Classification of signals
 □ Computation of power of a signal - simulation and verification
Polynomials, convolution and Toeplitz matrices
 □ Polynomial functions
 □ Representing single variable polynomial functions
 □ Multiplication of polynomials and linear convolution
 □ Toeplitz matrix and convolution
Methods to compute convolution
 □ Method 1: Brute-force method
 □ Method 2: Using Toeplitz matrix
 □ Method 3: Using FFT to compute convolution
 □ Miscellaneous methods
Analytic signal and its applications
 □ Analytic signal and Fourier transform
 □ Extracting instantaneous amplitude, phase, frequency
 □ Phase demodulation using Hilbert transform
Choosing a filter : FIR or IIR : understanding the design perspective
 □ Design specification
 □ General considerations in design

Books by the author


Wireless Communication Systems in Matlab
Second Edition(PDF)

Note: There is a rating embedded within this post, please visit this post to rate it.
Checkout Added to cart

Digital Modulations using Python
(PDF ebook)

Note: There is a rating embedded within this post, please visit this post to rate it.
Checkout Added to cart

Digital Modulations using Matlab
(PDF ebook)

Note: There is a rating embedded within this post, please visit this post to rate it.
Checkout Added to cart
Hand-picked Best books on Communication Engineering
Best books on Signal Processing