Help to convert formula into matlab code - QAM channel capacity

5 views (last 30 days)
Hi!
I'm trying to compute the Channel capacity of a QAM input x with M-Size-alphabet and vector-probability p over a AWGN channel in order to compare channel capacity with respect to input symbol distribution and SNR.
Essentially I'm having trouble in Matlab implementing the formula that allows the capacity calculation, and the implementation of the equations:
where the conditional probability for 2-dimensional constellation over awgn channel is:
The problem is in 2D, so all the integral are double (over dy1 and dy2) and x is an array of complex value.
It doesn't work because I got confused between symbolic expression and function handle and the integral is giving me problems.
x is a vector of M complex value (QAM-symbols) while p is the probability of each of these symbols.
SNR_dB = 10;
QAM_sym16 = qammod(alphabet16, 16);
QAM_normSym16 = QAM_sym16 / sqrt(1/16*norm(QAM_sym16, 'fro')^2);
% Maxwell-Boltzmann distribution (normalized)
energy16 = real(QAM_normSym16).^2 + imag(QAM_normSym16).^2;
MB_distr16 = exp(-lamba*energy16)/sum(exp(-lamba*energy16));
pd16 = makedist('Multinomial', 'Probabilities', MB_distr16);
QAM_MB_distr16 = pd16.Probabilities;
capacity_16 = awgn_channel_capacity(SNR_dB, QAM_normSym16, QAM_MB_distr16);
The "awgn_channel_capacity: function should compute the channel capacity (a single value) having as input.
integral2 is receiving in input a wrong mixture of symbolic expression and function-handle. I'm sure that the biggest problem is sum_den.
function capacity = awgn_channel_capacity(SNR, x, p)
syms y1 y2
snr = 10^(SNR/10) ; %linear
H_X = -sum(p .* log2(p)); %input-entropy
cond_H=0;
for i = 1 : length(x)
sum_den = 0;
for j = 1 : length(x)
sum_den = sum_den + p(j) * 1/pi * exp( -(abs(sqrt(snr)*x(j) - (y1 + 1i*y2)).^2));
end
cond_P = p(i) * 1/pi * exp( -(abs(sqrt(snr)*x(i) - (y1 + 1i*y2)).^2));
integrand = @(y1,y2) cond_P*log2(cond_P/sum_den(y1,y2));
cond_H = cond_H + integral2(integrand, -10,10,-10,10);
end
capacity = cond_H - H_X;
end

Answers (1)

UDAYA PEDDIRAJU
UDAYA PEDDIRAJU on 3 Nov 2023
Hi Nicolò,
I understand that you are having trouble implementing the formula for calculating the channel capacity of a QAM input over an AWGN channel in MATLAB. The issue seems to be with the integral and confusion between symbolic expressions and function handles.
To resolve this, you can modify your code as follows:
  1. Remove the “syms y1 y2” declaration since it is not needed.
  2. Define the integrand as a function handle that takes “y1” and “y2” as input arguments.
  3. Update the “integral” line to use “integral2” with the integrand function handle and the desired integration limits.
Here are the modified lines of code:
function capacity = awgn_channel_capacity(SNR, x, p)
% Your existing code...
for i = 1:length(x)
% Your existing code...
integrand = @(y1, y2) cond_P * log2(cond_P / sum_den);
cond_H = cond_H + integral2(integrand, -10, 10, -10, 10);
end
% Your existing code...
end
By making these changes, you should be able to successfully calculate the channel capacity.
I hope this resolves your issue!
  1 Comment
Nicolò Pisanu
Nicolò Pisanu on 4 Nov 2023
This is the modification you suggest to me.
%syms y1 y2
cond_H=0;
for i = 1 : length(x)
sum_den = 0;
for j = 1 : length(x)
sum_den = sum_den + p(j) * 1/pi * exp( -(abs(sqrt(snr)*x(j) - (y1 + 1i*y2)).^2));
end
cond_P = p(i) * 1/pi * exp( -(abs(sqrt(snr)*x(i) - (y1 + 1i*y2)).^2));
integrand = @(y1, y2) cond_P * log2(cond_P / sum_den);
%integrand = @(y1,y2) cond_P*log2(cond_P/sum_den(y1,y2));
cond_H = cond_H + integral2(integrand, -10, 10, -10, 10);
end
capacity = cond_H - H_X;
But unfortunately this is not yet the solution because I need y1 and y2 in the definition of sum_den, the denominator in the logarithm argument you can see in the main formula.
Unrecognized function or variable 'y1'
The problem is having to carry out a summation of integrals with a summation of its own inside!

Sign in to comment.

Products


Release

R2021b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!