alamouti_2xnr.m

基于MIMO通信系统中的各种空时编译码算法

% BPSK modulation, 2 X Nr antenna array
clear all;
%% initialztion
Nr=2; % Modify here!
M=2;% M-QAM,can't be modified
data_length=500000;% there are data_length*Nt bits at all
Nt=2; % can't be modified!
EbNo=(0:3:18);
BER=zeros(size(EbNo));
fr_length=data_length/(log2(M));
h=waitbar(0,'Percentage Completed');
set(h,'name','Alamouti STBC 2XNr...');
step_wb=100/(length(EbNo)-1);
wb=0;
%% QAM modulation
data=randint(data_length,Nt);
X=[];
for n=1:Nt
    temp=qammod(data(:,n),M); %M-QAM
    X=[X temp];
end
clear temp;
constellation=qammod((0:M-1),M);
%%  Alamouti STBC
x0=X;
x1(:,1)=-conj(X(:,2));
x1(:,2)=conj(X(:,1));
% x1 and x2 are consecutive symbols
%% channel matrix
for n=1:Nr  % rx_antennas
    H(n,:,:)=(randn(fr_length,Nt)+j*randn(fr_length,Nt))/sqrt(2); %Raylaigh channel
end
%% transmit and receive
SNR=EbNo;
for ii=1:length(SNR)
    for n=1:Nr  %received signal per receiver antenna  
        Hn=reshape(H(n,:,:),fr_length,Nt);
        Habs(:,n)=sum(abs(Hn).^2,2);
        r0(:,n)= sum(Hn.*x0,2)/sqrt(Nt)  +  sqrt(0.5/(10^(SNR(ii)/10)))*(randn(fr_length,1)+j*randn(fr_length,1));% AWGN
        r1(:,n)= sum(Hn.*x1,2)/sqrt(Nt)  +  sqrt(0.5/(10^(SNR(ii)/10)))*(randn(fr_length,1)+j*randn(fr_length,1));% Alamouti's paper (14)
        % combine the received signals
        z0(:,n)= r0(:,n).*conj(Hn(:,1))  +  conj(r1(:,n)).*Hn(:,2);
        z1(:,n)= r0(:,n).*conj(Hn(:,2))  -  conj(r1(:,n)).*Hn(:,1);   % Alamouti's paper (15)
    end
    Xe(:,1)=sum(z0,2);
    Xe(:,2)=sum(z1,2); % esitmation of X
    % Maximum Likelihood Decision
    for p=1:M
        d0(:,p)=(sum(Habs,2)-1)*abs(constellation(p))^2  +  abs(sum(z0,2)-constellation(p)).^2;
        d1(:,p)=(sum(Habs,2)-1)*abs(constellation(p))^2  +  abs(sum(z1,2)-constellation(p)).^2; % Alamouti's paper (7)
    end
    [y0,i0]=min(d0,[],2);
    y0=constellation(i0);
    [y1,i1]=min(d1,[],2);
    y1=constellation(i1);
    Xd=[y0.' y1.'];
    % QAM demodulation
    Xd_bit=[];
    for n=1:Nt
        temp=qamdemod(Xd(:,n),M); %M-QAM 
        Xd_bit=[Xd_bit temp];
    end
    [errnum,BER(ii)]=biterr(data,Xd_bit);
    %update the waitbar
    str_bar=[num2str(wb) '% Completed'];
    waitbar(wb/100,h,str_bar);
    wb=wb+step_wb;
end
close(h);
semilogy(EbNo,BER,'b-v');
xlabel('Eb/No');
ylabel('Bit Error Rate');
grid on;