I found a few issues in the code that caused a mismatch between the theoretical and simulated Bit Error Rate values.
Here are the issues and their corresponding fixes:
1. Directly used binary bits (data) with "pskmod", but "pskmod" requires integer symbols.
- Perform bit-to-symbol conversion using "bi2de" to correctly map input bits to QPSK symbols.
2. Used 1/2 instead of the correct factor 
for Rayleigh fading, resulting in incorrect channel power.
for Rayleigh fading, resulting in incorrect channel power.- Normalize the channel coefficient using
to ensure an average power of 1.
3. The code only used "abs(h)", ignoring phase shifts, which QPSK is sensitive to. Additionally, no equalization was performed, leading to channel distortion.
- Use the full complex channel coefficient h for both magnitude and phase effects, and apply equalization by dividing the received signal by h.
4. Theoretical BER Calculation:
- For QPSK in AWGN, use
. - For Rayleigh fading, use
5. Errors were averaged per iteration using only Nb, which was incorrect.
- Accumulate total bit errors across all iterations and divide by the total number of transmitted bits (Nb * N_iterations) for accurate BER estimation.
Here is the full corrected code:
%% Clean Up
clc;
clear;
close all;
%% Parameters Definition
Eb_No = 0:20; % [dB]
Eb_No_linear = db2pow(Eb_No); % Linear scale
N = 64; % number of symbols
N_iterations = 1000;
M = 4; % QPSK
Im = log2(M); % Bits per symbol
Nb = N * Im; % Total number of bits per iteration
% Initialize error arrays
error_awgn = zeros(length(Eb_No), N_iterations);
error_rayleigh_awgn = zeros(length(Eb_No), N_iterations);
BER_awgn = zeros(size(Eb_No));
BER_rayleigh_awgn = zeros(size(Eb_No));
%% Calculating The Theoretical Results
BER_theoretical_awgn = qfunc(sqrt(2 * Eb_No_linear)); % Corrected QPSK BER for AWGN
BER_theoretical_fading = 0.5 * (1 - sqrt(Eb_No_linear ./ (1 + Eb_No_linear))); % Corrected BER for Rayleigh fading
%% Calculating The Simulation Results
for n = 1:length(Eb_No)
Es_No = Eb_No(n) + 10*log10(Im); % Es/No in dB
total_bit_errors_awgn = 0;
total_bit_errors_rayleigh = 0;
for k = 1:N_iterations
% Generate random binary data
data = randi([0 1], Nb, 1);
% Map bits to symbols (0, 1, 2, 3 for QPSK)
data_symbols = bi2de(reshape(data, Im, []).', 'left-msb');
% QPSK Modulation with Gray coding
x = pskmod(data_symbols, M, pi/M, 'gray');
% Rayleigh Fading Channel Coefficient
h = 1/sqrt(2) * (randn(N, 1) + 1i * randn(N, 1)); % Normalized
% AWGN Channel
y = awgn(x, Es_No, "measured");
% Rayleigh + AWGN Channel
xf = h .* x; % Apply full Rayleigh fading (magnitude and phase)
yf = awgn(xf, Es_No, "measured");
% Equalize Rayleigh channel effects
yf_equalized = yf ./ h; % Equalize by dividing with channel coefficient
% Demodulation
y_est = pskdemod(y, M, pi/M, 'gray');
yf_est = pskdemod(yf_equalized, M, pi/M, 'gray');
% Convert received symbols back to bits
y_est_bits = de2bi(y_est, Im, 'left-msb').';
yf_est_bits = de2bi(yf_est, Im, 'left-msb').';
% Reshape bits for comparison
y_est_bits = y_est_bits(:);
yf_est_bits = yf_est_bits(:);
% Count bit errors
total_bit_errors_awgn = total_bit_errors_awgn + nnz(data ~= y_est_bits);
total_bit_errors_rayleigh = total_bit_errors_rayleigh + nnz(data ~= yf_est_bits);
end
% Calculate BER for current Eb/No
BER_awgn(n) = total_bit_errors_awgn / (Nb * N_iterations);
BER_rayleigh_awgn(n) = total_bit_errors_rayleigh / (Nb * N_iterations);
end
%% Plotting The Theoretical and Simulation Results
figure(1);
semilogy(Eb_No, BER_theoretical_awgn, 'r-', 'LineWidth', 1.5); hold on;
semilogy(Eb_No, BER_awgn, 'r-o', 'LineWidth', 1);
semilogy(Eb_No, BER_theoretical_fading, 'b-', 'LineWidth', 1.5);
semilogy(Eb_No, BER_rayleigh_awgn, 'b-o', 'LineWidth', 1);
grid on;
ylim([10e-5 0.5]);
hold off;
ylabel('BER');
xlabel('E_b/N_o');
title('QPSK BER curves for AWGN and Rayleigh');
legend('AWGN_{theory}', 'AWGN_{simulation}', 'Rayleigh_{theory}', 'Rayleigh_{simulation}', 'location', 'best');
Refer to the following documentation for more details about the function "pskomd": https://www.mathworks.com/help/comm/ref/pskmod.html
Hope this solves your query.
