%function Pb=rayl_theo(M,K,rho_in_db)
M=64;
K=2;
rho_in_db=8;
rho=10^(rho_in_db/10);
t=2^(K-1)/(2^K-1);
part=0;
for i=1:M-1
    comb=factorial(M-1)/(factorial(i)*factorial(M-1-i));
    part1=((-1)^(i+1))/(1+i+i*K*rho);
    part=part+comb*part1;
    a(i)=part;
end
Pb=t*part;