How to solve this ‘’’Input to EIG must not contain NaN or Inf’’’. problem
5 min readOct 19, 2023
This is the code:
% 0: optimal BF
% 0q: optimal BF
clearvars;
close all;
filename='stat';
seed=1970;
rng(seed);
% rng(107);
cl=3;
ray=10;
deg=5;
path=cl*ray;
t=16; %number of transmit antennas
r=16; %number of receive antennas
dt=1/2*2*pi; %normalized transmit antenna spacing
dr=1/2*2*pi; %normalized receive antenna spacing
phi_itr=100;mu_itr=100; itr=phi_itr*mu_itr;
ULA=0; %1 for ULA and 0 for UPAsector=4
if sector==1
Ta=[0 2*pi];Te=[0 2*pi];
elseif sector==2
Ta=[0 pi];
elseif sector==3
Ta=[pi/6 5*pi/6];
elseif sector==4
Ta=[0 2*pi];Te=[0 2*pi];
elseif sector==5
Ta=[pi/3 2*pi/3]; Te=[pi/2, pi/2+pi/4];
end
num=1;
if num==1
range=[0:3:15]; % range of angular spread in deg
%range=3; % range of angular spread in deg
disp('~~SNR vs angular spread');
elseif num==2
range=1:8; % number of cluster
disp('~~SNR vs number of cluster');
elseif num==3
range=1:10; % number of ray per cluster
disp('~~SNR vs number of rays (per cluster)');
elseif num==4
range=2:7; %number of antenna
disp('~~SNR vs number of antennas');
elseif num==5
range=1; %
disp('~~BER vs Pt/N0');
snr_range=-25:2:5;
%snr_range=0:2:25;
qam=4;
pt=10.^(snr_range/10);
pt=pt*3/(qam-1);qam_const=1;
ber0=zeros(size(snr_range));ber0q=ber0;ber1=ber0;ber1w=ber0;ber1q=ber0;ber_ay=ber0;ber_ay1=ber0;elseif num==6
range=0;
disp('~~capacity~~');
pt=10.^(range/10);
cp0=zeros(1,itr);cp0q=cp0;
cp1=cp0;cp1q=cp0;
cp_ay=cp0;cp_ay1=cp0;elseif num==7
range=0;
disp('~~transmission rate~~');
snr_range=-25:2:5;
pt=10.^(snr_range/10);
Rb0=zeros(size(snr_range));Rb0q=Rb0;Rb1=Rb0;Rb2=Rb0;Rb1q=Rb0;Rb1s=Rb0; Rb_ay=Rb0;Rb_ay_bf=Rb0;
end
Nth=8; Bf=4; Nrf=1;
Rh=zeros(t,t);asnr0=zeros(size(range)); asnr0q=asnr0; asnr1=asnr0;asnr2=asnr0; asnr2q=asnr0; asnr10=asnr0; asnr1q=asnr0;asnr1a=asnr0;asnr1b=asnr0;if ULA==1
tW=1;
tH=t;
else
if ceil(sqrt(t))==sqrt(t)
tW=sqrt(t);
tH=tW;
else
tW=sqrt(2*t);
tH=tW/2;
end
end
if ceil(sqrt(r))==sqrt(r)
rW=sqrt(r);
rH=rW;
else
rW=sqrt(2*r);
rH=rW/2;
end
index=0;
for nn=1:length(range)
if num==1
deg=range(nn); % angular spread in deg
disp(strcat('~~~~~angular spread=',num2str(deg)));
elseif num==2
cl=range(nn); % number of cluster
elseif num==3
ray=range(nn); % number of ray per cluster
elseif num==4
t=2^range(nn); %number of antenna
r=t;
if ceil(sqrt(t))==sqrt(t)
tW=sqrt(t);tH=tW;
else
tW=sqrt(2*t);tH=tW/2;
end
if ceil(sqrt(r))==sqrt(r)
rW=sqrt(r);rH=rW;
else
rW=sqrt(2*r);rH=rW/2;
end
end
sigma_ta=deg*ones(1,cl)/180*pi;% standard deviation of phi_t
sigma_te=sigma_ta;% standard deviation of phi_t
sigma_ra=sigma_ta;
sigma_re=sigma_ta;
sig_ta=deg*ones(1,path)/180*pi;
sig_te=deg*ones(1,path)/180*pi;
for n_mu=1:mu_itr
mu_ta=Ta(1)+(Ta(2)-Ta(1))*rand(cl,1);
mu_ra=2*pi*rand(cl,1);
mu_te=Te(1)+(Te(2)-Te(1))*rand(cl,1);
mu_re=2*pi*rand(cl,1); alpha=ones(1,path);
b=stat_Rh(sector,ULA,t,r,cl,ray,dt,tH,tW,alpha,mu_ta,mu_te,sig_ta,sig_te,Ta,Te);
[V,D]=eig(b);
[la1,i0]=max(diag(D));
f1=V(:,i0);
if ULA
At=exp(1i*dt*(0:t-1)'*cos(mu_ta)');
[f10,~]=eigs(At*At',1);
end
A=zeros(t,t);
f2=zeros(t,1);
for i_cl=1:cl
temp_t=exp(1i*dt*(0:tH-1)'*cos(mu_te(i_cl)))*exp(1i*dt*(0:tW-1)*sin(mu_te(i_cl))*sin(mu_ta(i_cl)));
temp_t=reshape(temp_t,t,1);
f2=f2+temp_t;
A=A+temp_t*temp_t';
end
[f10,~]=eigs(A,1);
f2=f2/norm(f2);
asnr1a(nn)=asnr1a(nn)+la1;
Rh=zeros(t,t);
for n_itr=1:phi_itr
index=index+1;
if sector==1
[H, At, Ar]=stat_mmCh(ULA,dr,dt,t,r,tW,tH,rW,rH, path, cl,ray,mu_ta, mu_te, mu_ra, mu_re, sigma_ta,sigma_te,sigma_ra,sigma_re);
else
[H, At, Ar,spc_phi_ta,spc_phi_te]=stat_mmCh_sector(Ta,Te,ULA,dr,dt,t,r,tW,tH,rW,rH, path, cl,ray,mu_ta, mu_te, mu_ra, mu_re, sigma_ta,sigma_te,sigma_ra,sigma_re);
end
[V,D]=eig(H'*H);
[la,i0]=max(diag(D));
f0=V(:,i0);
Rh=Rh+H'*H;
snr0=la;
% f0q=exp(1i*round(angle(f0)/2/pi*Nth)*2*pi/Nth);
f0=f0/max(abs(f0))*2; f0q=q1(f0,Nth,Nth);
g0q=H*f0q;
% g0q=exp(1i*round(angle(g0q)/2/pi*Nth)*2*pi/Nth);
g0q=q1(g0q/max(abs(g0q))*2,Nth,Nth);
f1q=q1(f1/max(abs(f1))*2,Nth,Nth);
% f1q=exp(1i*round(angle(f1)/2/pi*Nth)*2*pi/Nth);
g1q=H*f1q;
% g1q=exp(1i*round(angle(g1q)/2/pi*Nth)*2*pi/Nth);
g1q=q1(g1q/max(abs(g1q))*2,Nth,Nth);
snr0q=norm(g0q'*H*f0q,2)^2/norm(g0q,2)^2/norm(f0q,2)^2; %using quanized optimal bf snr1=norm(H*f1,2)^2/norm(f1,2)^2;;
snr1q=norm(g1q'*H*f1q,2)^2/norm(g1q,2)^2/norm(f1q,2)^2; %quantized snr10=norm(H*f10,2)^2;
snr2=norm(H*f2,2)^2;
if num<=4
asnr0(nn)=asnr0(nn)+snr0;
asnr0q(nn)=asnr0q(nn)+snr0q;
asnr1(nn)=asnr1(nn)+snr1;
asnr2(nn)=asnr2(nn)+snr2;
asnr10(nn)=asnr10(nn)+snr10;
asnr1q(nn)=asnr1q(nn)+snr1q; end
if num==5
[snr_ay,snr_ay1]= Ayach(1,t,dt,tH,tW,Nrf,H,ULA,At,Ar,f0,Ta,Te, 4, Nth,spc_phi_ta,spc_phi_te);
[snr_ay_bf,snr_ay1]= Ayach(1,t,dt,tH,tW,Nrf,H,ULA,At,Ar,f0,Ta,Te, 6, Nth,spc_phi_ta,spc_phi_te);
ber0=ber0+qam_const*Q_fun(sqrt(pt*snr0));
ber0q=ber0q+qam_const*Q_fun(sqrt(pt*snr0q));
ber1=ber1+qam_const*Q_fun(sqrt(pt*snr1));
ber2=ber2+qam_const*Q_fun(sqrt(pt*snr2));
% ber1w=ber1w+qam_const*Q_fun(sqrt(pt*snr1w));
ber1q=ber1q+qam_const*Q_fun(sqrt(pt*snr1q));
% ber1s=ber1s+qam_const*Q_fun(sqrt(pt*snr1s));
% ber1sq=ber1sq+qam_const*Q_fun(sqrt(pt*snr1sq)); ber_ay=ber_ay+qam_const*Q_fun(sqrt(pt*snr_ay));
ber_ay1=ber_ay1+qam_const*Q_fun(sqrt(pt*snr_ay1));
end
if num==6
[snr_ay,snr_ay1]= Ayach(1,t,dt,tH,tW,Nrf,H,ULA,At,Ar,f0,Ta,Te, Bf, Nth,spc_phi_ta,spc_phi_te); cp0(index)=snr0;
cp0q(index)=snr0q;
cp1(index)=snr1;
cp1q(index)=snr1q;
cp_ay(index)=snr_ay;
cp_ay1(index)=snr_ay1;end
if num==7
[snr_ay,snr_ay1]= Ayach(1,t,dt,tH,tW,Nrf,H,ULA,At,Ar,f0,Ta,Te, Bf, Nth,spc_phi_ta,spc_phi_te);
[snr_ay_bf,snr_ay1]= Ayach(1,t,dt,tH,tW,Nrf,H,ULA,At,Ar,f0,Ta,Te, Bf+1, Nth,spc_phi_ta,spc_phi_te);
Rb0=Rb0+log2(1+pt*snr0);
Rb0q=Rb0q+log2(1+pt*snr0q);
Rb1=Rb1+log2(1+pt*snr1);
Rb2=Rb2+log2(1+pt*snr2);
Rb1q=Rb1q+log2(1+pt*snr1q);
Rb_ay=Rb_ay+log2(1+pt*snr_ay);
Rb_ay_bf=Rb_ay_bf+log2(1+pt*snr_ay_bf);
Rb1s=Rb1s+log2(1+pt*la1);
end end Rh=Rh/phi_itr;
temp=eig(Rh);
asnr1b(nn)=asnr1b(nn)+max(temp); end
end%%
if num<=4
asnr0=10*log10(asnr0/itr);
asnr0q=10*log10(asnr0q/itr);
asnr1=10*log10(asnr1/itr);
asnr1q=10*log10(asnr1q/itr);
asnr2=10*log10(asnr2/itr);
asnr2q=10*log10(asnr2q/itr);
asnr10=10*log10(asnr10/itr);
asnr1b=10*log10(asnr1b/mu_itr);
asnr1a=10*log10(asnr1a/mu_itr);
%% spawc
figure;
plot(range, asnr0,'k--o'); hold on;
%plot(range, asnr0q,'k->');
%plot(range, asnr1b,'b-s');
plot(range, asnr1,'r--d');
%plot(range, asnr1q,'r--*');
%plot(range, asnr2q,'g--+');%legend('full CSI','full CSI,quantized','opt stat','stat','stat,quantized');if num==1
xlabel('Angular spread (deg)');% range of angular spread in deg
elseif num==2
xlabel('Number of clusters');% number of cluster
elseif num==3
xlabel('Number of rays per cluster'); % number of ray per cluster
else
xlabel('Number of antennas'); %number of antenna
end
%axis([0 15 18 25]); %for t=16
axis([0 15 22 32]); %t=64, r=16
%axis([0 15 28 35]); %t=r=64
ylabel('Average SNR (dB)');
grid on;
%%
figure;
plot(range, asnr0,'k'); hold on;
plot(range, asnr1,'g--d');
plot(range, asnr2,'r--+');
plot(range, asnr10,'g--o');
plot(range, asnr1a,'r-*');
plot(range, asnr1b,'b-s');
legend('optimal','statistical','multibeam','stat,zero var','theoretical','opt stat');if num==1
xlabel('Angular spread (deg)');% range of angular spread in deg
elseif num==2
xlabel('Number of clusters');% number of cluster
elseif num==3
xlabel('Number of rays per cluster'); % number of ray per cluster
else
xlabel('Number of antennas'); %number of antenna
end
ylabel('Average SNR (dB)');
grid on; temp=clock;
st=strcat(num2str(temp(2)),'.',num2str(temp(3)),'.',num2str(temp(4)),'.',num2str(temp(5)));
st=strcat(st,',stat-new.m, BER, t=',num2str(t),',r=',num2str(r));
st=strcat(st,',cl=',num2str(cl));
st=strcat(st,',ray=',num2str(ray(1)));
st=strcat(st,',deg=',num2str(sigma_ta(1)/pi*180));
st1=strcat('Ta=',num2str(Ta/pi),'\pi');
st1=strcat(st1,',Te=',num2str(Te/pi),'\pi');
st1=strcat(st1,',mu_itr=',num2str(mu_itr),',phi_itr=',num2str(phi_itr));
if ULA==1
st1=strcat(st1,',ULA');
else
st1=strcat(st1,',UPA');
end
st1=strcat(st1,',Bf=',num2str(Bf),',Nth=',num2str(Nth));
title({st,st1});
grid onend
%%
if num==4
figure;
range=2.^range;
plot(range, asnr0,'k-.o'); hold on;
plot(range, asnr0q,'k->');
plot(range, asnr1,'r-.d');
plot(range, asnr1q,'r-*');
legend('FULL CSI','FULL CSI,quantized','stat','stat,quantized');
xlabel('Number of antennas');
ylabel('Average SNR (dB)');
grid on;end
%%
if num==5
ber0=ber0/itr;
ber0q=ber0q/itr;
ber1=ber1/itr;
ber2=ber2/itr;
ber1q=ber1q/itr;
ber_ay=ber_ay/itr;
ber_ay1=ber_ay1/itr;
%%
figure;
semilogy(snr_range, ber0, 'k-');hold on;
semilogy(snr_range, ber0q, 'k-*');
semilogy(snr_range, ber1, 'r-s');
semilogy(snr_range, ber1q, 'r-*');
semilogy(snr_range, ber_ay, 'g-d');
axis([min(snr_range) max(snr_range) 10^(-5) 1]);
legend('0','0q','1','1q','ay');
temp=clock;
st=strcat(num2str(temp(2)),'.',num2str(temp(3)),'.',num2str(temp(4)),'.',num2str(temp(5)));
st=strcat(st,',stat-new.m, BER, t=',num2str(t),',r=',num2str(r));
st=strcat(st,',cl=',num2str(cl));
st=strcat(st,',ray=',num2str(ray(1)));
st=strcat(st,',deg=',num2str(sigma_ta(1)/pi*180));
st1=strcat('Ta=',num2str(Ta/pi),'\pi');
st1=strcat(st1,',Te=',num2str(Te/pi),'\pi');
st1=strcat(st1,',mu_itr=',num2str(mu_itr),',phi_itr=',num2str(phi_itr));
if ULA==1
st1=strcat(st1,',ULA');
else
st1=strcat(st1,',UPA');
end
st1=strcat(st1,',Bf=',num2str(Bf),',Nth=',num2str(Nth));
title({st,st1});
grid on
end%% plot cdfif num==6 cp0=log2(1+cp0*pt);
cp0q=log2(1+cp0q*pt);
cp1=log2(1+cp1*pt);
cp1q=log2(1+cp1q*pt);
cp_ay=log2(1+cp_ay*pt);
cp_ay1=log2(1+cp_ay1*pt);
[pdf0,x0]=hist(cp0,100);
[pdf0q,x0q]=hist(cp0q,100);
[pdf1,x1]=hist(cp1,100);
[pdf1q,x1q]=hist(cp1q,100);
[pdf_ay,x_ay]=hist(cp_ay,100);
[pdf_ay1,x_ay1]=hist(cp_ay1,100);
pdf0=pdf0/itr;
pdf0q=pdf0q/itr;
pdf1=pdf1/itr;
pdf1q=pdf1q/itr;
pdf_ay=pdf_ay/itr;
pdf_ay1=pdf_ay1/itr;
cdf0=zeros(size(pdf0));cdf0q=cdf0;cdf1=cdf0; cdf1q=cdf0; cdf_ay=cdf0;cdf_ay1=cdf0;
cdf0(1)=pdf0(1);
cdf0q(1)=pdf0q(1);
cdf1(1)=pdf1(1);
cdf1q(1)=pdf1q(1);
cdf_ay(1)=pdf_ay(1);cdf_ay1(1)=pdf_ay1(1);
for ii=2:100
cdf0(ii)=cdf0(ii-1)+pdf0(ii-1);
cdf0q(ii)=cdf0q(ii-1)+pdf0q(ii-1);
cdf1(ii)=cdf1(ii-1)+pdf1(ii-1);
cdf1q(ii)=cdf1q(ii-1)+pdf1q(ii-1);
cdf_ay(ii)=cdf_ay(ii-1)+pdf_ay(ii-1);
cdf_ay1(ii)=cdf_ay1(ii-1)+pdf_ay1(ii-1);
end
%%
figure;
plot(x0,cdf0,'k'); hold;
plot(x0q,cdf0q,'k:');
plot(x1,cdf1,'r');
plot(x1q,cdf1q,'--');
plot(x_ay,cdf_ay,'-.m');
plot(x_ay1,cdf_ay1,'--g');
legend('cdf0','cdf0q','cdf1','cdf1q','SPC, N_{rf}=1','SPC, N_{rf}=2','Location','NorthWest');
temp=clock;
st=strcat(num2str(temp(2)),'.',num2str(temp(3)),'.',num2str(temp(4)),'.',num2str(temp(5)));
st=strcat(st,filename,'cdf, t=',num2str(t),',r=',num2str(r));
st=strcat(st,',cl=',num2str(cl));
% st=strcat(st,',ray=',num2str(ray));
st=strcat(st,',deg=',num2str(sigma_ta(1)/pi*180));
st=strcat(st,',pt(dB)=',num2str(range));
st=strcat(st,',Bf=',num2str(Bf));
st=strcat(st,',Nth=',num2str(Nth));
st=strcat(st,',Ta=',num2str(Ta/pi),'\pi');
st=strcat(st,',Te=',num2str(Te/pi),'\pi');
if ULA==1
st1=strcat('ULA,','itr=');
else
st1=strcat('UPA,','itr=');
end st1=strcat(st1,num2str(itr));
title({st,st1}); grid on
axis([2 10 0 1.1]);
grid on;
end
if num==7
Rb0=Rb0/itr;
Rb0q=Rb0q/itr;
Rb1s=Rb1s/itr;
Rb1=Rb1/itr;
Rb2=Rb2/itr;
Rb1q=Rb1q/itr;
Rb_ay=Rb_ay/itr;
Rb_ay_bf=Rb_ay_bf/itr;
%%
figure;
plot(snr_range, Rb0, 'k-.');hold on;
plot(snr_range, Rb0q, 'k-*');
plot(snr_range, Rb1, 'r-d');
plot(snr_range, Rb2, 'r-+');
plot(snr_range, Rb1q, 'r-s');
plot(snr_range, Rb1s, 'b-x')
plot(snr_range, Rb_ay, 'm->'); legend('optimal','0q','1','2multibeam','1q','1s','ay');
%axis([min(snr_range) max(snr_range) 10^(]);
temp=clock;
st=strcat(num2str(temp(2)),'.',num2str(temp(3)),'.',num2str(temp(4)),'.',num2str(temp(5)));
st=strcat(st,'stat-new.m, Rb, t=',num2str(t),',r=',num2str(r));
st=strcat(st,',cl=',num2str(cl));
st=strcat(st,',ray=',num2str(ray));
st=strcat(st,',deg=',num2str(sigma_ta(1)/pi*180));
if ULA==1
st1=strcat('ULA,','itr=');
else
st1=strcat('UPA,','itr=');
end st1=strcat(st1,num2str(itr));
title({st,st1});
grid on
%%
figure;
plot(snr_range, Rb0, 'k-.o');hold on;
plot(snr_range, Rb0q, 'k->');
plot(snr_range, Rb_ay_bf, 'm:o');
plot(snr_range, Rb1, 'r-.d');
plot(snr_range, Rb1q, 'r-*');
plot(snr_range, Rb_ay, 'b-+');
legend('Full CSI','Full CSI, quantized','SPC, 10 bits','stat','stat, quantized','SPC, 8bits');
axis([-20 0 1 10]) ;
grid on
ylabel('Transmission rate (bits)');
xlabel('P_t/N_0');end
NOTE:-
Matlabsolutions.com provide latest MatLab Homework Help,MatLab Assignment Help , Finance Assignment Help for students, engineers and researchers in Multiple Branches like ECE, EEE, CSE, Mechanical, Civil with 100% output.Matlab Code for B.E, B.Tech,M.E,M.Tech, Ph.D. Scholars with 100% privacy guaranteed. Get MATLAB projects with source code for your learning and research.
num=1;
if num==1
range=[0:3:15]; % range of angular spread in degTherefore range will start with 0.
for nn=1:length(range)
if num==1
deg=range(nn); % angular spread in degnum is still 1, so the angular speed case applies. range(nn) will be accessed, and as indicated earlier that starts with 0. So deg will be assigned 0.
sig_ta=deg*ones(1,path)/180*pi;
sig_te=deg*ones(1,path)/180*pi;with deg being 0, this makes sig_ta an sig_te all 0.
b=stat_Rh(sector,ULA,t,r,cl,ray,dt,tH,tW,alpha,mu_ta,mu_te,sig_ta,sig_te,Ta,Te);and there they are passed to stat_Rh.
Inside stat_Rh
SEE COMPLETE ANSWER CLICK THE LINK
