% nTx = 2; nRx = 2;
clear all;
close all;
NFFT = 32;%1024; FFT size
EbNo = 0:5:25; %0:5:20;
M = 2;
Samp = 1e3; %1000;
nTx = 2;
nRx = 2;
for ii = 1:length(EbNo)
err = 0;
for value = 1:Samp
%*******Generate QPSK data***************
%ip = 2*(rand(nTx, NFFT)>0.5)-1 + j*(2*(rand(nTx, NFFT)>0.5)-1); % 0-->-1, 1-->1
%plot(real(ip),imag(ip),'o');
ip = rand(nTx, NFFT)>0.5; % BPSK
%*******Generate QPSK data***************
%sMod = ip;
sMod = 2*ip-1; %BPSK 0-->-1, 1-->1
%*******Generate square matrix H***************
CIR = 1/sqrt(2)*(randn(nRx, nTx, NFFT)+j*randn(nRx, nTx, NFFT));
%*******Noise power***************
%Noise generation
npow = 10^(-EbNo(ii)/10);
attn = 10^(-EbNo(ii)/20);
ipHat = zeros(nTx, NFFT);
for iii = 1:NFFT
n = (1/sqrt(2))*(randn(nRx, 1)+j*randn(nRx, 1));
%*******Receiver***************
y = CIR(:,:,iii)*sMod(:,iii)+attn*n;
%*******Equalizer ML: J = |y-Hx_estimated|***************
% Maximum Likelihood Receiver
% if [s1 s2 ] = [+1;+1 ] Demap: [1;1]
s_1_1 = [1;1];
J1 = sum(abs(y-CIR(:,:,iii)*s_1_1));
% if [s1 s2 ] = [+1;-1 ] Demap: [1;0]
s_1_1 = [1;-1];
J2 = sum(abs(y-CIR(:,:,iii)*s_1_1));
% if [s1 s2 ] = [-1;+1 ] Demap: [0;1]
s_1_1 = [-1;1 ];
J3 = sum(abs(y-CIR(:,:,iii)*s_1_1));
% if [s1 s2 ] = [-1;-1 ] Demap: [0;0]
s_1_1 = [-1;-1 ];
J4 = sum(abs(y-CIR(:,:,iii)*s_1_1));
% finding the minimum from the four alphabet combinations
mid = [J1,J2,J3,J4];
[value index] = min(mid); % min(A,[],B): if A is a matrix, B=1 means min A in the row and B = 2,...in the column
%*******Demodulation***************
% mapping the minima to bits
ref = [1 1 0 0; 1 0 1 0 ];
ipHat(:, iii) = ref(:, index);
end
%*******Calculate errors***************
err = err + size(find([reshape(ip,1,[])- reshape(ipHat,1,[])]),2); % counting the number of errors
end
nErr(ii) = err;
end
sErr = nErr/NFFT/nTx/Samp;
ebno = 10.^(EbNo/10);
theory_QPSK = 0.5*(1-sqrt(ebno./(ebno+1)));
%theoryBer_nRx1 = 0.5.*(1-1*(1+1./ebno).^(-0.5));
p = 1/2 - 1/2*(1+1./ebno).^(-1/2);
theoryBerMRC_nRx2 = p.^2.*(1+2*(1-p));
pAlamouti = 1/2 - 1/2*(1+2./ebno).^(-1/2);
theoryBerAlamouti_nTx2_nRx1 = pAlamouti.^2.*(1+2*(1-pAlamouti));
semilogy(EbNo, theory_QPSK, 'bs-', EbNo, theoryBerMRC_nRx2, 'kd-', EbNo, theoryBerAlamouti_nTx2_nRx1, 'c+-', EbNo, sErr, 'mo-', 'linewidth', 2);
legend('theory (nTx=1,nRx=1)', 'theory (nTx=1,nRx=2, MRC)', 'theory (nTx=2, nRx=1, STBC)', 'sim (nTx=2, nRx=2, STBC)');
axis([0 20 10^-5 1])
grid on
%legend('theory: 1x1 QPSK-MIMO', 'simulation: nTx x nRx QPSK-MIMO');
xlabel('Eb/No, dB')
ylabel('BER')
title('BER curve for BPSK')