how to calculate duration of peaks

4 visualizaciones (últimos 30 días)
Abdelhakim Souidi
Abdelhakim Souidi el 9 de Jul. de 2022
Editada: Star Strider el 11 de Jun. de 2023
I have an envelope signal attached y.mat and time array t.mat
I want to calculate the duration in second between the peaks (mentionned in red color in the figure below)
How I can reach this task

Respuesta aceptada

Star Strider
Star Strider el 9 de Jul. de 2022
Editada: Star Strider el 9 de Jul. de 2022
Try this —
LD1 = load('t.mat');
ty = LD1.ty.';
LD2 = load('y.mat');
sh1 = LD2.sh1;
y = sh1.';
[pks,plocs] = findpeaks(y, 'MinPeakProminence',0.5); % Desired Peak & Index
[vys,vlocs] = findpeaks(-y, 'MinPeakProminence',0.0000001); % Valleys & Indices
for k = 1:numel(pks)
vlcs(1,k) = max(vlocs(vlocs<=plocs(k))); % First Valley Index (Beginning Of Spike)
vlcs(2,k) = min(vlocs(vlocs>=plocs(k))); % Last Valley Index (End Of Spike)
end
deltat = -diff(ty(vlcs));
figure
plot(ty,sh1)
hold on
plot(ty(plocs), pks, '^r')
plot(ty(vlcs), sh1(vlcs),'.r')
hold off
grid
xlabel('Time')
ylabel('Amplitude')
text(mean(ty(vlcs)),zeros(size(plocs)), compose('\\leftarrow \\Delta t = %.6f',deltat), 'Horiz','left', 'Vert','middle', 'Rotation',60)
producing —
EDIT — (9 Jul 2022 at 17:24)
Minor edit to my original code to make it more robust and more efficient. It produces the same reault. I have no idea how to make it more robust with respect to other data sets, since I have only one data set to work with here, and it works correctly. It will be necessary to experiment with it to get it to work with other data sets. One possibiity is to experiment with the name-value pair arguments in the second findpeaks call to be certain that it detects the beginning and end ‘valleys’.for each peak, and rejects any noise that may be in the data. My code requires that each peak has valleys on either side of it, since they are in the posted data.
EDIT — (9 Jul 2022 at 22:58)
This slightly revised code works for both sets of files —
% LD1 = load('Abdelhakim Souidi t.mat');
LD1 = load('Abdelhakim Souidi (2) t.mat');
ty = LD1.ty.';
% LD2 = load('Abdelhakim Souidi y.mat');
LD2 = load('Abdelhakim Souidi(2) y.mat');
sh1 = LD2.sh1;
y = sh1.';
[pks,plocs] = findpeaks(y, 'MinPeakProminence',0.1, 'MinPeakDistance',250); % Desired Peak & Index
[vys,vlocs] = findpeaks(-y, 'MinPeakProminence',0.0001, 'MinPeakDistance',500); % Valleys & Indices
for k = 1:numel(pks)
vlcs(1,k) = max(vlocs(vlocs<=plocs(k))); % First Valley Index
vlcs(2,k) = min(vlocs(vlocs>=plocs(k))); % Last Valley Index
end
deltat = diff(ty(vlcs));
figure
plot(ty,sh1)
hold on
plot(ty(plocs), pks, '^r')
plot(ty(vlcs), sh1(vlcs),'.r')
hold off
grid
xlabel('Time')
ylabel('Amplitude')
text(mean(ty(vlcs)),zeros(size(plocs)), compose('\\leftarrow \\Delta t = %.6f',deltat), 'Horiz','left', 'Vert','middle', 'Rotation',60)
producing —
It matters not which peak of the double-peak envelopes are identified as the actual peak, since the baseline duration is the desired result, and that simply depends on the relation of individual ‘plocs’ indices and the nearest ‘vlocs’ indices.
.
  2 comentarios
Lam Ha
Lam Ha el 11 de Jun. de 2023
I follow your code and apply it on my data but I got the problems like it. It marks the higest point with x=y. I'm looking forward to hear suggestion from u.
This is my data
Star Strider
Star Strider el 11 de Jun. de 2023
Editada: Star Strider el 11 de Jun. de 2023
My code does exactly what it is supposed to do with your data.
Please post this as a new Question, specifically with a description of what you want to do, and include the URL of this thread.
I will not respond to it further here, since it is not directly relevant to this thread.

Iniciar sesión para comentar.

Más respuestas (2)

Image Analyst
Image Analyst el 9 de Jul. de 2022
Try this:
clc; % Clear the command window.
close all; % Close all figures (except those of imtool.)
clear; % Erase all existing variables. Or clearvars if you want.
workspace; % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 15;
st = load('t.mat')
t = st.ty;
sy = load('y.mat')
sh1 = sy.sh1;
y = sh1;
plot(t, y, 'b-', 'LineWidth', 2);
grid on;
hold on;
xlabel('t', 'FontSize', fontSize);
ylabel('y', 'FontSize', fontSize);
[peakValues, indexesOfPeaks] = findpeaks(y, 'MinPeakProminence', 0.25); % Desired Peak & Index
% Plot peaks
plot(t(indexesOfPeaks), peakValues, 'rv')
numPeaks = numel(peakValues)
% Fall down the left sides
leftIndexes = zeros(numPeaks, 1);
tLeft = zeros(numPeaks, 1);
for k = 1 : numPeaks
for k2 = indexesOfPeaks(k) : -1 : 1
thisValue = y(k2);
if thisValue < 0.1 && thisValue > y(k2+1)
leftIndexes(k) = k2;
tLeft(k) = t(leftIndexes(k));
% Put a green line at the start of the peak.
xline(tLeft(k), 'Color', 'g', 'LineWidth', 2);
break;
end
end
end
% Fall down the right sides
rightIndexes = zeros(numPeaks, 1);
tRight = zeros(numPeaks, 1);
for k = 1 : numPeaks
for k2 = indexesOfPeaks(k) : length(y)
thisValue = y(k2);
if thisValue < 0.1 && thisValue > y(k2-1)
rightIndexes(k) = k2;
tRight(k) = t(rightIndexes(k));
% Put a red line at the end of the peak.
xline(tRight(k), 'Color', 'r', 'LineWidth', 2);
break;
end
end
end
% Get the widths
peakWidths = tRight - tLeft
The left side is tLeft and is indicated by the green line.
The right side is tRight and is indicated by the red line.
The widths are given by peakWidths.
peakWidths =
0.146280579131303
0.157264103844234
0.155516724912631
0.159510733899151
0.170993509735397
0.159011482775836
0.177983025461807
0.179980029955067
0.150773839241138
0.159760359460809
  2 comentarios
Abdelhakim Souidi
Abdelhakim Souidi el 11 de Jul. de 2022
@Image Analyst you are the best. If you don't mind please how can I implement findchgpts to detect the boundaries
thank you
Image Analyst
Image Analyst el 11 de Jul. de 2022
@Abdelhakim Souidi I don't use it that much. All I'd do is look at the documentation and the examples there and try to apply it to your data, something you can do as well as me.

Iniciar sesión para comentar.


Image Analyst
Image Analyst el 9 de Jul. de 2022
Editada: Image Analyst el 9 de Jul. de 2022
  7 comentarios
Abdelhakim Souidi
Abdelhakim Souidi el 9 de Jul. de 2022
@Image Analyst here is an example
Image Analyst
Image Analyst el 11 de Jun. de 2023
Try adjusting the 'MinPeakDistance' to make sure your valleys are not too close together.

Iniciar sesión para comentar.

Productos


Versión

R2022a

Community Treasure Hunt

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

Start Hunting!

Translated by