%% (Internal) Delayed covariance matrix calculation
% Function for calculating covariance matrix in the vicinity of
% peaks with k_hb heartbeats of delay
%
% Example
%
% A = DelayedCovMat(x, peaks, win_size, nsamp, k_hb)
%
% Arguments:
% +x: signal matrix
%
% +peaks: QRS complexes location
%
% +win_size: HALF of the size of the win around each QRS complex to
% compute te cov matrix.
%
% +nsamp: size(x,1)
%
% +k_hb: delay in number of heartbeats to calculate the cov mat.
%
% Output:
% + A: The computed cov mat.
%
% See also PiCA
%
% Author: Mariano Llamedo Soria llamedom@electron.frba.utn.edu.ar
% Version: 0.1 beta
% Last update: 14/5/2014
% Birthdate : 23/4/2013
% Copyright 2008-2015
function A = DelayedCovMat(x, peaks, win_size, nsamp, k_hb)
peaks = peaks( (peaks - win_size) > 0 & (peaks + win_size ) <= nsamp );
aux_seq_idx = 1:length(peaks);
T0 = cell2mat(arrayfun(@(a)( peaks(a) - win_size:peaks(a) + win_size ), aux_seq_idx(1:end-k_hb), 'UniformOutput', false));
T1 = cell2mat(arrayfun(@(a)( peaks(a) - win_size:peaks(a) + win_size ), aux_seq_idx(1+k_hb:end), 'UniformOutput', false));
x_T0_centered = bsxfun(@minus,x(T0,:),mean(x(T0,:))); % Remove mean
x_T1_centered = bsxfun(@minus,x(T1,:),mean(x(T1,:))); % Remove mean
A = 1/length(T0)*(x_T0_centered'*x_T1_centered);
A = (A+A')/2;