📄 cslmmse.asv
字号:
clc
clear all
Nt = 1; % the number of transmitted antennas
Nr = 1; % the number of transmitted antennas
T=2; % the number of symbol periods
Tp=Nt;
Td=T-Tp;
SNR = 0:5:70; % the system signal-to-noise ratio with dB
snr = 10.^(0.1*SNR); % the signal-to-noise with normal scale
snrp=snr;
MM=100; % 1000 channel
NN = 100;% using 10000 Monte-Carlo runs
hm=zeros(Nt*Nr,1); % channel mean
hh=eye(Nt*Nr); % channel covariance
vv=eye(Tp*Nr);
nn=eye(Td*Nr); % noise covariance
nn_ray=eye(T*Nr); % rayleigh fading with perfect csi and without training symbol
capa_ray = zeros(1,length(SNR)); % rayleigh fading with perfect csi and without training symbol
A_ray=zeros(T*Nr,T*Nr); % rayleigh fading with perfect csi and without training symbol
capa1 = zeros(1,length(SNR)); % perfect csi
A10=zeros(Td*Nr,Td*Nr);
capa2 = zeros(1,length(SNR)); % worst case noise lower bound
A20=zeros(Td*Nr,Td*Nr);
AA20=zeros(Td*Nr,Td*Nr);
capa3 = zeros(1,length(SNR)); % inequality lower bound
capa30=zeros(1,length(SNR));
capa300=zeros(1,length(SNR));
A30=zeros(Td*Nr,Td*Nr);
AA30=zeros(Td*Nr,Td*Nr);
aa300=0;
for ii = 1:length(SNR)
waitbar(ii/length(SNR));
% LMMSE estimator
p=eye(Nt*Nr); % training symbol
pp=inv(p'*p)*p'; % wei ni
F=sqrt(snrp(ii))*hh*p'*inv(snrp(ii)*p*hh*p'+vv');
f0=(eye(size(F*pp))-sqrt(snrp(ii))*F*pp)*hm;
for jj=1:MM
Hw=sqrt(1/2)*(randn(Nr*Nt,1) + j*randn(Nr*Nt,1)); % channel
v=sqrt(1/2)*(randn(Nr*Tp,1) + j*randn(Nr*Tp,1)); % training noise
z=sqrt(snrp(ii))*p*Hw+v; % receive
g=F*z+f0; % estimation
B=eye(Nt*Nr);
C=inv(eye(Nt*Nr)+snrp(ii)*hh*p'*inv(vv)*p);
d=B*g+C*hm-B*f0;
for kk =1:NN
X_ray =sqrt(1/(2*Nt))* (randn(Nt,T + j* randn(Nt,T)); % rayleigh fading with perfect csi and without training symbol
XX_ray=kron(X_ray.',eye(Nr));
X =sqrt(1/(2*Nt))* (randn(Nt,Td) + j* randn(Nt,Td));
XX=kron(X.',eye(Nr));
A_ray=A_ray+XX_ray*Hw*Hw'*XX_ray'; % rayleigh fading with perfect csi and without training symbol
A10=A10+XX*Hw*Hw'*XX'; % perfect csi
A20=A20+XX*d*d'*XX'; % worst case noise lower bound
AA20=AA20+snrp(ii)*XX*C*hh*XX'+nn;
A30=A30+XX*Hw*Hw'*XX'; % inequality lower bound
aa300=aa300+1/T*log2(det(eye(Nt*Nr)+snr(ii)*XX'*inv(nn)*XX*C*hh));
end
A_ray=1/NN*A_ray; % rayleigh fading with perfect csi and without training symbol
capa_ray(ii)=capa_ray(ii)+1/T(ii)*log2(det(eye(T*Nr)+snr*A_ray*inv(nn_ray)));
A10=1/NN*A10; % perfect csi
capa1(ii)=capa1(ii)+1/T*log2(det(eye(Td*Nr)+snr(ii)*A10*inv(nn)));
A20=1/NN*A20; % worst case noise lower bound
AA20=1/NN*AA20;
capa2(ii)=capa2(ii)+1/T*log2(det(eye(Td*Nr)+snr(ii)*A20*inv(AA20)));
A30=1/NN*A30; % inequality lower bound
capa30(ii)=capa30(ii)+1/T*log2(det(eye(Td*Nr)+snr(ii)*A30*inv(nn)));
aa300=1/NN*aa300;
capa300(ii)=aa300+capa300(ii);
end
capa_ray(ii)=1/MM*capa_ray(ii); % rayleigh fading with perfect csi and without training symbol
capa1(ii)=1/MM*capa1(ii);
capa2(ii)=1/MM*capa2(ii);
capa3(ii)=1/MM*(capa30(ii)-capa300(ii));
ww= waitbar(ii/length(SNR));
end
close(ww)
% figure
plot(SNR,capa_ray,'ko-')
hold on
plot(SNR,capa1,'ks-')
hold on
plot(SNR,capa2,'kv-')
hold on
plot(SNR,capa3,'k*-')
hold on
grid on
title('The Ergodic Capacity of MIMO Systems')
xlabel('SNR (dB)')
ylabel('Ergodic Capacity (bit/channel use)')
legend('Rayleigh','perfect csi','worst case','inequality')
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -