QPSK modulation for cooperative communications

hi,
i had done plot graph BER cooperative communications with using QAF (Quantize-Amplify and forward) protocol..i use bpsk modulation and 4 level quantization..please tell me how to use qpsk modulation with 16 level quantization..the matlab code as below..actually i dont know this code right or not, so please correct and help me..thank you


close all

sum_aqf = 0;
N_bits = 1000; %Number of data bits
BERthreshold=1e-3;
N_iter = 100;%number of iterations

for iter = 1:N_iter

SNRdB = 0:30;%Range SNR

for k = 1:length(SNRdB)

SNR = 10^(SNRdB(k)/10); %convert SNRdB to linear value SNR

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

%channel coding using 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 channel response for the relay channel:
%Source to Destination
noise_d = 1/sqrt(2) * (randn(2 * N_bits,1) + j * randn(2 * N_bits,1));
h_d = 1/sqrt(2) * (randn(2* N_bits,1) + j * randn(2 * N_bits,1));

%Source to Relay
noise_r1 = 1/sqrt(2) * (randn(2 * N_bits,1) + 1i * randn(2 * N_bits,1));
h_r1 = 1/sqrt(2) * (randn(2 * N_bits,1) + 1i * randn(2 * N_bits,1));

%Relay to Destination
noise_r2 = 1/sqrt(2) * (randn(2 * N_bits,1) + j * randn(2 * N_bits,1));
h_r2 = 1/sqrt(2) * (randn(2 * N_bits,1) + j * randn(2 * N_bits,1));

ftx_r1 = sqrt(SNR) * h_r1 .* tx + noise_r1;

%%Relay %%
%%Apmlify_Quantize and Forwar%%
beta = sqrt(1./((SNR * abs(h_r1).^2) + 1));
adc_resolution = 2; % Or Quantizer resolution
no_quantiz_levels = power(2, adc_resolution);

%Quantize%
ftx_quant = 2*(floor((no_quantiz_levels-1)*ftx_r1)/(no_quantiz_levels/2))-1;

%amplify:
ftx_ampq = ftx_quant .* beta;

%Amplify_Quantize and Forward
ftx_aqf_r2 = sqrt(SNR) * ftx_ampq .* h_r2 + noise_r2 ;

%Direct Link
ftx_d = sqrt(SNR)* tx .* h_d + noise_d;

%%Destination%%
%%Quantize and Forward%%%
%MRC%%
R_aqf = ftx_aqf_r2 .* conj(h_r2) .* conj(h_r1) + ftx_d .* conj(h_d);
dec_com_aqf = sign(real(R_aqf));

%BER%
%at destination:
err_com1(k) = sum(abs(dec_com_aqf - tx)/2);
simber_com1(k) = err_com1(k) / (2 * N_bits);

end
sum_aqf = sum_aqf + simber_com1;
if (sum_aqf > BERthreshold) sum_aqf = sum_aqf+1;
else sum_aqf = sum_aqf; end
end

avgber_aqf = sum_aqf/N_iter;

figure
semilogy(SNRdB,avgber_aqf,'go-','LineWidth',2);
axis([0 30 10^-5 0.5]);
legend('AQF');
xlabel('SNR dB');
ylabel('BER');
title(['Curve of BER']);
grid on