Help to convert formula into matlab code - QAM channel capacity

5 次查看（过去 30 天）
Nicolò Pisanu
Nicolò Pisanu 2023-10-24
评论： Nicolò Pisanu 2023-11-4
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

请先登录，再回答此问题。

回答（1 个）

UDAYA PEDDIRAJU
UDAYA PEDDIRAJU 2023-11-3
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 个评论
Nicolò Pisanu
Nicolò Pisanu 2023-11-4
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!

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Formula Manipulation and Simplification 的更多信息

标签

产品


版本

R2021b

Community Treasure Hunt

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

Start Hunting!

Translated by