Need Help Please! Simulate Quantize-and-Forward Multi-hop Relay (BER vs SNR)

[Isizaki]

Hi guys, i'm trying to simulate Multi-hop Relay cooperative communication using Quantize-and-Forward to find BER vs SNR. The result isn't correct yet, and i'm really confused right now. I really need your help guys.

Here's the Matlab code i've made (just try to run it on your Matlab):

%BER vs SNR for Direct and QF multi-hop relay%
%using BPSK and rayleigh fading channel + AWGN%
%using Convolutional code 1/2%
%using MRC%
%Isyqi Harzaki

clc
clear all
close all

sum_qf = 0;
sum_qf2 = 0;
sum_qf3 = 0;
sum_qf4 = 0;

sum_d = 0;


N_bits = 1000; % bit data
BERthreshold=1e-5;
N_iter = 10;%iteration


for iter = 1:N_iter
SNRdB = 0:30;%range of SNR

for k = 1:length(SNRdB)
SNR = 10^(SNRdB(k)/10); %convert SNRdB to linear SNR


%Source
data = round(rand(N_bits,1));%random data bit

%kanal koding dengan rate 1/2 convolutional code:
trellis = poly2trellis(3,[7 5]); %trellis structure
c_data = convenc(data,trellis); %encoding data

tx = 2*c_data - 1;% BPSK modulation


%additive noise and fading channel at relays channel:
%Source to Destination
noise_sd = 1/sqrt(2) * (randn(2 * N_bits,1) + j * randn(2 *N_bits,1));
h_sd = 1/sqrt(2) .* (randn(2 * N_bits,1) + j * randn(2 * N_bits,1));

%Source to Relay1
noise_sr1 = 1/sqrt(2) * (randn(2 * N_bits,1) + j * randn(2 *N_bits,1));
h_sr1 = 1/sqrt(2) * (randn(2 * N_bits,1) + j * randn(2 * N_bits,1));

%Relay1 to Relay2
noise_r1r2 = 1/sqrt(2) * (randn(2 * N_bits,1) + j * randn(2 *N_bits,1));
h_r1r2 = 1/sqrt(2) * (randn(2 * N_bits,1) + j * randn(2 * N_bits,1));

%Relay2 to Relay3
noise_r2r3 = 1/sqrt(2) * (randn(2 * N_bits,1) + j * randn(2 *N_bits,1));
h_r2r3 = 1/sqrt(2) * (randn(2 * N_bits,1) + j * randn(2 * N_bits,1));

%Relay3 to Destination
noise_r3d = 1/sqrt(2) * (randn(2 * N_bits,1) + j * randn(2 *N_bits,1));
h_r3d = 1/sqrt(2) * (randn(2 * N_bits,1) + j * randn(2 * N_bits,1));


%the signal from Source to Relay1%
tx_sr1 = sqrt(SNR) * tx .* h_sr1 + noise_sr1;

%quantize
L=2^7;
Xmax=1;
Xmin=-1;
delta=(Xmax-Xmin)/L;
ih1=round((tx_sr1-Xmin)/delta);
if (ih1==L)
ih1=ih1-1;
end
if ih1<0;
ih1=0;
end
tx_qf1=Xmin+ih1*delta;


%the signal from Relay 1 to Relay 2%
tx_r1r2 = sqrt(SNR) * tx_qf1 .* h_r1r2 + noise_r1r2;
%Quantize:
ih2=round((tx_r1r2-Xmin)/delta);
if (ih2==L)
ih2=ih2-1;
end
if ih2<0;
ih2=0;
end
tx_qf2=Xmin+ih2*delta;


%the signal from Relay 2 to Relay 3%
tx_r2r3 = sqrt(SNR) * tx_qf2 .* h_r2r3 + noise_r2r3;
%Quantize:
ih3=round((tx_r2r3-Xmin)/delta);
if (ih3==L)
ih3=ih3-1;
end
if ih3<0;
ih3=0;
end
tx_qf3=Xmin+ih3*delta;


%the signal from Relay 3 to Destination%
tx_r3d = sqrt(SNR) * tx_qf3 .* h_r3d + noise_r3d;

%%the signal from Source to Destination%%
tx_sd = sqrt(SNR) * tx .* h_sd + noise_sd;


%Destination
%%Quantize and Forward%%%
%MRC%%
R_qf = tx_qf1 .* conj(h_sr1).* conj(h_r1r2);
dec_com_qf = sign(real(R_qf));

%BER calculations Hop 1%
%at destination:
err_com1(k) = sum(abs(dec_com_qf - tx_qf1)/2);
simber_com1(k) = err_com1(k) / (2 * N_bits);


R_qf2 = tx_qf2 .* conj(h_r1r2).* conj(h_r2r3);
dec_com_qf2 = sign(real(R_qf2));

%BER calculations Hop 2%
%at destination:
err_com2(k) = sum(abs(dec_com_qf2 - tx_qf1)/2);
simber_com2(k) = err_com2(k) / (2 * N_bits);


R_qf3 = tx_qf3 .* conj(h_r2r3).* conj(h_r3d);
dec_com_qf3 = sign(real(tx_qf3));
%BER calculations Hop 3%
%at destination:
err_com3(k) = sum(abs(dec_com_qf3 - tx_r2r3)/2);
simber_com3(k) = err_com3(k) / (2 * N_bits);


R_qf4 = tx_r3d .* conj(h_r3d).* conj(h_sd) + tx_sd .* conj(h_sd);
dec_com_qf4 = sign(real(R_qf4));
%BER calculations Hop 4%
%at destination:
err_com4(k) = sum(abs(dec_com_qf4 - tx_r3d)/2);
simber_com4(k) = err_com3(k) / (2 * N_bits);

%theoratical rayleigh
direct(k) = 0.5 * (1 - sqrt(SNR./(1 + SNR)));


end
sum_qf = sum_qf + simber_com1 ;
sum_qf2 = sum_qf2 + simber_com2 ;
sum_qf3 = sum_qf3 + simber_com3 ;
sum_qf4 = sum_qf4 + simber_com4 ;
sum_d = direct;

if (sum_qf > BERthreshold) sum_qf = sum_qf
else sum_qf = sum_qf;
if (sum_qf2 > BERthreshold) sum_qf2 = sum_qf2+1
else sum_qf2 = sum_qf2;
if (sum_qf3 > BERthreshold) sum_qf3 = sum_qf3+1
else sum_qf3 = sum_qf3;
if (sum_qf4 > BERthreshold) sum_qf4 = sum_qf4+1
else sum_qf4 = sum_qf4;
if (sum_d > BERthreshold) sum_d = sum_d+1;
else sum_d = sum_d;
end
end
end
end
end
end


avgber_qf = sum_qf/N_iter;
avgber_qf2 = sum_qf2/N_iter;
avgber_qf3 =sum_qf3/N_iter;
avgber_qf4 =sum_qf4/N_iter;

avgber_d = sum_d;


semilogy(SNRdB,avgber_d,'s-r','LineWidth',2)
hold on
semilogy(SNRdB,avgber_qf, '-b','LineWidth',2)
hold on
semilogy(SNRdB,avgber_qf2, '-k','LineWidth',2)
hold on
semilogy(SNRdB,avgber_qf3, 'o-k','LineWidth',2)
hold on
semilogy(SNRdB,avgber_qf4, '-y','LineWidth',2)
hold on

axis([0 30 10^-4 10^0]);
xlabel('SNR dB');
ylabel('BER');
legend('Direct','Hop 1','Hop 2','Hop 3','Hop 4')
grid on

Here's the result i've got so far, and still incorrect


Hope you will help me solve the problem, Thanks

(sorry for my bad english, i'm newbie here)