ctlmmse.m

分别给出单天线和多天线下的MIMO－仿真

clc
clear all
Nt = 1; % the number of transmitted antennas
Nr = 1; % the number of transmitted antennas

T=2:2:20; % the number of symbol periods
Tp=Nt;
Td=T-Tp;

SNR = 42; % 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);       % noise covariance
% nn=eye(Td*Nr);       % noise covariance

capa_ray = zeros(1,length(T));  % rayleigh fading with perfect csi and without training symbol
capa1 = zeros(1,length(T));  % perfect csi
capa2 = zeros(1,length(T));    % worst case noise lower bound
capa3 = zeros(1,length(T));    % inequality lower bound
capa30=zeros(1,length(T));
capa300=zeros(1,length(T));

for ii = 1:length(T)
waitbar(ii/length(T));

A_ray=zeros(T(ii)*Nr,T(ii)*Nr);   % rayleigh fading with perfect csi and without training symbol

% capa1 = zeros(1,length(T));  % perfect csi
A10=zeros(Td(ii)*Nr,Td(ii)*Nr);

%capa2 = zeros(1,length(T));    % worst case noise lower bound
A20=zeros(Td(ii)*Nr,Td(ii)*Nr);
AA20=zeros(Td(ii)*Nr,Td(ii)*Nr);


% capa3 = zeros(1,length(T));    % inequality lower bound
% capa30=zeros(1,length(T));
% capa300=zeros(1,length(T));
A30=zeros(Td(ii)*Nr,Td(ii)*Nr);
AA30=zeros(Td(ii)*Nr,Td(ii)*Nr);
aa300=0;

      % LMMSE estimator
p=eye(Nt*Nr);       % training symbol
pp=inv(p'*p)*p';    % wei ni 
F=sqrt(snrp)*hh*p'*inv(snrp*p*hh*p'+vv');
f0=(eye(size(F*pp))-sqrt(snrp)*F*pp)*hm;
nn=eye(Td(ii)*Nr);       % noise covariance
nn_ray=eye(T(ii)*Nr);     % rayleigh fading with perfect csi and without training symbol

    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)*p*Hw+v;                                        % receive
    g=F*z+f0;                                                   % estimation

    % B=zeros(Nt*Nr);
    B=eye(Nt*Nr);
    C=inv(eye(Nt*Nr)+snrp*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(ii)) + j* randn(Nt,T(ii)));   % rayleigh fading with perfect csi and without training symbol
             XX_ray=kron(X_ray.',eye(Nr)); 
    
             X =sqrt(1/(2*Nt))* (randn(Nt,Td(ii)) + j* randn(Nt,Td(ii))); 
             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+snr*XX*C*hh*XX'+nn;
    
             A30=A30+XX*Hw*Hw'*XX';  % inequality lower bound
             aa300=aa300+1/T(ii)*log2(det(eye(Nt*Nr)+snr*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(ii)*Nr)+snr*A_ray*inv(nn_ray)));

    A10=1/NN*A10;     % perfect csi
    capa1(ii)=capa1(ii)+1/T(ii)*log2(det(eye(Td(ii)*Nr)+snr*A10*inv(nn)));

    A20=1/NN*A20;      % worst case noise lower bound
    AA20=1/NN*AA20;
    capa2(ii)=capa2(ii)+1/T(ii)*log2(det(eye(Td(ii)*Nr)+snr*A20*inv(AA20)));

    A30=1/NN*A30;        % inequality lower bound
    capa30(ii)=capa30(ii)+1/T(ii)*log2(det(eye(Td(ii)*Nr)+snr*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(T));
    end
    close(ww)
% figure
plot(T,capa_ray,'ko-')
hold on
plot(T,capa1,'ks-')
hold on
plot(T,capa2,'kv-')
hold on
plot(T,capa3,'k*-')
hold on
grid on
title('The Ergodic Capacity of MIMO Systems')
xlabel('T [channel use]')
ylabel('Ergodic Capacity [bit/channel use]') 

legend('Rayleigh','perfect csi','worst case','inequality')