clc;
clear;
close all;
%
% Quaternary Phase Shift Keying (QPSK) in MIMO
%
Iter=100000; % Number of iterations
SNR_dB=0:10:30; % Signal-to-Noise Ratio (SNR) in dB scale
SNR_linear=10.^(SNR_dB/10); % Signal-to-Noise Ratio (SNR) in linear scale
error_ML=zeros(length(SNR_dB),1);
error_ZF=zeros(length(SNR_dB),1);
error_MMSE=zeros(length(SNR_dB),1);
error_ZF_DF=zeros(length(SNR_dB),1);
error_MMSE_DF=zeros(length(SNR_dB),1);
Eb=1;
k=1;
M=2; % No of transmit antennas
N=2; % No of receive antenna
data=[];
for p=-1:2:1
for q=-1:2:1
for s=-1:2:1
for t=-1:2:1
temp=[p+sqrt(-1)*q;s+sqrt(-1)*t];
data=[data temp];
end
end
end
end
for snr=SNR_linear,
N0=Eb/snr; % snr=Eb/N0
for i=1:Iter,
data1=sign(randn(M,1)); % data is +1 or -1
data2=sign(randn(M,1)); % data is +1 or -1
X=sqrt(Eb)*(data1+sqrt(-1)*data2);
H=(randn(N,M)+sqrt(-1)*randn(N,M))*sqrt(0.5);
n=sqrt(N0/2)*(randn(N,1)+sqrt(-1)*randn(N,1));
r=H*X+n; % Received signal
%%%%% ML Detection %%%%%%%
temp = zeros(1,16);
for b = 1 : 16
temp(b) = norm(r-H*data(:,b));
end
[Dist_min index]=min(temp);
detected_signal_ML1(1)=real(data(1,index));
detected_signal_ML1(2)=real(data(2,index));
detected_signal_ML2(1)=imag(data(1,index));
detected_signal_ML2(2)=imag(data(2,index));
%%%%% ZF Detection %%%%%%%
C_ZF = inv(H'*H)*H';
y_ZF = C_ZF*r;
detected_signal_ZF1(1) = sign(real(y_ZF(1,1)));
detected_signal_ZF1(2) = sign(real(y_ZF(2,1)));
detected_signal_ZF2(1) = sign(imag(y_ZF(1,1)));
detected_signal_ZF2(2) = sign(imag(y_ZF(2,1)));
% MMSE DETECTION %
I = eye(2,2);
C_MMSE = inv(H'*H + N0*I)*H';
y_MMSE = C_MMSE*r;
detected_signal_MMSE1(1) = sign(real(y_MMSE(1,1)));
detected_signal_MMSE1(2) = sign(real(y_MMSE(2,1)));
detected_signal_MMSE2(1) = sign(imag(y_MMSE(1,1)));
detected_signal_MMSE2(2) = sign(imag(y_MMSE(2,1)));
%%%%% ZF -DF DETECTION %%%%%%
W1_ZF = (H'*H)\H';
a1_zf = W1_ZF(1,:)*r;
a1_ZF = sign(real(a1_zf)) + sqrt(-1)*sign(imag(a1_zf));% detect
H_1 = H(:,1);
r2_ZF = r - a1_ZF*H_1 ;
H_2 = H(:,2);
h_2 = inv(H_2'*H_2);
W2_ZF = h_2*H_2';
a2_zf = W2_ZF*r2_ZF;
a2_ZF = sign(real(a2_zf)) + sqrt(-1)*sign(imag(a2_zf));% detect
y_zfdf = [a1_ZF ; a2_ZF];
detected_signal_ZF_DF1(1) = sign(real(y_zfdf(1,1)));
detected_signal_ZF_DF1(2) = sign(real(y_zfdf(2,1)));
detected_signal_ZF_DF2(1) = sign(imag(y_zfdf(1,1)));
detected_signal_ZF_DF2(2) = sign(imag(y_zfdf(2,1)));
% MMSE-DF DETECTION %
W1_MMSE = (H'*H + N0*I)\H';
a1_mmse = W1_MMSE(1,:)*r;
a1_MMSE = sign(real(a1_mmse))+sqrt(-1)*sign(imag(a1_mmse));% detect
r2_MMSE = r - a1_MMSE*H_1;
h_2_mmse = inv(H_2'*H_2 + N0);
W2_MMSE = h_2_mmse*H_2';
a2_mmse = W2_MMSE*r2_MMSE;
a2_MMSE = sign(real(a2_mmse))+sqrt(-1)*sign(imag(a2_mmse));% detect
y_mmse_df = [ a1_MMSE ; a2_MMSE];
detected_signal_MMSE_DF1(1) = sign(real(y_mmse_df(1,1)));
detected_signal_MMSE_DF1(2) = sign(real(y_mmse_df(2,1)));
detected_signal_MMSE_DF2(1) = sign(imag(y_mmse_df(1,1)));
detected_signal_MMSE_DF2(2) = sign(imag(y_mmse_df(2,1)));
%%%%% Detection %%%%%%%
for m=1:2,
if detected_signal_ML1(m) ~= data1(m)
error_ML(k)=error_ML(k)+1;
end
if detected_signal_ML2(m) ~= data2(m)
error_ML(k)=error_ML(k)+1;
end
if detected_signal_ZF1(m) ~= data1(m)
error_ZF(k)=error_ZF(k)+1;
end
if detected_signal_ZF2(m) ~= data2(m)
error_ZF(k)=error_ZF(k)+1;
end
if detected_signal_MMSE1(m) ~= data1(m)
error_MMSE(k)=error_MMSE(k)+1;
end
if detected_signal_MMSE2(m) ~= data2(m)
error_MMSE(k)=error_MMSE(k)+1;
end
if detected_signal_ZF_DF1(m) ~= data1(m)
error_ZF_DF(k)=error_ZF_DF(k)+1;
end
if detected_signal_ZF_DF2(m) ~= data2(m)
error_ZF_DF(k)=error_ZF_DF(k)+1;
end
if detected_signal_MMSE_DF1(m) ~= data1(m)
error_MMSE_DF(k)=error_MMSE_DF(k)+1;
end
if detected_signal_MMSE_DF2(m) ~= data2(m)
error_MMSE_DF(k)=error_MMSE_DF(k)+1;
end
end
end
k=k+1;
end
% Divide the number of error by the total number of transmitted bits
error_ML=error_ML/(4*Iter);
error_ZF=error_ZF/(4*Iter);
error_MMSE=error_MMSE/(4*Iter);
error_ZF_DF=error_ZF_DF/(4*Iter);
error_MMSE_DF=error_MMSE_DF/(4*Iter);
figure;
semilogy(SNR_dB, error_ZF,'rs-')
hold on;
semilogy(SNR_dB, error_MMSE,'bs-')
semilogy(SNR_dB, error_ZF_DF,'r*-')
semilogy(SNR_dB, error_MMSE_DF,'b*-')
semilogy(SNR_dB, error_ML,'g^-')
legend('ZF','MMSE','ZF-DF','MMSE-DF','ML');
xlabel('E_b/N_0'); % label of x axis
ylabel('BER'); % label of y axis
grid on;