How to get the indices corresponding to the width of a peak?

조회 수: 12 (최근 30일)
Nut
Nut 2017년 12월 19일
댓글: Damien Bestard 2023년 4월 7일
Hi,
I'm using the findpeaks function to find the peaks of an array. After this, I should get the indices of the array corresponding to the width of the peak, both if calculated referring to the half-heigth or to the half-prominence of the peak.
Please, look at the attached images (taken from the findpeaks reference): in other words, I'm interested to find the indices corresponding to the horizontal yellow segments showed into the images.
Is there a way to do this? I'm looking for an existent function or method, but also any suggestion about the approach is welcome.
Finally, I'm also wondering about the possibility of extending this concept, i.e. to get the indices corresponding to a quote given as input, not only at 0.5*height or 0.5*prominence (as example, the width at 0.7*prominence, or at 0.3*height). Any suggestions?
Thank you so much for answering.
Kind regards.
  댓글 수: 1
Damien Bestard
Damien Bestard 2023년 4월 7일
Hi,
I was looking for the exact same thing (I want to integrate the signal covered by the horizontal yellow segments, to calculate the corresponding percentage of the total energy of the signal contained in each peak). And I find a way to do exactly this (skip the bullet points if you just want the result without the investigations) :
  • Open the findpeaks function (using debug mode and ctrl+D with the text cursor somewhere on the function name in my code) and search for "annotate" with ctrl+F. It will quickly lead you to lines 187-190 :
hAxes = plotSignalWithPeaks(x,y,iPk);
if strcmp(annotate,'extents')
plotExtents(hAxes,x,y,iPk,bPk,bxPk,byPk,wxPk,refW);
end
which seems to be where it plots the so called "horizontal yellow segments".
  • It seems that the parameter wxPk contains what we want. The function returns wPk as 3-rd output, which is pretty similar. To verify, you can place a breakpoint on the line 189 (the one calling plotExtents() ), and lauch the findpeaks function to get paused on this line. Then copy it in your command window, the figure will appear and the values in wxPk will describe exactly the bounds of each yellow segment.
Once this investigation done, it turns out that all you need to do is the following :
_______________________________
Since you aren't allowed to modify directly this function to add wkPk as an output, just copy its code in a new function file that you can rename (my_findpeaks.m for example). There you just have to modify the header :
function [Ypk,Xpk,Wpk,Ppk] = my_findpeaks(Yin,varargin)
in :
function [Ypk,Xpk,Wpk,Ppk,wxPk] = my_findpeaks(Yin,varargin)
This allows you to get the wanted parameter when calling the function in your code :
[pks,loc,w,p,wxPk] = my_findpeaks(...).
And wxPk(n,:) will give you the two boundaries of the n-th peak.
_______________________________
Done !
Regards.

댓글을 달려면 로그인하십시오.

채택된 답변

Csaba
Csaba 2017년 12월 20일
I do not know why these 0.5*height indices are important since they actually do not have any meaning. When the peaks are overwhelming 0.5*height can be practically anywhere for a given height. Especially if the vector contains noisy data. The correct way would be to fit peaks.
But if you want to find the half value, you can use the find function and look for the 0.5*height values in the vector. Then you can choose the closest or any of the found indices as you wish. Like:
indices=find(X<height/2)
where X is your vector (or part of your vector).if you want to find the left 0.5*height value then
indice=find(X(1:peakindice)<height/2,'last')
for right indice=find(X(peakindice:end)<height/2,'first')
But if your vector is rather noisy it will not be a correct value anyway.
  댓글 수: 2
Nut
Nut 2017년 12월 20일
Very good answer, I think it will be useful for me. It is also very simple, I don't know why I didn't think it before.
Thank you very much for your help.
Csaba
Csaba 2017년 12월 20일
You are welcome.

댓글을 달려면 로그인하십시오.

추가 답변 (0개)

카테고리

Help CenterFile Exchange에서 Descriptive Statistics에 대해 자세히 알아보기

Community Treasure Hunt

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

Start Hunting!

Translated by