v1.3.0

See Changelog for v1.3.0.  Added ability to filter out P and T waves when detecting R peaks.  Updated Threshold figure.  Fixed bug that prevented loading .scp files if there were errors in non-ECG data (demographics etc).

AnnoResult.m

@@ -106,7 +106,7 @@
                     obj.missing_lead = num2cell(missing_lead);
 
                     % VERSION MANUALLY UPDATED HERE
-                    obj.version = {'1.2.2'};
+                    obj.version = {'1.3.0'};
 
 			end
 		end

Annoparams.m

@@ -40,6 +40,7 @@
 		% Load in ECG/peak detection
 		maxBPM = 150;                    % Sets window for detecting R peaks
 		pkthresh = 95;                   % Percentile of ECG signal above which a peak can be found
+        pkfilter = 0;                    % Additional filtering prior to peak detection on/off
 		lowpass = 1;                     % Low pass wavelet filter on/off
 		highpass = 1;                    % High pass wavelet filter on/off
 		wavelet_level_lowpass = 1;       % LPF freq is > samp freq/2^(wavelet_level_lowpass + 1)

VCG.m

@@ -110,7 +110,7 @@
 
         function l = length(obj); l = length(obj.VM); end
 
-        function QRS = peaks(obj,a); QRS = findpeaksecg(obj.VM, a.maxBPM, obj.hz, a.pkthresh); end
+        function QRS = peaks(obj,a); QRS = findpeaksecg(obj.VM, a.maxBPM, obj.hz, a.pkthresh, a.pkfilter); end
 
         function e = ecg(obj, transform_matrix)
             switch transform_matrix

braveheart_batch_windows.prj

@@ -8,13 +8,13 @@
       <file>${PROJECT_ROOT}\braveheart_batch_windows_resources\icon_24.png</file>
       <file>${PROJECT_ROOT}\braveheart_batch_windows_resources\icon_16.png</file>
     </param.icons>
-    <param.version>1.2.2</param.version>
+    <param.version>1.3.0</param.version>
     <param.authnamewatermark>Hans F. Stabenau and Jonathan W. Waks</param.authnamewatermark>
     <param.email>braveheart.ecg@gmail.com</param.email>
     <param.company />
     <param.summary />
     <param.description>BRAVHEART ECG/VCG Analysis Software -- Command Line Version
-Version 1.2.2</param.description>
+Version 1.3.0</param.description>
     <param.screenshot>${PROJECT_ROOT}\braveheart_splash.png</param.screenshot>
     <param.guid />
     <param.installpath.string>\braveheart_batch\</param.installpath.string>
@@ -50,7 +50,7 @@ Version 1.2.2</param.description>
 	 varargin  input arguments 
  Description 
  	 BRAVHEART ECG/VCG Analysis Software -- Command Line Version
-Version 1.2.2</param.help.text>
+Version 1.3.0</param.help.text>
     <unset>
       <param.company />
       <param.summary />
@@ -86,6 +86,7 @@ Version 1.2.2</param.help.text>
       <file>${PROJECT_ROOT}\ecg_axis.m</file>
       <file>${PROJECT_ROOT}\elecpos.mat</file>
       <file>${PROJECT_ROOT}\findHighpassLvl.m</file>
+      <file>${PROJECT_ROOT}\findpeakswavelet.m</file>
       <file>${PROJECT_ROOT}\getfiletype.m</file>
       <file>${PROJECT_ROOT}\HRV_Calc.m</file>
       <file>${PROJECT_ROOT}\leadidcodexyz.m</file>

braveheart_gui.m

@@ -44,7 +44,7 @@
 
 % Edit the above text to modify the response to help braveheart_gui
 
-% Last Modified by GUIDE v2.5 20-Jun-2024 10:57:50
+% Last Modified by GUIDE v2.5 11-Sep-2024 13:17:11
 
 % Update the current L&F for mac button issues...
 % Windows will use the normal Windows theme
@@ -411,7 +411,7 @@ function about_help_button_Callback(hObject, eventdata, handles)
 {'\fontsize{14}\it\bf \color[rgb]{0.09,0.078,0.377}';...
 'BRAVE\color[rgb]{0.89,0.016,0.016}H\fontsize{11}EART';...
 '\fontsize{8}\rm\color{black}(Beth Israel Analysis of Vectors of the Heart)';...
-'Version 1.2.2' ; 
+'Version 1.3.0' ; 
 ' ' ;...
 'Copyright 2016-2024  Hans F. Stabeneau and Jonathan W. Waks' ;...
 ' ' ;...
@@ -439,7 +439,7 @@ function about_help_button_Callback(hObject, eventdata, handles)
 {'\fontsize{18}\it\bf \color[rgb]{0.09,0.078,0.377}';...
 'BRAVE\color[rgb]{0.89,0.016,0.016}H\fontsize{14}EART\fontsize{11}';...
 '\rm\color{black}(Beth Israel Analysis of Vectors of the Heart)';...
-'Version 1.2.2' ; 
+'Version 1.3.0' ; 
 ' ' ;...
 'Copyright 2016-2024  Hans F. Stabeneau and Jonathan W. Waks' ;...
 ' ' ;...
@@ -4807,7 +4807,7 @@ function threshold_fig_button_Callback(hObject, eventdata, handles)
 % handles    structure with handles and user data (see GUIDATA)
 
 aps = pull_guiparams(hObject, eventdata, handles);
-threshold_figure_gui(handles.vcg.VM, aps.pkthresh, aps.maxBPM, handles.vcg.hz);
+threshold_figure_gui(handles.vcg, aps);
 
 
 % --- Executes on button press in pushbutton150.
@@ -5449,3 +5449,12 @@ function BRAVEHEART_GUI_WindowKeyPressFcn(hObject, eventdata, handles)
 end
 
 
+
+
+% --- Executes on button press in pkfilter_checkbox.
+function pkfilter_checkbox_Callback(hObject, eventdata, handles)
+% hObject    handle to pkfilter_checkbox (see GCBO)
+% eventdata  reserved - to be defined in a future version of MATLAB
+% handles    structure with handles and user data (see GUIDATA)
+
+% Hint: get(hObject,'Value') returns toggle state of pkfilter_checkbox

braveheart_gui_windows.prj

@@ -8,13 +8,13 @@
       <file>${PROJECT_ROOT}\braveheart_gui_windows_resources\icon_24.png</file>
       <file>${PROJECT_ROOT}\braveheart_gui_windows_resources\icon_16.png</file>
     </param.icons>
-    <param.version>1.2.2</param.version>
+    <param.version>1.3.0</param.version>
     <param.authnamewatermark>Hans F. Stabenau and Jonathan W. Waks</param.authnamewatermark>
     <param.email>bravheart.ecg@gmail.com</param.email>
     <param.company />
     <param.summary />
     <param.description>BRAVEHEART ECG/VCG Analysis Software -- GUI Version
-Version 1.2.2</param.description>
+Version 1.3.0</param.description>
     <param.screenshot>${PROJECT_ROOT}\braveheart_splash.png</param.screenshot>
     <param.guid />
     <param.installpath.string>\braveheart_gui\</param.installpath.string>
@@ -50,7 +50,7 @@ Version 1.2.2</param.description>
 	 varargin  input arguments 
  Description 
  	 BRAVEHEART ECG/VCG Analysis Software -- GUI Version
-Version 1.2.2</param.help.text>
+Version 1.3.0</param.help.text>
     <unset>
       <param.company />
       <param.summary />
@@ -95,6 +95,7 @@ Version 1.2.2</param.help.text>
       <file>${PROJECT_ROOT}\ecg_formats.csv</file>
       <file>${PROJECT_ROOT}\elecpos.mat</file>
       <file>${PROJECT_ROOT}\findHighpassLvl.m</file>
+      <file>${PROJECT_ROOT}\findpeakswavelet.m</file>
       <file>${PROJECT_ROOT}\getfiletype.m</file>
       <file>${PROJECT_ROOT}\HRV_Calc.m</file>
       <file>${PROJECT_ROOT}\leadidcodexyz.m</file>

findpeaksecg.m

@@ -20,10 +20,26 @@
 % This software is for research purposes only and is not intended to diagnose or treat any disease.
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
-function QRS = findpeaksecg(VM, maxbpm, freq, pkthresh)
+function QRS = findpeaksecg(VM, maxbpm, freq, pkthresh, filter)
+
+if filter == 1
+    [QRS, ~, ~] = findpeakswavelet(VM, maxbpm, freq, pkthresh, 'sym4');
+else
+    q = quantile(abs(VM), 100);
+    [~, QRS] = findpeaks(abs(VM), ...
+        'MinPeakHeight', q(pkthresh), 'MinPeakDistance', round(freq*60/maxbpm));
+end
+
+% QRS = QRS';
+% 
+% QRS
+%QRS2
+
+% figure
+% hold on
+% plot(VM)
+% scatter(QRS,VM(QRS))
+% scatter(QRS2,VM(QRS2))
 
-q = quantile(abs(VM), 100);
-[~, QRS] = findpeaks(abs(VM), ...
-    'MinPeakHeight', q(pkthresh), 'MinPeakDistance', round(freq*60/maxbpm));
 
 end

findpeakswavelet.m

@@ -0,0 +1,128 @@
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% BRAVEHEART - Open source software for electrocardiographic and vectorcardiographic analysis
+% findpeakswavelet.m -- Find R peaks using wavelet filtered signals 
+% Copyright 2016-2024 Hans F. Stabenau and Jonathan W. Waks
+% 
+% Source code/executables: https://github.com/BIVectors/BRAVEHEART
+% Contact: braveheart.ecg@gmail.com
+% 
+% BRAVEHEART is free software: you can redistribute it and/or modify it under the terms of the GNU 
+% General Public License as published by the Free Software Foundation, either version 3 of the License, 
+% or (at your option) any later version.
+%
+% BRAVEHEART is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; 
+% without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. 
+% See the GNU General Public License for more details.
+% 
+% You should have received a copy of the GNU General Public License along with this program. 
+% If not, see <https://www.gnu.org/licenses/>.
+%
+% This software is for research purposes only and is not intended to diagnose or treat any disease.
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+function [QRS, QRS_wavelet, R] = findpeakswavelet(VM, maxbpm, freq, pkthresh, wavelet)
+
+% Max level of wavelet decomposition based on signal length
+max_lvl = floor(log2(length(VM)));
+
+% Get freq bands for decomposition at each level (col 1 lower bound,
+% col2 upper bound)
+fc = zeros(max_lvl,2);
+
+for i = 1:max_lvl
+    fc(i,2) = freq / (2^i);
+    fc(i,1) = freq / (2^(i+1));
+end
+
+fc(max_lvl,1) = 0;
+
+% Want to keep frequencies 10-40 Hz which should capture most of QRS energy
+% without T wave (<10 Hz)
+Fa = 40;
+Fb = 10;
+
+% Find closest value to 40 in upper bound of freq cutoffs
+U = fc(:,2) - Fa;
+
+% Take smallest positive value so dont undercut the freq band
+U(U<0) = nan;
+[~,idxU] = min(U);
+
+% Find closest value to 10 in lower bound of freq cutoffs
+L = fc(:,1) - Fb;
+
+% Take smallest postive value so dont undercut the freq band
+L(L<0) = nan;
+[~,idxL] = min(L);
+
+% Decompose signal to max level using MODWT
+% Signal is mirrored to reduce edge effects
+M = modwt(mirror(VM), max_lvl, wavelet);
+
+% Reconstruct using only the coefficients in levels between idxL and idxU
+% and then square the resultant signal (R)
+Mrecon = zeros(size(M));
+Mrecon(idxU:idxL,:) = M(idxU:idxL,:);
+R = imodwt(Mrecon,wavelet);
+R = abs(R).^2;
+R = middlethird(R);
+
+% Find peaks in the reconstructed signal which will be estimates of the
+% actual R peaks.
+
+q = quantile(R, 100);
+
+[~, QRS_wavelet] = findpeaks(R, ...
+    'MinPeakHeight', q(pkthresh), 'MinPeakDistance', round(freq*60/maxbpm));
+
+% Now have to use the peaks in the wavelet signal as an estimate of where
+% to search for true peaks in VM since the max in the wavelet
+% reconstruction is not always going to correlate exactly with the max in
+% each QRS complex in VM
+
+% Open a window of 100 ms on each side of each peak detected in the wavelet
+% reconstruction signal
+win_len = 100;
+w = round(win_len * freq / 1000);
+QRS = zeros(1,length(QRS_wavelet));
+
+for i = 1:length(QRS_wavelet)
+    ll = QRS_wavelet - w;
+    ul = QRS_wavelet + w;
+
+    % Prevent indices from being out of bounds for QRS complexes near
+    % start/end of signal
+
+    ll(ll<0) = 1;
+    ul(ul>length(VM)) = length(VM);
+
+    % Just take maximum as should only be 1 peak
+    [~,QRS(i)] = max(VM(ll(i):ul(i)));
+
+    % Correct indexing
+    QRS(i) = QRS(i) + ll(i) - 1;
+
+end
+
+debug = 0;
+if debug
+    q1 = quantile(abs(VM), 100);
+    q2 = quantile(abs(R), 100);
+
+    s = max(VM)/max(R);
+    R=R*s;
+
+    figure
+    plot(R)
+    hold on
+    plot(VM)
+    scatter(QRS_wavelet,R(QRS_wavelet))
+    scatter(QRS, VM(QRS))
+
+    line([0 length(VM)],[q(pkthresh)*s q(pkthresh)*s])
+    %line([0 5000],[q2(95) q2(95)])
+
+end
+
+
+end

pull_guiparams.m

@@ -43,6 +43,13 @@
 
 % Peak Threshold for QRS detection
     aps.pkthresh = str2double(get(handles.pkthresh, 'String'));
+
+% Peak Threshold filtering
+    if get(handles.pkfilter_checkbox, 'Value') == 1
+        aps.pkfilter = 1;
+    else
+        aps.pkfilter = 0;
+	end
 
 % Lowpass Wavelet filtering        
     if get(handles.wavelet_filter_box, 'Value') == 1

push_guiparams.m

@@ -38,6 +38,9 @@ function push_guiparams(aps, hObject, eventdata, handles)
 
 % Peak Threshold for QRS detection
     set(handles.pkthresh, 'String', aps.pkthresh);
+
+% Peak Threshold filtering
+    set(handles.pkfilter_checkbox, 'Value', aps.pkfilter);
 
 % Lowpass Wavelet filtering        
    set(handles.wavelet_filter_box, 'Value', aps.lowpass);

scpopen.m

@@ -142,6 +142,10 @@
                                 Sect1Len = Sect1Len - 3 - len; 
 %% [tag,len,Sect1Len],          %% DEBUGGING information
                                 field = fread(fid,[1,len],'uchar');
+
+                            % Disable if field is not read correctly
+                            if ~isempty(field)
+
                                 if tag == 0,	
                                         ListOfRecommendedTags(ListOfRecommendedTags==tag)=[];
                                         HDR.Patient.Name = char(field);  %% LastName
@@ -195,8 +199,8 @@
                                 elseif tag == 9,
                                         HDR.Patient.Race = field;
                                 elseif tag == 10,
-field,
-					if (field(1)~=0)
+                                        field,
+					if (field(1)~=0) 
 	                                        HDR.Patient.Medication = field;
 					else	
 	                                        HDR.Patient.Medication.Code = {field(2:3)};
@@ -274,10 +278,12 @@
                                 else
                                         fprintf(HDR.FILE.stderr,'Warning SCOPEN: unknown tag %i (section 1)\n',tag);
                                 end;
+
+                            end
                         end;
                         if ~isempty(ListOfRequiredTags)
-                                fprintf(HDR.FILE.stderr,'Warning SCPOPEN: the following tags are required but missing in file %s\n',HDR.FileName);
-                                disp(ListOfRequiredTags);
+                                %fprintf(HDR.FILE.stderr,'Warning SCPOPEN: the following tags are required but missing in file %s\n',HDR.FileName);
+                                %disp(ListOfRequiredTags);
                         end;
                         if ~isempty(ListOfRecommendedTags)
                                 %fprintf(HDR.FILE.stderr,'Warning SCPOPEN: the following tags are recommended but missing in file %s\n',HDR.FileName);

test_braveheart.m

@@ -50,6 +50,7 @@
     ap = Annoparams();                  % Blank Annoparams class
     ap.maxBPM = 150;                    % Sets window for detecting R peaks
     ap.pkthresh = 95;                   % Percentile of ECG signal above which a peak can be found
+    ap.pkfilter = 0;                    % Additional filtering prior to peak detection on/off
     ap.lowpass = 1;                     % Low pass wavelet filter on/off
     ap.highpass = 1;                    % High pass wavelet filter on/off
     ap.wavelet_level_lowpass = 2;       % Lvl 2 at 500 Hz and Lvl 3 at 997 Hz is 62.5 Hz low-pass filter