16-QAM...problem in code

qam problem ber

hi

i am facing problem in showing BER curve.... my code work correctly upto bit recieved but facing problem in calculating errors due to which my BER curve appears to be straight line while it should be decreasing with growing SNR... also my code takes above 1 hr to compute the result .... code is given below
kindly suggest me where i am wrong and help me out

thanks


clear all;
N=10^4 %symbols
M=16; % for 16 qam
error=0;
count=1;
SNR1=[0:20];
a=SNR1/10;
SNR=10.^a;
% EbN0=20;
for i=1:length(SNR)
EbN0=SNR(i);

L=log2(M);% number of bits
bits=1000;
b=0;
while (a<bits)
b=1000*count;
count=count+1;
data=((sign(randn(N*L,1)+1))>0.5).'; % bit to send


tbl=[0;1];
for i=1L-1) tbl=[tbl zeros(length(tbl),1); flipud(tbl) ones(length(tbl),1)]; end;

tmp = sortrows([bi2de(tbl)'; [0 : M-1]]');
Reordering = tmp,2)';

Vals = [0:M-1];

Y= fix(Vals./sqrt(M)); % y=0:3 => row
% % counting forward for odd rows, backward for even rows
x = -sign(mod(Y,2)-0.5).*(bitand(Vals,sqrt(M)-1)*2+1-sqrt(M)); % x = +/- (0:3)*2+1-4
y = sqrt(M)-Y*2-1; % y = 4 - y*2-1
constellation = x + j*y; % constellation points % construct result
% scatterplot(constellation);

constellationpower = sqrt(mean(abs(constellation).^2)); % calculate average power
modulationconstellation= constellation(Reordering+1)/constellationpower;



EsN0 = EbN0 * L;


% needed to perform bi2de conversion
bpwr = (2.^(0L-1)));
% scale looked up modulated constellation index obtained by bi2de conversion

if ((mod(size(data,2),L)~=0) & (1))
data = [data zeros(1,L-mod(size(data,2),L))];
end;

z = sqrt(EsN0) * modulationconstellation( bpwr*reshape(data,L,size(data,2)/L) + 1 );
% noise
noise=1/sqrt(2)*(randn(1,length(z))+j*randn(1,length(z)));
h=1;
y=h*z+noise;%modulation

% demodulation
% getting index of a minimum value
for k=1:length
sig=y(k)*ones(M,1);
[m,index(k)]=min(sig-(h*sqrt(EsN0)*modulationconstellation.'));
end


% need to perform decimal to binary conversion

binary=(de2bi(index-1)).';
bitrcv=binary); % bit recieved
data=data);



% error=error+e;
end
%actually facing problem in this part
e(i) =length(find([data- bitrcv])); % finding errors...
end

BER=e/b;
close all;
figure

semilogy(SNR,BER,'mx-','LineWidth',2);
xlabel('SNR');
ylabel('Bit Error Rate');
title('Bit error probability curve for 16-QAM modulation');

16 qam

oh so sad ...there is no one to help me