classdef BioSigKit < handle
%% ============== BioSigKit ==================
% BioSigKit is a set of useful signal processing tools that are either
% developed by me personally or others in different fields of biosignal
% processing. BioSigKit is a wrapper with a simple visual interface that
% gathers this tools under a simple easy to use platform. This tool might
% be only used for non-commercial, academic, research and learning
% purposes.
%
%% ============== How to start ===================
% example:: obj = BioSigKit(Sig, Fs, Gr)
% Where:::
% Sig: is the signal
% Fs : Sample Rate
% Gr : Flag for showing the interface or not
% Then you can call any subroutine
%% ====List Of Subroutines that you can call for QRS detection =========
% ------- Algorithm -------------------- How to Call ---------------
% (1) Pan-Tompkins Algorithm : obj.PanTompkins
% (2) Phase Space Reconstruction : obj.PhaseSpaceAlg
% (3) RST State-Machine : obj.StateMachine
% (4) Filter Bank: obj.FilterBankQRS
% (5) MTEO qrstAlg: obj.MTEO_qrstAlg
% (6) AMPD: obj.AMPD_PAlg
%% ==== List of all subroutines for ACC, EMG and etc processing ========
% (7) Activity Detection Hilbert: obj.Env_hilbert
% ---------------- Inputs ------------------ %
% Smooth_window : Length of smoothing window in nr of sample
% threshold_style : Set 0 for Automatic, Set 1 for Manual
% DURATION : The number of samples for signal to be above
% threshold to be considered active
% ---------------- Output ------------------ %
% alarm : Pattern of activities
% ---------------- Demo -------------------- %
% v = repmat([.1*ones(200,1);ones(100,1)],[10 1]); % generate true variance profile
% obj.sig = sqrt(v).*randn(size(v));
% obj.Env_hilbert;
%-----------------
% (8) Comp Mobility and Complexity: obj.ComputeHjorthP
%-----------------
% (9) Posture detection 3 Chann ACC: obj.ACC_Act
% ------------- Inputs ----------------------- %
% obj.sig : 3 axis Accelerometer signal where, each row is an axis
% and each column a sample (e.g. (3,:))
% obj.Fs : Sampling frequency of the Accelerometer
% ------------- Output ----------------------- %
% output : adaptively filtered ACC channels based on activity
% state : activity level (0:steady,1:walking,2:joggin)
% EE : Energy Expenditure over 1 min (or length sig)
% F : Bandpass filter in Hz
% SMA : Signal Magnitude area
%------------------
% (10)PsuedoCorr template matching : obj.TemplateMatch
% ------------------- Inputs --------------------- %
% template : A template in the form of a vector, the length of
% the template should be smaller than the signal.
% lag : a lag in terms of nr of samples to move the template,
% it should be smaller than the length of the template
% ------------------- Outputs -------------------- %
% PsC_s : template matching score in range [0,1].
% best_lag: the lag that gave the highest correlation score.
%-----------------
% (11)ECG derived respiration : obj.EDR_comp
%-----------------
% (12)ACC derived respiration : obj.ADR_comp
%% ===== General Projective, linear and nonlinear filterings ============ %%
% (13)Real-time neural PCA: obj.neural_pca
% ----------------- Inputs ------------------- %
% X : Multi-channel signal, each row represents a channel and
% each column a sample.
% nPCA: Number of PCAs to extract
% nit: Number of iterations to go through the whole signal
% T : Learning rate in range [0,1], default:0.9
% ---------------- Outputs ------------------- %
% EigVec: Eigen vectors
% PC : PCs
% Eigval : Eigenvalues
%-----------------
% (14)Adaptive Filtering: obj.adaptive_filter
% * RLS : Recursive Least Squares Filter
% * ALE : Adaptive Line Enhancer (Delayed Filter)
% * VLALE : Variable Leaky Adaptive Line Enhancer
% * NLMS_ecg : Normalized Least Mean Squares filter for artficat
% removal in ECG based on 3 channel Accelerometer recordings
%-------------------- Inputs ---------------------------%
% type: type of the filter (numeric):
% (1) RLS : Recursive Least Squares Filter
% (2) ALE : Adaptive Line Enhancer (Delayed Filter)
% (3) VLALE : Variable Leaky Adaptive Line Enhancer
% (4) NLMS_ecg : Normalized Least Mean Squares filter for
% artficat removal in ECG based on 3 channel Accelerometer
% recordings
% ref: Reference signal:
% * For RLS filter it is single channel (1*N)
% * For NLMS_ECG filter it should be 3 Channel
% Accelerometer (3*N)
% obj.Sig: input signal to clean (single channel vector)
% order : order of the filter (for VLALE and ALE also delay)
% lambda : learning rate(0<= lambda <=1, usually close to 1)
%-------------------- Outputs --------------------------%
% output: cleaned signal
% error_sig : error signal
%-----------------
% (15)Nonlinear phasespace filtering: obj.nonlinear_phase_filt
%-------------------- Method ---------------------------%
% Employs nonlinear phasespace filter to clean up the signal. This
% method is very strong and even able to extract foetal ecg from
% single channel maternal recordings. Please refer to examples of
% BioSigKit for further details.
%-------------------- Inputs ----------------------------%
% sig : Signal to be analyzed (single channel)
% t : Number of samples for computing delayed phase space (def: 1)
% d : Embedding dimension to consider (def: 50)
% m : dimension of null space (def: 49)
% r : number of nearest neighbors to consider
% (normally a large number def: length(sig)/4)
% --------------------- Output ------------------------ %
% output : Cleaned Signal
% --------------------- example ----------------------- %
% output = projective(foetal_ecg(:,2), 1, round(Fs/5), round(Fs/6.25), 1500);
% (16)Teager-Keiser energy operator: obj.TK_comp
%% ============== Licensce ========================================== %%
% THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
% "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
% LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
% FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
% OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
% SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
% TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
% PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
% LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
% NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
% SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
% Author :
% Hooman Sedghamiz, Feb, 2018
% MSc. Biomedical Engineering, Linkoping University
% Email : Hooman.sedghamiz@gmail.com
%% ============== End Helper ================= %%
%------------ Class Definition of BioSigKit Toolbox -------------- %
properties(GetAccess = 'public', SetAccess = 'private')
panel = []; % Main Visualization Panel
SigView = []; % Signal View Axis
FreqValHolder; % Holds SampFreq from interface
Alg; % Alg type from interface
Status; % Loader
LoadedSig; % Loaded Sig
path_alg = addpath(genpath(fileparts(fileparts(which(mfilename)))));% Path of Algorithm
slashchar = char('/'*isunix + '\'*(~isunix)); % OS dependent /
end
%------------ Interface and private methods -------------%
methods(Access = protected)
BioSigKitPanel(obj); % Setup the interface
end
%------------ Properties for BioSigKit ------------------%
properties(GetAccess = 'public', SetAccess = 'public')
Fs; % Sample Freq
Sig; % Signal
Gr = 0; % Flag for showing the interface
Results; % Struct holding results
statsC; % Holding Stats
PlotResult = 0; % Verbose results of each algorithm
PhasePeriod = 0.020; % Phase Period Used for Phase Space
ScalogramL; % Scalogram Length
complexity; % Complexity based on Hjorth Param
mobility; % mobility based on Hjorth
ACC; % 3 channel ACC
end
%% ------------------ Methods ----------------------------- %%
methods
%% --------------------- Constructor -------------------------- %%
function obj = BioSigKit(Sig, Fs, Gr)
% --------- Check Inputs -------------%
if nargin < 3
Gr = 0;
if nargin < 2
Fs = 250;
if nargin < 1
BioSigKitPanel(obj);
return;
end
Gr = 1;
end
end
%--------- Pack Necessary Parameters -------------%
obj.Sig = Sig;
obj.Gr = Gr;
obj.Fs = Fs;
obj.ScalogramL = ceil(length(Sig)/2) - 1;
%--------- Create the interface of Flag is on ------------- %
if obj.Gr
%------ Initialize Figure------%
BioSigKitPanel(obj);
%------ Update Figure -----------%
UpdateFig(obj);
end
end
%% ----------------------Run Algorithms method-------------------- %%
function RunAlg(obj)
if isempty(obj.Sig)
msgbox('First Load a Signal!');
return;
end
[~,ind] = min(size(obj.Sig));
if ind == 1
obj.Sig = obj.Sig(1,:);
else
obj.Sig = obj.Sig(:,1);
end
% ------------------ Set the Values --------------------------%
obj.Fs = str2double(obj.FreqValHolder.String);
% ------------------ Start Loader ----------------------------%
obj.Status.start;
obj.Status.setBusyText('Processing...');
pause(0.01);
UpdateFig(obj);
obj.Results = [];
ResetStats(obj);
% --------------- Run the appropriate method ------------------%
try
switch obj.Alg.Value
% ----- Pan-Tompkins -------%
case 1
PanTompkins(obj);
% ----- Phase Space --------%
case 2
PhaseSpaceAlg(obj);
% ----- RST State-Machine ------%
case 3
StateMachine(obj);
% ----- QRS Filter-Bank -------- %
case 4
FilterBankQRS(obj);
% ---------- QRST MTEO ------------ %
case 5
MTEO_qrstAlg(obj);
%----------- AMPD Peak Detector --------%
case 6
AMPD_PAlg(obj);
otherwise
return;
end
% -------------------- Update Figure --------------------- %
visualizeResults(obj);
catch ME
msgbox(ME.message);
end
% ------------------- Update Loader ---------------------- %
obj.Status.setBusyText('Done!');
pause(1);
obj.Status.stop;
end
%% ---------------------- Import Sig ------------------------------ %%
function ImportSig(obj)
[FileName,Path] = uigetfile(fullfile(pwd,...
'*.mat'),...
'Select Your Signal');
if ~FileName
return;
end
data = load([Path,obj.slashchar,FileName]);
FN = fields(data);
if ~isempty(FN)
obj.Sig = data.(FN{1});
else
msgbox('Bad Input Signal! Input should be a vector!');
return;
end
if ~isempty(obj.LoadedSig)
obj.LoadedSig.String = [Path,FileName];
% ----------------- Update Figure ------------------ %
UpdateFig(obj);
obj.Results = [];
end
end
%% ================= Phase Space ===================== %%
%----- Employs Phase Space for QRS detection ------ %
% Call : obj.PhaseSpaceAlg OR PhaseSpaceAlg(obj)
% Outputs :
% R : R Beats,
% V : Processed area computed from Phase Space
function [R,V] = PhaseSpaceAlg(obj)
obj.Results = [];
[R,V] = PhaseSpaceQRS(obj.Sig,obj.Fs,obj.PhasePeriod,...
obj.PlotResult);
obj.Results.R(1,:) = R(:,1);
end
%% ================= MTEO QRST Delineation ============== %%
%------ Emlpoys MultiValued MTEO to delineate ECG---------- %
% Call : obj.MTEO_qrstAlg OR MTEO_qrstAlg(obj)
% Outputs: Delineated Waves (QRSTP)
function [R,Q,S,T,P_w] = MTEO_qrstAlg(obj)
obj.Results = [];
[R,Q,S,T,P_w] = MTEO_qrst(obj.Sig,obj.Fs,obj.PlotResult);
obj.Results.R(1,:) = R(:,1);
obj.Results.Q(1,:) = Q(:,1);
obj.Results.S(1,:) = S(:,1);
obj.Results.T(1,:) = T(:,1);
obj.Results.P(1,:) = P_w(:,1);
end
%% =================== AMPD Peak Detection ================== %%
% ------------- Employs AMPD Algorithm for Peak detection -----------%
% Inputs :
% Sig
% L : Max sample number for Scalogram computation(Leave free if have not idea)
function R = AMPD_PAlg(obj)
obj.Results = [];
if isempty(obj.ScalogramL)
obj.ScalogramL = ceil(length(obj.Sig)/2) - 1;
end
try
R = AMPD_P(obj.Sig,obj.ScalogramL,obj.PlotResult);
obj.Results.R(1,:) = R(:,1);
catch ME
msgbox(ME.message);
end
end
%% ================== Filter Bank QRS Detector ================== %%
function R = FilterBankQRS(obj)
obj.Results = [];
R = nqrsdetect(obj.Sig,obj.Fs);
obj.Results.R(1,:) = R(:);
end
%% ==================== Pan-Tompkins QRS Detector ================ %%
function [R_amp,R_ind] = PanTompkins(obj)
obj.Results = [];
[R_amp,R_ind]=pan_tompkin(obj.Sig,obj.Fs,obj.PlotResult);
obj.Results.R(1,:) = R_ind(:);
end
%% ===================== Simple State-Machine ===================== %%
% Outputs : R_i : Index of R peaks, R_a : Amplitude
function [R_i,R_a,S_i,S_a,T_i,T_a,Q_i,Q_a] = StateMachine(obj)
[R_i,R_a,S_i,S_a,T_i,T_a,Q_i,Q_a] = SimpleRST(obj.Sig,obj.Fs,...
obj.PlotResult);
obj.Results.R(1,:) = R_i;
obj.Results.Q(1,:) = Q_i;
obj.Results.S(1,:) = S_i;
obj.Results.T(1,:) = T_i;
end
%% ===============Update the plot in interface ========== %%
function visualizeResults(obj)
% --------------- First the Figure ------------------%
if ~isempty(obj.Results)
PP = fields(obj.Results);
Index = [5,10,15;3,8,13;2,7,12;1,6,11;4,9,14];
colors = distinguishable_colors(length(PP),{'k','y'});
for i= 1: length(PP)
% ----------- Update Signal ------------------ %
line(repmat(obj.Results.(PP{i})(1,:),[2 1]),...
repmat([min(obj.Sig-mean(obj.Sig))/2; max(obj.Sig-mean(obj.Sig))/2],...
size(obj.Results.(PP{i})(1,:))),'LineWidth',1.5,...
'LineStyle','-.','color',colors(i,:),'Parent',obj.SigView);
% ----------- Update Stats -------------------- %
intervals = round((diff(obj.Results.(PP{i})(1,:)))./obj.Fs,3);
obj.statsC.Children(Index(i,1)).String = mat2str(max(intervals));
obj.statsC.Children(Index(i,2)).String = mat2str(round(mean(intervals),3));
obj.statsC.Children(Index(i,3)).String = mat2str(length(obj.Results.(PP{i})(1,:)));
end
end
end
%% ================= Update VisualInterface Axis ================= %%
function UpdateFig(obj)
cla(obj.SigView);
plot(obj.Sig - mean(obj.Sig),'LineWidth',2,'Parent',obj.SigView,'color','yellow');
set(obj.SigView,'XGrid','on',...
'YGrid','on','XMinorGrid','on','YMinorGrid','on',...
'Color',[0,0,0],'YColor',[1,1,1],'XColor',[1,1,1]);
axis(obj.SigView,'tight');
end
%% ================ Update Stats ================================== %%
function ResetStats(obj)
for j = 1: 15
obj.statsC.Children(j).String = '\bf --';
end
end
end
%% ============ Experimental Methods not in GUI Yet ================== %%
methods
%% ------ Use Hilber Tr to detect Activity in signals --------- %%
function alarm = Env_hilbert(obj,Smooth_window,threshold_style,...
DURATION)
% ---------------- Inputs ------------------ %
% Smooth_window : Length of smoothing window in nr of sample
% threshold_style : Set 0 for Automatic, Set 1 for Manual
% DURATION : The number of samples for signal to be above
% threshold to be considered active
% ---------------- Output ------------------ %
% alarm : Pattern of activities
% ---------------- Demo -------------------- %
% v = repmat([.1*ones(200,1);ones(100,1)],[10 1]); % generate true variance profile
% obj.sig = sqrt(v).*randn(size(v));
% obj.Env_hilbert;
% ---------------- Input Handling ------------------%
if nargin < 4
DURATION = 20; % default
if nargin < 3
threshold_style = 1; % default 1 , means it is done automatic
if nargin < 2
Smooth_window = 20; % default for smoothing length
end
end
end
alarm = envelop_hilbert(obj.Sig,...
Smooth_window,...
threshold_style,DURATION,obj.PlotResult);
end
%% ------ Employ Hjorth Parameters to analyze a signal --------- %%
function [mobility,complexity] = ComputeHjorthP(obj)
% ---------------- Hjorth Parameters ------------------- %
% ---------------- Output ----------------------%
% mobility : Central frequency
% complexity : bandwith of the signal
[mobility,complexity] = simpleHjorth(obj.Sig,obj.Fs);
obj.mobility = mobility;
obj.complexity = complexity;
end
%% ------ Activity Detection in ACC signals ---------------- %%
function [output,state,EE,F,SMA] = ACC_Act(obj)
% ------------- Inputs ----------------------- %
% obj.sig : 3 axis Accelerometer signal where, each row is an axis
% and each column a sample (e.g. (3,:))
% obj.Fs : Sampling frequency of the Accelerometer
% ------------- Output ----------------------- %
% output : adaptively filtered ACC channels based on activity
% state : activity level (0:steady,1:walking,2:joggin)
% EE : Energy Expenditure over 1 min (or length sig)
% F : Bandpass filter in Hz
% SMA : Signal Magnitude area
if isempty(obj.Sig)
error('Please first load your ACC to obj.ACC!');
end
[output,state,EE,F,SMA] = ACC_Activity(obj.Sig,obj.Fs);
end
%% ------------- Template Matching with Psuedo Corr--------------------- %%
function [PsC_s,best_lag] = TemplateMatch(obj,template,lag)
% ------------------- Inputs --------------------- %
% template : A template in the form of a vector, the length of
% the template should be smaller than the signal.
% lag : a lag in terms of nr of samples to move the template,
% it should be smaller than the length of the template
% ------------------- Outputs -------------------- %
% PsC_s : template matching score in range [0,1].
% best_lag: the lag that gave the highest correlation score.
if nargin < 3
[PsC_s,best_lag] = PsC(template,obj.Sig);
else
[PsC_s,best_lag] = PsC(template,obj.Sig,lag);
end
end
%% ---------------- Real time Neural PCA Filter ----------------- %%
function [EigVec,PC,Eigval] = neural_pca(obj,nPCA,nit)
% ----------------- Inputs ------------------- %
% X : Multi-channel signal, each row represents a channel and
% each column a sample.
% nPCA: Number of PCAs to extract
% nit: Number of iterations to go through the whole signal
% T : Learning rate in range [0,1], default:0.9
% ---------------- Outputs ------------------- %
% EigVec: Eigen vectors
% PC : PCs
% Eigval : Eigenvalues
% ---------------- Check Inputs -------------- %
if size(obj.Sig,1) < 2
fprintf('|#| BioSigKit> Input signals should be channles*N.\n');
end
[EigVec,PC,Eigval] = RTpca(obj.Sig,nPCA,nit);
end
%% ----------------- Adaptive Filters --------------------- %%
function [output,error_sig] = adaptive_filter(obj,type,ref,...
order,lamda)
%-------------------- Inputs ---------------------------%
% type: type of the filter (numeric):
% (1) RLS : Recursive Least Squares Filter
% (2) ALE : Adaptive Line Enhancer (Delayed Filter)
% (3) VLALE : Variable Leaky Adaptive Line Enhancer
% (4) NLMS_ecg : Normalized Least Mean Squares filter for
% artficat removal in ECG based on 3 channel Accelerometer
% recordings
% ref: Reference signal:
% * For RLS filter it is single channel (1*N)
% * For NLMS_ECG filter it should be 3 Channel
% Accelerometer (3*N)
% obj.Sig: input signal to clean (single channel vector)
% order : order of the filter (for VLALE and ALE also delay)
% lambda : learning rate(0<= lambda <=1, usually close to 1)
%-------------------- Outputs --------------------------%
% output: cleaned signal
% error_sig : error signal
% ------------------- Process Inputs -------------------%
output= [];
error_sig = [];
if nargin < 2
fprintf('|#| BioSigKit> Please select a filter type! (1-4).\n');
fprintf('|#| BioSigKit> Exiting.');
return;
end
if nargin < 5
lamda = 1;
if nargin < 4
order = 1*obj.Fs;
if nargin < 3
ref= [];
end
end
end
% ------------------- Run the Filters-------------------%
switch type
% ------ Recursive Least Squars (RLS) ------ %
case 1
if isempty(ref) || size(ref,1) > 1
fprintf('|#| BioSigKit> Ref signal should be 1*N!\n');
return;
end
fprintf('|#| BioSigKit> Adaptive RLS Filter...\n');
[output,error_sig] = RLS(ref,obj.Sig,order,lamda);
% ------ Adaptive Line Enhancer (ADLE) ------ %
case 2
fprintf('|#| BioSigKit> Adaptive Delayed Line Enhancer...\n');
[output,error_sig] = ALE_imp(obj.Sig,obj.Fs,order);
% ------ Variable Leaky ALE (VLALE) ------ %
case 3
fprintf('|#| BioSigKit> Variable Leaky ALE Running...\n');
[output,error_sig] = VLALE_imp(obj.Sig,obj.Fs,order);
% ------ ECG artifact removal with 3 channel ACC ------ %
case 4
if isempty(ref) || size(ref,1) < 3 || size(ref,1) > 3
fprintf('|#| BioSigKit> Accelerometer signals should be 3*N!\n');
return;
end
fprintf('|#| BioSigKit> ECG Artifact removal with NLMS...\n');
output = NLMS_ecg(ref,obj.Sig,order,obj.PlotResult);
otherwise
fprintf('|#| BioSigKit> Filter type not recognized!\n');
end
fprintf('|#| Done!\n');
end
%% ----------------- Non-Linear PhaseSpace Filter -------------- %%
function output = nonlinear_phase_filt(obj,t,d,m,r)
%-------------------- Method ---------------------------%
% Employs nonlinear phasespace filter to clean up the signal. This
% method is very strong and even able to extract foetal ecg from
% single channel maternal recordings. Please refer to examples of
% BioSigKit for further details.
%-------------------- Inputs ----------------------------%
% sig : Signal to be analyzed (single channel)
% t : Number of samples for computing delayed phase space (def: 1)
% d : Embedding dimension to consider (def: 50)
% m : dimension of null space (def: 49)
% r : number of nearest neighbors to consider
% (normally a large number def: length(sig)/4)
% --------------------- Output ------------------------ %
% output : Cleaned Signal
% --------------------- example ----------------------- %
% output = projective(foetal_ecg(:,2), 1, round(Fs/5), round(Fs/6.25), 1500);
%---------------------- Input Handling ----------------------%
if nargin < 5
r = round(length(obj.Sig)/4);
if nargin < 4
m = 49;
if nargin < 3
d = 50;
if nargin < 2
t = 1;
end
end
end
end
% ------------------- Run filter ------------------- %
fprintf('|#| BioSigKit> Running Non-linear Filter...\n');
output = projective (obj.Sig, t, d, m, r);
fprintf('|#| Done!\n');
end
%% =================== ECG derived Respiration =============== %%
function EDR = EDR_comp(obj)
% --------------------- Method ------------------------ %
% 1) Use Pan-tompkins for R peak detection
% 2) Extract templates
% 3) Use neural PCA to extract the eigenvectors
% 4) construct EDR signal.
% ---------------- Compute R peaks --------------------- %
obj.PanTompkins();
% --------------- Segment window size ------------------ %
Win = round(0.356*obj.Fs);
if ~mod(Win,2)
L = Win/2;
Win = Win + 1;
else
L = (Win-1)/2;
end
% --------------- Segment the signal --------------------- %
R = zeros(Win,length(obj.Results.R));
for i = 1 : length(obj.Results.R)
if obj.Results.R(i)-L < 1
dummy = Win - length(obj.Sig(1:obj.Results.R(i)+L));
R(:,i) = [ones(dummy,1)*mean(obj.Sig(1:obj.Results.R(i)+L));...
obj.Sig(1:obj.Results.R(i)+L)];
elseif obj.Results.R(i)+L > length(obj.Sig)
dummy = Win - length(obj.Sig(obj.Results.R(i)-L:end));
R(:,i) = [obj.Sig(obj.Results.R(i)-L:end);...
ones(dummy,1)*mean(obj.Sig(obj.Results.R(i)-L:end))];
else
R(:,i) = obj.Sig(obj.Results.R(i)-L:obj.Results.R(i) + L);
end
end
% ------------------- Real Time PCA ------------------------ %
temp_ecg = obj.Sig;
obj.Sig = R';
PC = neural_pca(obj,size(R,2),2);
obj.Sig = temp_ecg;
PC = PC(:,1);
% ------------------- Build EDR Signal --------------------- %
xx = 1 : length(temp_ecg);
interp1 = [0 obj.Results.R length(temp_ecg)];
PC = [0 PC' 0];
EDR = spline(interp1,PC,xx);
EDR = detrend(EDR);
end
%% ================ 3 Channel ACC derived Respiration ========== %%
function ADR = ADR_comp(obj)
%------------ Method -----------------------------%
% (1) Adaptively filter the ACC signals
% (2) Real Time PCA
% (3) Mix PCs and reconstruct the signal
%------------------- Input handling ----------------%
if size(obj.Sig,1) < 3
fprintf('|#| BioSigKit> ACC signal should be 3*n!');
return;
end
output = ACC_Activity(obj.Sig,obj.Fs);
% ------------------- Real Time PCA ------------------------ %
PC = neural_pca(obj,3,2);
PC = PC(:,1);
ADR = output'*PC;
ADR = detrend(ADR);
end
%% =================== Teager-Keiser Energy Op =========== %%
function [ey,ex] = TK_comp(obj)
% --------------- Outputs ----------------%
% ey: energy operator
% ex: Teager operator
[ey,ex]=energyop(obj.Sig,obj.PlotResult);
end
end
end