how to plot in matlab a formula with multiple summation ?

4 Ansichten (letzte 30 Tage)
Fares Zaritt
Fares Zaritt am 19 Mai 2023
Kommentiert: Fares Zaritt am 11 Jun. 2023
i have the following Spectral efficiency in line of sight (SE_NLOS) formula :
i have K angles Q0 random value between (0;2Pi) we'll refer to them as Q in my code ,and K random angles Q1 random in value between (0;2PI) we'll refer to them as W in my code,
in the first G(.,.) , G1( Q(k), Q(i) ) where i and k are the i_th and k_th Q angles
and the second G(.,.), G2( Q(k), W(i) ) where i is the i_th angle of W, and k is the k_th Q angles
due do the random nature of generating the random angles Q and W, i calculate the mean of G over many realizations for accuracy
beta_bar and SNR0 are given in dB we convert them to linear scale (and to avoid any errors in the case of 1/SNR0)
the plotted results should be similar to that present in the figure below
i only attempted to simulate SE_LOS for M=100 so it's the upper most black line that we need to pay attention to
(i'am yet to attemp recreating the code for the NLOS case (blue line) but considering my evident incompetence when it comes to loops and summations it won't be long before i post a question about it)
my best attempt code is below:
% parameters
dH = 0.5;
M = 100;
K_UE = 1:20;
SNR0_dB = 0;
SNR0 = 10^(0/10);
beta_bar_dB = -10;
beta_bar = 10^(-10/10);
num_realizations = 1000;
% Compute SE_LOS for each value of K
% SE_LOS = zeros(length(K));
% SE_LOS=0;
SE_LOS = 0;
G1_sum = 0;
G2_sum = 0;
TOT_SUM = 0;
for kk = 1:numel(K_UE)
K = K_UE(kk);
for k = 1:K % for the innitial sum k to K
G1 = zeros(1, num_realizations);
G2 = zeros(1, num_realizations);
for i = 1:K % for the G sums i to K per each k
for j = 1:num_realizations % we calculate the mean of multiple realizations of G1 & G2 for more accuracy
Q = 2*180*rand(1:K); % Q angles of each intra_cell UE
W = 2*180*rand(1:K); % W angles of each inter_cell UE
if sind(Q(k)) ~= sind(W(i))
G2(j) = (sind(180*dH*M*(sind(Q(k))-sind(W(i)))))^2/ ...
(M*(sind(180*dH*(sind(Q(k))-sind(W(i)))))^2);
else
G2(j) = M;
end
if i ~= k
G1(j) = (sind(pi*dH*M*(sind(Q(k))-sind(Q(i)))))^2/ ...
(M*(sind(pi*dH*(sind(Q(k))-sind(Q(i)))))^2);
else
G1(j) = 0;
end
end
G2_sum = G2_sum + mean(G2);
G1_sum = G1_sum + mean(G1);
end
TOT_SUM = TOT_SUM + mean( log2(1 + M ./ (G1_sum + beta_bar * G2_sum + 1./SNR0)) );
G2_sum = 0;
G1_sum = 0;
end
SE_LOS(kk) = TOT_SUM
TOT_SUM=0;
end
% Plot SE_LOS as function of K_UE
plot(K_UE, SE_LOS);
xlabel('Number of UEs (K)');
ylabel('Sum SE (bits/s/Hz)');
when it comes to matlab summations and loops continue to be the bane of my existence when it comes to matlab
of coure all assistance is greatly appreciated!
  3 Kommentare
VBBV
VBBV am 20 Mai 2023
May be you can do something like this
% parameters
dH = 0.5;
M = 100;
K_UE = 1:20;
SNR0_dB = 0;
SNR0 = 10^(0/10);
beta_bar_dB = -10;
beta_bar = 10^(-10/10);
num_realizations = 100;
% Compute SE_LOS for each value of K
% SE_LOS = zeros(length(K));
% SE_LOS=0;
SE_LOS = 0;
G1_sum = 0;
G2_sum = 0;
TOT_SUM = 0;
for kk = 1:numel(K_UE)
K = K_UE(kk);
for k = 1:K % for the innitial sum k to K
G1 = zeros(1, num_realizations);
G2 = zeros(1, num_realizations);
Q1 = 2*180*rand(K,1); %
for i = 1:K % for the G sums i to K per each k
Q2 = 2*180*rand(K,1); %
Q = 2*180*rand(K,1);
W = 2*180*rand(K,1); % W angles of each inter_cell UE
for j = 1:num_realizations % we calculate the mean of multiple realizations of G1 & G2 for more accuracy
if sind(Q(k)) ~= sind(W(i))
G2(j) = (sind(180*dH*M*(sind(Q(k))-sind(W(i)))))^2/ ...
(M*(sind(180*dH*(sind(Q(k))-sind(W(i)))))^2);
else
G2(j) = M;
end
if i~=k
G1(j) = (sind(180*dH*M*(sind(Q1(k))-sind(Q2(i)))))^2/ ...
(M*(sind(180*dH*(sind(Q1(k))-sind(Q2(i)))))^2);
else
G1(j) = 0;
end
end
G2_sum = G2_sum + mean(G2);
G1_sum = G1_sum + mean(G1);
end
TOT_SUM = TOT_SUM + ( log2(1 + M ./ (G1_sum + beta_bar * G2_sum + (1./SNR0))) );
G2_sum = 0;
G1_sum = 0;
end
SE_LOS(kk) = TOT_SUM;
TOT_SUM=0;
end
G1, G2
G1 = 1×100
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
G2 = 1×100
1.0e-06 * 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086 0.6086
% Plot SE_LOS as function of K_UE
plot(K_UE, SE_LOS);
grid
xlabel('Number of UEs (K)');
ylabel('Sum SE (bits/s/Hz)');
Fares Zaritt
Fares Zaritt am 11 Jun. 2023
from my experience trying to to saparating the Q(i) and Q(k) and having Q(k) outside the
"for j = 1:num_realizations" always results in a jagged plot, although it seems to be converging towards a value similar to that of the original figure

Melden Sie sich an, um zu kommentieren.

Antworten (1)

VBBV
VBBV am 20 Mai 2023
% parameters
dH = 0.5;
M = 100;
K_UE = 1:20;
SNR0_dB = 0;
SNR0 = 10^(0/10);
beta_bar_dB = -10;
beta_bar = 10^(-10/10);
num_realizations = 100;
% Compute SE_LOS for each value of K
% SE_LOS = zeros(length(K));
% SE_LOS=0;
SE_LOS = 0;
G1_sum = 0;
G2_sum = 0;
TOT_SUM = 0;
for kk = 1:numel(K_UE)
K = K_UE(kk);
for k = 1:K % for the innitial sum k to K
G1 = zeros(1, num_realizations);
G2 = zeros(1, num_realizations);
for i = 1:K % for the G sums i to K per each k
for j = 1:num_realizations % we calculate the mean of multiple realizations of G1 & G2 for more accuracy
Q = 2*180*rand(K,1); % Q angles of each intra_cell UE
W = 2*180*rand(K,1); % W angles of each inter_cell UE
if sind(Q(k)) ~= sind(W(i))
G2(j) = (sind(180*dH*M*(sind(Q(k))-sind(W(i)))))^2/ ...
(M*(sind(180*dH*(sind(Q(k))-sind(W(i)))))^2);
else
G2(j) = M;
end
%
% if i ~= k
% G1(j) = (sind(pi*dH*M*(sind(Q(k))-sind(Q(i)))))^2/ ...
% (M*(sind(pi*dH*(sind(Q(k))-sind(Q(i)))))^2);
% else
% G1(j) = 0;
% end
end
G2_sum = G2_sum + mean(G2);
G1_sum = G1_sum + mean(G1);
end
TOT_SUM = TOT_SUM + ( log2(1 + M ./ (G1_sum + beta_bar * G2_sum + (1./SNR0))) );
G2_sum = 0;
G1_sum = 0;
end
SE_LOS(kk) = TOT_SUM;
TOT_SUM=0;
end
% Plot SE_LOS as function of K_UE
plot(K_UE, SE_LOS);
grid
xlabel('Number of UEs (K)');
ylabel('Sum SE (bits/s/Hz)');
  2 Kommentare
VBBV
VBBV am 20 Mai 2023
% if i ~= k
% G1(j) = (sind(pi*dH*M*(sind(Q(k))-sind(Q(i)))))^2/ ...
% (M*(sind(pi*dH*(sind(Q(k))-sind(Q(i)))))^2);
% else
% G1(j) = 0;
% end
The above code may be superfluous , otherwise code looks ok except for these lines below
Q = 2*180*rand(1:K); % Q angles of each intra_cell UE
W = 2*180*rand(1:K);
K is a scalar, when you use it as 1:K it becomes a vector which is interpreted incorrectly by rand function. Change it to
Q = 2*180*rand(K,1); % Q angles of each intra_cell UE
W = 2*180*rand(K,1);
Fares Zaritt
Fares Zaritt am 11 Jun. 2023
hello thank you for the answer! i also apologies for the super late response, some more urgent problems popped up and i had to abandone this for a while, and then forgot i even asked on the forum,
/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
%
% if i ~= k
% G1(j) = (sind(pi*dH*M*(sind(Q(k))-sind(Q(i)))))^2/ ...
% (M*(sind(pi*dH*(sind(Q(k))-sind(Q(i)))))^2);
% else
% G1(j) = 0;
% end
unfortunatly can't be a superfluoues segment, since its a relevent part of the equation,
i took into considuration
Q = 2*180*rand(K,1); % Q angles of each intra_cell UE
W = 2*180*rand(K,1);
and after fixing some other errors that were in the code (mainly using the cammand sinD while still using Pi instead of 180) , with added code of NLOS LOWER BOUND to serve as a point of reference, and the resulting code is:
(note that K_UE = 1:10; instead of K_UE = 1:20; because setting K_UE to 1:20 takes way too long to compute, so i set it this way to conserve on some time)
% SE_LOS FOR multiple values of M, FORMULA ON PAGE (196), FIGURE ON PAGE (199)
% parameters
dH = 0.5;
M_values = [10 100];
K_UE = 1:10;
SNR0_dB = 0;
SNR0 = 10^(0/10);
beta_bar_dB = -10;
beta_bar = 10^(-10/10);
num_realizations = 8000;
% Compute SE_LOS and SE_NLOS_LOWER_BOUND for each value of K and M
SE_NLOS_LOWER_BOUND = zeros(numel(M_values), numel(K_UE));
SE_LOS = zeros(numel(M_values), numel(K_UE));
for mm = 1:numel(M_values)
M = M_values(mm);
for kk = 1:numel(K_UE)
K = K_UE(kk);
TOT_SUM = 0;
for k = 1:K % for the initial sum k to K
G1_sum = 0;
G2_sum = 0;
for i = 1:K % for the G sums i to K per each k
G1 = zeros(1, num_realizations);
G2 = zeros(1, num_realizations);
for j = 1:num_realizations % we calculate the mean of multiple realizations of G1 & G2 for more accuracy
W = 2*180*rand(K,1); % W angles of each inter_cell UE
Q = 2*180*rand(K,1); % Q angles of each intra_cell UE
if sind(Q(k)) ~= sind(W(i))
G2(j) = (sind(180*dH*M*(sind(Q(k))-sind(W(i)))))^2./ ...
(M*(sind(180*dH*(sind(Q(k))-sind(W(i)))))^2);
else
G2(j) = M;
end
if i ~= k && sind(Q(k)) ~= sind(Q(i))
G1(j) = (sind(180*dH*M*(sind(Q(k))-sind(Q(i)))))^2./ ...
(M*(sind(180*dH*(sind(Q(k))-sind(Q(i)))))^2);
else
G1(j) = 0;
end
end
G2_sum = G2_sum + mean(G2);
G1_sum = G1_sum + mean(G1);
end
TOT_SUM = TOT_SUM + ( log2(1 + M ./ (G1_sum + beta_bar .* G2_sum + (1./SNR0))) );
end
SE_LOS(mm, kk) = TOT_SUM;
SE_NLOS_LOWER_BOUND(mm, kk) = K * log2(1 + ((M - 1) ./ ((K - 1) + K * beta_bar + (1/SNR0))));
end
% Plot SE_LOS and SE_NLOS_LOWER_BOUND for the current M value
plot(K_UE, SE_LOS(mm, :), 'LineWidth', 1.5);
hold on;
plot(K_UE, SE_NLOS_LOWER_BOUND(mm, :),'red', 'LineWidth', 1.5);
end
hold off;
grid on;
xlabel('Number of UEs (K)');
ylabel('Sum SE (bits/s/Hz)');
legend('M = 10 (LOS)', 'M = 10 (NLOS lowerbound)', 'M = 100 (LOS)', 'M = 100 (NLOS lowerbound)','Location', 'northwest');
the higher the numaber realizations the smoother the LOS curves becomes, although i still dont understand why the resulting curve gives the wrong values and doesnt follow the the same path as the one in figure i provided in my origin post,
at the point K_UE =10, SE should've reached the value around 9 and 45 respectively, but my sim gives at that point 7.48 and 26.38 respectively,not to mention LOS values should be luch greater than that of NLOS LOWER BOUND
may have to repost the question, again but regardless thank you for giving it your time i greatly appretiate it!

Melden Sie sich an, um zu kommentieren.

Kategorien

Mehr zu Signal Processing finden Sie in Help Center und File Exchange

Produkte


Version

R2014a

Community Treasure Hunt

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

Start Hunting!

Translated by