BPSK Demodulation

I found out the reason behind that error... NCO_model produces the same sample value for a single sample of phase input, I need to pass a vector instead.

Code now is as follows:

load LPF.mat
N = length(rf_signal);
phi_inc = zeros(1, N);
bb_sig_I = zeros(1, N);
bb_sig_Q = zeros(1, N);
I_f = zeros(1, N);
Q_f = zeros(1, N);
error_sig = zeros(1, N);
err = zeros(1, N);
% NCO
fs = 1.0E8;
f = 5.0E6;
phi_inc = ((2^32)*f/fs).*ones(1, N);
for i = 1:N
     = NCO_model(phi_inc);
    % Downconverting to Baseband
    bb_sig_I(i) = rf_signal(i)*cos_out(i);
    bb_sig_Q(i) = rf_signal(i)*sin_out(i);
    % Filtering
    I_f(i) = filter(LPF, bb_sig_I(i));
    I_f(i) = filter(LPF, I_f(i));
    Q_f(i) = filter(LPF, bb_sig_Q(i));
    Q_f(i) = filter(LPF, Q_f(i));
    % Error
    error_sig(i) = I_f(i).*Q_f(i);
    % Loop Filter (Leaky Integrator)
    alpha = 0.1;
    err(i) = filter(alpha, , error_sig(i));
    phi_inc = phi_inc - err(i);
end

phi_inc is a vector now, passed to NCO_model on each iteration, and produces new vectors for sin & cos. I only multiply with the index i of both sin_out & cos_out.

However, as you stated... filtering is not correct.