function [fqrs, ch_ref, X, Xcl] = MyFBeatDect(ECG, Fs, bs, mqrs, bPlot)
%
% function: computes beat detection
%
% IN:
% ECG: raw data. N x 4 matrix
% Fs: sampling rate
% bs: annotated fqrs series (for debug only)
%
% OUT:
% fqrs: foetal QRS annotation. Lq x 1 matrix
% ch_ref: reference channel (used for fqrs)
% X: feature vector
% Xcl: classif outcome vector (set A)
%
fqrs = [];
N = 48+1;
b = ones(N, 1)/N; a = 1;
for k=1:size(ECG,2)
D = diff(ECG(:,k)).^2;
Sk = filtfilt(b, a, D);
S(:,k) = Sk/std(Sk);
[Tq, std_rr(k) rr_iqr(k) rr_med(k) ampl(k)] = MyBuildSeq(S(:,k), Fs, bs);
end
[~, ch_ref] = min(std_rr);
fqrs = MyBuildSeq(S(:,ch_ref), Fs, bs);
for k=1:size(ECG,2)
if( bPlot )
figure(1)
set(gcf, 'Name', sprintf('rECG - QRS CH: %d', ch_ref));
t = (1:size(ECG,1))/Fs;
subplot(2,2,k)
plot(t, ECG(:,k))
hold on
plot(t(bs), ECG(bs,k)+1, '.r'),
if( isempty(fqrs)==0 )
plot(t(fqrs), ECG(fqrs,k), '.k'),
end
set(gca, 'TickLength', [0 0])
xlim([6 10])
hold off
figure(2)
subplot(2,2,k)
plot(t(1:end-1), S(:,k))
hold on
plot(t(bs), S(bs,k), '.r'),
if( isempty(fqrs)==0 )
plot(t(fqrs), S(fqrs,k), '.k'),
end
set(gca, 'TickLength', [0 0])
xlim([6 10])
hold off
end
end