In the MATLAB code below, a comparison between frequency-selective fading and flat fading is shown.
In frequency-selective fading, multipath propagation causes multiple delayed copies of the signal to arrive at the receiver. When the channel delay spread exceeds the symbol duration, these delayed components overlap, resulting in inter-symbol interference (ISI).
In flat fading, ISI does not occur because the signal bandwidth is much smaller than the channel’s coherence bandwidth. Therefore, the channel response remains approximately constant across the signal bandwidth, and all symbols experience the same fading.
MATLAB Code for frequency selective fading channel
% OFDM over frequency selective Rayleigh fading channel
clc; clearvars; close all;
% Simulation parameters
nSym = 10^4; % Number of OFDM symbols
EbN0dB = 0:2:20; % Eb/N0 range
MOD_TYPE = 'MPSK'; % 'MPSK' or 'MQAM'
M = 4; % QPSK
N = 64; % Total number of subcarriers
Ncp = 16; % Cyclic prefix length
L = 10; % Number of channel taps
k = log2(M); % Bits per symbol
EsN0dB = 10*log10(k*(N/(N+Ncp))) + EbN0dB; % account CP loss
errors = zeros(1,length(EsN0dB));
% Monte Carlo simulation
for i = 1:length(EsN0dB)
for j = 1:nSym
% -------- Transmitter --------
d = ceil(M*rand(1,N)); % Random symbols
[X , ref] = modulation_mapper(MOD_TYPE,M,d);
x = ifft(X,N); % IDFT
s = add_cyclic_prefix(x,Ncp); % Add CP
% -------- Channel (Freq-selective Rayleigh) --------
h = 1/sqrt(2)*(randn(1,L)+1i*randn(1,L)); % CIR
H = fft(h,N); % Freq response
hs = conv(h,s); % Channel filtering
r = awgn(hs,EsN0dB(i)); % Add noise
% -------- Receiver --------
y = remove_cyclic_prefix(r,Ncp,N); % Remove CP
Y = fft(y,N); % DFT
V = Y ./ H; % Equalization
[~, dcap] = iqOptDetector(V,ref);
errors(i) = errors(i) + sum(d ~= dcap);
end
end
% SER results
simulatedSER = errors/(nSym*N);
theoreticalSER = ser_rayleigh(EbN0dB,MOD_TYPE,M);
% Plot
semilogy(EbN0dB,simulatedSER,'ko'); hold on;
title(['Performance of ',num2str(M),'-',MOD_TYPE,...
' OFDM over Freq Selective Rayleigh channel']);
xlabel('Eb/N0 (dB)');
ylabel('Symbol Error Rate');
function [X, ref] = modulation_mapper(MOD_TYPE, M, d)
% Modulation mapper for OFDM transmitter
% MOD_TYPE - 'MPSK' or 'MQAM'
% M - Modulation order (BPSK=2, QPSK=4, 16-QAM=16, etc.)
% d - Data symbols drawn from {1,2,...,M}
%
% Outputs:
% X - Modulated complex symbols
% ref - Ideal constellation points (for IQ detector)
if strcmpi(MOD_TYPE, 'MPSK')
[X, ref] = mpsk_modulator(M, d); % M-PSK modulation
elseif strcmpi(MOD_TYPE, 'MQAM')
[X, ref] = mqam_modulator(M, d); % M-QAM modulation
else
error('Invalid modulation type. Use ''MPSK'' or ''MQAM''.');
end
end
function [s,ref]=mpsk_modulator(M,d)
%Function to MPSK modulate the vector of data symbols - d
%[s,ref]=mpsk_modulator(M,d) modulates the symbols defined by the
%vector d using MPSK modulation, where M specifies the order of
%M-PSK modulation and the vector d contains symbols whose values
%in the range 1:M. The output s is the modulated output and ref
%represents the reference constellation that can be used in demod
ref_i= 1/sqrt(2)*cos(((1:1:M)-1)/M*2*pi);
ref_q= 1/sqrt(2)*sin(((1:1:M)-1)/M*2*pi);
ref = ref_i+1i*ref_q;
s = ref(d); %M-PSK Mapping
end
function s = add_cyclic_prefix(x,Ncp)
%function to add cyclic prefix to the generated OFDM symbol x that
%is generated at the output of the IDFT block
% x - ofdm symbol without CP (output of IDFT block)
% Ncp-num. of samples at x's end that will copied to its beginning
% s - returns the cyclic prefixed OFDM symbol
s = [x(end-Ncp+1:end) x]; %Cyclic prefixed OFDM symbol
end
function y = remove_cyclic_prefix(r,Ncp,N)
%function to remove cyclic prefix from the received OFDM symbol r
% r - received ofdm symbol with CP
% Ncp - num. of samples at beginning of r that need to be removed
% N - number of samples in a single OFDM symbol
% y - returns the OFDM symbol without cyclic prefix
y=r(Ncp+1:N+Ncp);%cut from index Ncp+1 to N+Ncp
end
function [idealPoints,indices]= iqOptDetector(received,ref)
%Optimum Detector for 2-dim. signals (MQAM,MPSK,MPAM) in IQ Plane
%received - vector of form I+jQ
%ref - reference constellation of form I+jQ
%Note: MPAM/BPSK are one dim. modulations. The same function can be
%applied for these modulations since quadrature is zero (Q=0).
x=[real(received); imag(received)]';%received vec. in cartesian form
y=[real(ref); imag(ref)]';%reference vec. in cartesian form
[idealPoints,indices]= minEuclideanDistance(x,y);
end
function [idealPoints,indices]= minEuclideanDistance(x,y)
%function to compute the pairwise minimum Distance between two
%vectors x and y in p-dimensional signal space and select the
%vectors in y that provides the minimum distances.
% x - a matrix of size mxp
% y - a matrix of size nxp. This acts as a reference against
% which each point in x is compared.
% idealPoints - contain the decoded vector
% indices - indices of the ideal points in reference matrix y
[m,p1] = size(x);[n,p2] = size(y);
if p1~=p2
error('Dimension Mismatch: x and y must have same dimension')
end
X = sum(x.*x,2);
Y = sum(y.*y,2)';
d = X(:,ones(1,n)) + Y(ones(1,m),:) - 2*x*y';%Squared Euclidean Dist.
[~,indices]=min(d,[],2); %Find the minimum value along DIM=2
idealPoints=y(indices,:);
indices=indices.';
end
function [ser] = ser_rayleigh(EbN0dB,MOD_TYPE,M)
%Compute Theoretical Symbol Error rates for MPSK or MQAM modulations
%EbN0dB - list of SNR per bit points
%MOD_TYPE - 'MPSK' or 'MQAM'
%M - Modulation level for the chosen modulation
% - For MPSK M can be any power of 2
% - For MQAM M must be even power of 2 (square QAM only)
gamma_b = 10.^(EbN0dB/10); %SNR per bit in linear scale
gamma_s = log2(M)*gamma_b; %SNR per symbol in linear scale
switch lower(MOD_TYPE)
case {'bpsk'}
ser = 0.5*(1-sqrt(gamma_b/(1+gamma_b)));
case {'mpsk','psk'}
ser = zeros(size(gamma_s));
for i=1:length(gamma_s), %for each SNR point
g = sin(pi/M).^2;
fun = @(x) 1./(1+(g.*gamma_s(i)./(sin(x).^2))); %MGF
ser(i) = (1/pi)*integral(fun,0,pi*(M-1)/M);
end
case {'mqam','qam'}
ser = zeros(size(gamma_s));
for i=1:length(gamma_s) %for each SNR point
g = 1.5/(M-1);
fun = @(x) 1./(1+(g.*gamma_s(i)./(sin(x).^2)));%MGF
ser(i) = 4/pi*(1-1/sqrt(M))*integral(fun,0,pi/2)-4/pi*(1-1/sqrt(M))^2* ...
integral(fun,0,pi/4);
end
case {'mpam','pam'}
ser = zeros(size(gamma_s));
for i=1:length(gamma_s) %for each SNR point
g = 3/(M^2-1);
fun = @(x) 1./(1+(g.*gamma_s(i)./(sin(x).^2)));%MGF
ser(i) = 2*(M-1)/(M*pi)*integral(fun,0,pi/2);
end
end
end
Output
MATLAB Code for flat fading channel
% OFDM over flat Rayleigh fading channel
clc; clearvars; close all;
% Simulation parameters
nSym = 10^4; % Number of OFDM symbols
EbN0dB = 0:2:20; % Eb/N0 range
MOD_TYPE = 'MPSK'; % 'MPSK' or 'MQAM'
M = 4; % QPSK
N = 64; % Total number of subcarriers
Ncp = 16; % Cyclic prefix length
k = log2(M); % Bits per symbol
EsN0dB = 10*log10(k*(N/(N+Ncp))) + EbN0dB; % account CP loss
errors = zeros(1,length(EsN0dB));
% Monte Carlo simulation
for i = 1:length(EsN0dB)
for j = 1:nSym
% -------- Transmitter --------
d = ceil(M*rand(1,N)); % Random symbols
[X , ref] = modulation_mapper(MOD_TYPE,M,d);
x = ifft(X,N); % IDFT
s = add_cyclic_prefix(x,Ncp); % Add CP
% -------- Channel (Flat Rayleigh) --------
h = 1/sqrt(2)*(randn + 1i*randn); % Single-tap flat fading
H = h; % Frequency response
r = h*s; % Flat-fading multiplication
r = awgn(r,EsN0dB(i),'measured'); % Add AWGN noise
% -------- Receiver --------
y = remove_cyclic_prefix(r,Ncp,N); % Remove CP
Y = fft(y,N); % DFT
V = Y ./ H; % Equalization
[~, dcap] = iqOptDetector(V,ref);
errors(i) = errors(i) + sum(d ~= dcap);
end
end
% SER results
simulatedSER = errors/(nSym*N);
theoreticalSER = ser_rayleigh(EbN0dB,MOD_TYPE,M);
% Plot
semilogy(EbN0dB,simulatedSER,'ko'); hold on;
title(['Performance of ',num2str(M),'-',MOD_TYPE,...
' OFDM over Flat Rayleigh channel']);
xlabel('Eb/N0 (dB)');
ylabel('Symbol Error Rate');
%% ------------------------ Functions ------------------------
function [X, ref] = modulation_mapper(MOD_TYPE, M, d)
if strcmpi(MOD_TYPE, 'MPSK')
[X, ref] = mpsk_modulator(M, d);
elseif strcmpi(MOD_TYPE, 'MQAM')
[X, ref] = mqam_modulator(M, d);
else
error('Invalid modulation type. Use ''MPSK'' or ''MQAM''.');
end
end
function [s, ref] = mpsk_modulator(M,d)
ref_i= 1/sqrt(2)*cos(((1:M)-1)/M*2*pi);
ref_q= 1/sqrt(2)*sin(((1:M)-1)/M*2*pi);
ref = ref_i + 1i*ref_q;
s = ref(d);
end
function s = add_cyclic_prefix(x,Ncp)
s = [x(end-Ncp+1:end) x];
end
function y = remove_cyclic_prefix(r,Ncp,N)
y = r(Ncp+1:N+Ncp);
end
function [idealPoints,indices] = iqOptDetector(received,ref)
x = [real(received); imag(received)]';
y = [real(ref); imag(ref)]';
[idealPoints,indices] = minEuclideanDistance(x,y);
end
function [idealPoints,indices] = minEuclideanDistance(x,y)
[m,p1] = size(x); [n,p2] = size(y);
if p1 ~= p2, error('Dimension mismatch'); end
X = sum(x.*x,2);
Y = sum(y.*y,2)';
d = X(:,ones(1,n)) + Y(ones(1,m),:) - 2*x*y';
[~,indices] = min(d,[],2);
idealPoints = y(indices,:);
indices = indices.';
end
function ser = ser_rayleigh(EbN0dB,MOD_TYPE,M)
gamma_b = 10.^(EbN0dB/10); % SNR per bit linear
gamma_s = log2(M)*gamma_b; % SNR per symbol linear
switch lower(MOD_TYPE)
case {'bpsk'}
ser = 0.5*(1-sqrt(gamma_b./(1+gamma_b)));
case {'mpsk','psk'}
ser = zeros(size(gamma_s));
for i=1:length(gamma_s)
g = sin(pi/M)^2;
fun = @(x) 1./(1 + (g*gamma_s(i)./(sin(x).^2)));
ser(i) = (1/pi)*integral(fun,0,pi*(M-1)/M);
end
case {'mqam','qam'}
ser = zeros(size(gamma_s));
for i=1:length(gamma_s)
g = 1.5/(M-1);
fun = @(x) 1./(1 + (g*gamma_s(i)./(sin(x).^2)));
ser(i) = 4/pi*(1-1/sqrt(M))*integral(fun,0,pi/2) ...
- 4/pi*(1-1/sqrt(M))^2*integral(fun,0,pi/4);
end
end
end
Output
