Hey!
My qeustion is: I need to remove all the R peaks. I have all the Q and S junctions points and I thought about connecting them with interp1 but I doesn't seem to work.
How can I connect all of the Q and S points at once?
Thank you in advance!
This is my code so far: (of course not the best, but it works)
- I am using this Pan–Tompkins algorithm implemention
(I have added one ECG file as an example for what I use)
SE = strel('line',2*Fs,0); % 2 seconds ("fast timescale")
SE_1 = strel('line',3*Fs,0); % 3 seconds ("slow timescale")
signal_open = imopen(signal, SE); % Morphologically opens "image"
signal_close = imclose(signal_open, SE_1); % Morphologically closes "opened image"
signal_filtered = signal - signal_close; % Subtracts computed reference value
%% Step 2: R peak detection
%
% Pan-Tompkins Algorithm to detect R peaks and to calculate the RR Intervals
%[R_Peak_amp,R_Peak_ind] = pan_tompkin(signal_filtered,Fs,0);if(debug_flag) %TEST
N = numel(R_Peak_ind); Max_QRS_duration = 0.12; Min_QRS_duration = 0.08; QRS_Duration = (Min_QRS_duration-Max_QRS_duration).*rand(1,N) + Max_QRS_duration; return;endQ_j_ind = zeros(0,length(R_Peak_ind));S_j_ind = zeros(0,length(R_Peak_ind));%% Step 3: search for the junctions of Q wave and S wave
first_peak = R_Peak_ind(1);last_peak = R_Peak_ind(end-1);q = 1;s = 1;while s < param_limit_length_to_N_samples if(s == first_peak) R_Peak_ind = R_Peak_ind(2:end); end if(s == R_Peak_ind(2)) d = 0; Left_Peak = R_Peak_ind(1); Current_Peak = R_Peak_ind(2); Right_Peak = R_Peak_ind(3); Y_CurrentPeak = signal_filtered(Current_Peak); Y_LeftPeak = signal_filtered(Left_Peak); Y_RightPeak = signal_filtered(Right_Peak); delta_Y = signal_filtered(Current_Peak) - signal_filtered(Left_Peak); delta_X = Current_Peak - Left_Peak; for Pi = Current_Peak:-1:Left_Peak d = d + 1; % delta_Y = Y_R - Y_L;
D_Pi = abs(delta_Y * Pi - delta_X * signal_filtered(Pi) + Current_Peak * Y_LeftPeak - Y_CurrentPeak * Left_Peak)/sqrt(delta_Y^2 + delta_X^2); W_Pi = cos(2*pi*(Pi-Current_Peak)/(Current_Peak-Left_Peak)); WD_Pi_q = D_Pi * W_Pi; WD_Pi_amp_q(d) = WD_Pi_q; end half = floor((Current_Peak-Left_Peak)/2); [~,Q_j] = max(WD_Pi_amp_q(1:half)); Q_j_ind(q) = Current_Peak - Q_j + 1; delta_Y = signal_filtered(Right_Peak) - signal_filtered(Current_Peak); delta_X = Right_Peak - Current_Peak; d = 0; for Pi = Current_Peak:Right_Peak d = d + 1; % delta_Y = Y_R - Y_L; D_Pi = abs(delta_Y * Pi - delta_X * signal_filtered(Pi) + Right_Peak * Y_CurrentPeak - Y_RightPeak * Current_Peak)/sqrt(delta_Y^2 + delta_X^2); W_Pi = cos(2*pi*(Pi-Current_Peak)/(Right_Peak-Current_Peak)); WD_Pi_s = D_Pi * W_Pi; WD_Pi_amp_s(d) = WD_Pi_s; end half = floor((Right_Peak-Current_Peak)/2); [~,S_j] = max(WD_Pi_amp_s(1:half)); S_j_ind(q) = Current_Peak + S_j - 1; q = q + 1; R_Peak_ind = R_Peak_ind(2:end); end if(s == last_peak) break; end s = s+1; endsubplot(2,1,1)plot(t,signal_filtered,'-o','MarkerEdgeColor','b','MarkerIndices',Q_j_ind);subplot(2,1,2)plot(t,signal_filtered,'-o','MarkerEdgeColor','r','MarkerIndices',S_j_ind);
Best Answer