creating a 1D vector

10 visualizaciones (últimos 30 días)
Ken Chi
Ken Chi el 10 de En. de 2020
Comentada: Meg Noah el 12 de En. de 2020
I have been given a task to create a function. However i dont believe i have completed the whole task.
. It only gives me a figure with the plotted with all the maximum points. how would i alter this to function? Have been working on this code for many hours and any hints very much appreciated

Respuestas (2)

Meg Noah
Meg Noah el 10 de En. de 2020
Editada: Meg Noah el 11 de En. de 2020
Here's one of many possible solutions (edited 2020/01/11 for off-by-one error in indexing):
clc
close all
clear all
ECG = dlmread('ECG.csv');
nx = numel(ECG); x = 1:nx;
idx = find(ECG > 200);
mask = zeros(nx,1);
mask(idx) = 1;
labeled = bwlabel(mask);
stats = regionprops(labeled,ECG,'ALL');
ipeak = zeros(length(stats),1);
for iblob = 1:length(stats)
ipeak(iblob) = stats(iblob).PixelIdxList(1) - 1 + ...
find( ECG(stats(iblob).PixelIdxList) == stats(iblob).MaxIntensity);
end
figure('color','white','position',[70 100 600 300]);
plot(x,ECG);
hold on
title({'Heartbeat data'})
xlabel('x')
ylabel('ECG(x)')
plot(x(ipeak),ECG(ipeak),'+');

Meg Noah
Meg Noah el 11 de En. de 2020
Here's a solution that has a step by step process.
clc
close all
clear all
%% Step 0
% load the data
ECG = dlmread('ECG.csv');
nx = numel(ECG); x = 1:nx;
Fs = 350; % [samples/s] frequency
T = 1/Fs; % [s] sampling period
N = 3000; % [samples] Length of signal
t = (0:N-1)*T; % [s] Time vector
% visualization
figure('color','white','position',[70 100 900 500]);
plot(1e3*t,ECG,'k','DisplayName','ECG Amplitude')
title({'Heartbeat data'})
xlabel('t (milliseconds)')
ylabel('X(t)')
hold on
%% Step 1
% First, determine points of the ECG scan exceed a given threshold.
threshold = 100;
idx = find(ECG > threshold);
% visualization
mask = nan(size(ECG)); mask(idx) = 1;
plot(1e3*t,mask.*ECG,'b','DisplayName',['ECG > ' num2str(threshold)])
%% Step 2
% Next, determine which of the data-points found in the previous step
% correspond to the same heartbeat. (Hint: data-points corresponding
% to one heartbeat will have consecutive indices, whereas transitions
% between heartbeats correspond to discontinuities in the index. diff
% function might help.
idxStartNewHeartBeat = [idx(1) idx(find(diff(idx)>1)+1)']';
% visualization
plot(1e3*t(idxStartNewHeartBeat),ECG(idxStartNewHeartBeat),'+r', ...
'DisplayName',['Unique Heartbeat ECG > ' num2str(threshold)]);
%% Step 3
% For each of the groups of data-points (corresponding to different
% heartbeats), find the location of the element with the maximum amplitude,
% and store this location in a vector.
mask(isnan(mask)) = 0;
labeledHeartbeats = bwlabel(mask);
stats = regionprops(labeledHeartbeats,ECG,'MaxIntensity'); %#ok<MRPBW>
maxAmplitudeValues = [stats.MaxIntensity];
idxMaxIndex = zeros(size(maxAmplitudeValues));
for iheartbeat = 1:length(maxAmplitudeValues)
idxMaxIndex(iheartbeat) = find(labeledHeartbeats == iheartbeat & ...
ECG == maxAmplitudeValues(iheartbeat));
end
% visualization
plot(1e3*t(idxMaxIndex),ECG(idxMaxIndex),'or', ...
'DisplayName',['Maximum in Heartbeat']);
%% Step 4
% The output should be a 1D vector containing the indices (one for each
% heart beat) corresponding to the local maxima.
% visualization
fprintf(1,'Indices of the heartbeat local maxima\n');
for iheartbeat = 1:length(idxMaxIndex)
fprintf(1,'%d\n',idxMaxIndex(iheartbeat));
ht = text(1e3*t(idxMaxIndex(iheartbeat)), ...
ECG(idxMaxIndex(iheartbeat))+30, ...
['idx = ' num2str(idxMaxIndex(iheartbeat))]);
set(ht,'rotation',90);
end
ylim([-100 500]);
hl = legend('location','northoutside');
It produces this image
HeartBeatMaxima.png
However, the original solution had an off by one error. But that doesn't even matter. The solution works.
mask = zeros(nx,1);
mask(idx) = 1;
labeled = bwlabel(mask);
stats = regionprops(labeled,ECG,'ALL');
ipeak = zeros(length(stats),1);
for iblob = 1:length(stats)
ipeak(iblob) = stats(iblob).PixelIdxList(1) - 1 + ...
find( ECG(stats(iblob).PixelIdxList) == stats(iblob).MaxIntensity);
end
figure('color','white','position',[70 100 600 300]);
plot(x,ECG);
hold on
title({'Heartbeat data'})
xlabel('x')
ylabel('ECG(x)')
plot(x(ipeak),ECG(ipeak),'+')
  1 comentario
Meg Noah
Meg Noah el 12 de En. de 2020
Solution as function and a program to call it is attached.
If there is something else that this code needs to do, please explain.

Iniciar sesión para comentar.

Etiquetas

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by