Rayleigh flat fading simulation_ask for _help

I'd like to ask for your advice. You see, I simulated BPSK passband signal,
then I created flat Rayleigh fading and created multipath signal with no
delay of the 4 paths so that the fading would be flat.
But something's wrong with it though I based by simulation on theory.
Because when I use coherent correlator detector the BER I get is of different
nature than in theory for coherent BPSK flat Rayleigh fading case.

Could you have a look at the small piece of code I attached
BPSK_flat_Rayleigh_passband.m and see what might be wrong.

My detector that works perfectly in accordance with theory in AWGN channel
is BPSK_demod_BER.m
I use BPSK passband signal with 100 oversampling.
It's very important for me. If I learn how to properly do this I will
then shift to applying it to the GMSK passband signal which is what I have
to do for writing my diploma thesis.

Thank you in advance.

hi

see this:

clear
clc
N=1000000;
snr=0:2:14;
for m=1:size(snr,2)
efade1 = 0;
efade2 = 0;
egaus = 0;
No=10^(-snr(m)/10); % PSD of AWGN
sigma= sqrt(No/2); % Standard deviation
for n=1:N
s=(-1)^round(rand);
x1=randn/sqrt(2);
y1=randn/sqrt(2);
x2=randn/sqrt(2);
y2=randn/sqrt(2);
p1=x1^2 + y1^2; % Fading Power
p2=x2^2 + y2^2; % Fading Power
a1=sqrt(p1); % Fading Amplitude
a2=sqrt(p2); % Fading Amplitude
r1=a1*s + sigma*(randn); % faded signal + AWGN
r2=a2*s + sigma*(randn); % faded signal + AWGN

if (r1>0)
shat1=+1;
else
shat1=-1;
end

if abs(r1)>abs(r2)
rsc=r1;
else
rsc=r2;
end

if (rsc>0)
shat2=+1;
else
shat2=-1;
end

if (shat1 ~= s)
efade1 = efade1 + 1;
end

if (shat2 ~= s)
efade2 = efade2 + 1;
end

end
pe1(m)=efade1/N; % BER flat Rayleigh fading
pe2(m)=efade2/N; % BER flat Rayleigh fading

end
semilogy(snr,pe1,'ko-')
hold on
semilogy(snr,pe2,'bs-')

SNR=0:10;
BERtheory = berawgn(SNR,'psk',2,'non-diff');

semilogy(SNR,BERtheory,'k.--')
legend('Rayleigh channel','selection combinig','AWGN');
ylabel('Bit Error Rate');
xlabel('Average E_b/N_0 (dB)');
title('Perfomance of BPSK over flat-Rayleigh and AWGN channels');
axis([ snr(1) snr(end) 9.5e-6 1])
grid on;


plz press Helped me.