{{ message }}

Instantly share code, notes, and snippets.

# endolith/peakdet.m

Last active Oct 27, 2021
Peak detection in Python [Eli Billauer]
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters. Learn more about bidirectional Unicode characters
 function [maxtab, mintab]=peakdet(v, delta, x) %PEAKDET Detect peaks in a vector % [MAXTAB, MINTAB] = PEAKDET(V, DELTA) finds the local % maxima and minima ("peaks") in the vector V. % MAXTAB and MINTAB consists of two columns. Column 1 % contains indices in V, and column 2 the found values. % % With [MAXTAB, MINTAB] = PEAKDET(V, DELTA, X) the indices % in MAXTAB and MINTAB are replaced with the corresponding % X-values. % % A point is considered a maximum peak if it has the maximal % value, and was preceded (to the left) by a value lower by % DELTA. % Eli Billauer, 3.4.05 (Explicitly not copyrighted). % This function is released to the public domain; Any use is allowed. maxtab = []; mintab = []; v = v(:); % Just in case this wasn't a proper vector if nargin < 3 x = (1:length(v))'; else x = x(:); if length(v)~= length(x) error('Input vectors v and x must have same length'); end end if (length(delta(:)))>1 error('Input argument DELTA must be a scalar'); end if delta <= 0 error('Input argument DELTA must be positive'); end mn = Inf; mx = -Inf; mnpos = NaN; mxpos = NaN; lookformax = 1; for i=1:length(v) this = v(i); if this > mx, mx = this; mxpos = x(i); end if this < mn, mn = this; mnpos = x(i); end if lookformax if this < mx-delta maxtab = [maxtab ; mxpos mx]; mn = this; mnpos = x(i); lookformax = 0; end else if this > mn+delta mintab = [mintab ; mnpos mn]; mx = this; mxpos = x(i); lookformax = 1; end end end
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters. Learn more about bidirectional Unicode characters
 import sys from numpy import NaN, Inf, arange, isscalar, asarray, array def peakdet(v, delta, x = None): """ Converted from MATLAB script at http://billauer.co.il/peakdet.html Returns two arrays function [maxtab, mintab]=peakdet(v, delta, x) %PEAKDET Detect peaks in a vector % [MAXTAB, MINTAB] = PEAKDET(V, DELTA) finds the local % maxima and minima ("peaks") in the vector V. % MAXTAB and MINTAB consists of two columns. Column 1 % contains indices in V, and column 2 the found values. % % With [MAXTAB, MINTAB] = PEAKDET(V, DELTA, X) the indices % in MAXTAB and MINTAB are replaced with the corresponding % X-values. % % A point is considered a maximum peak if it has the maximal % value, and was preceded (to the left) by a value lower by % DELTA. % Eli Billauer, 3.4.05 (Explicitly not copyrighted). % This function is released to the public domain; Any use is allowed. """ maxtab = [] mintab = [] if x is None: x = arange(len(v)) v = asarray(v) if len(v) != len(x): sys.exit('Input vectors v and x must have same length') if not isscalar(delta): sys.exit('Input argument delta must be a scalar') if delta <= 0: sys.exit('Input argument delta must be positive') mn, mx = Inf, -Inf mnpos, mxpos = NaN, NaN lookformax = True for i in arange(len(v)): this = v[i] if this > mx: mx = this mxpos = x[i] if this < mn: mn = this mnpos = x[i] if lookformax: if this < mx-delta: maxtab.append((mxpos, mx)) mn = this mnpos = x[i] lookformax = False else: if this > mn+delta: mintab.append((mnpos, mn)) mx = this mxpos = x[i] lookformax = True return array(maxtab), array(mintab) if __name__=="__main__": from matplotlib.pyplot import plot, scatter, show series = [0,0,0,2,0,0,0,-2,0,0,0,2,0,0,0,-2,0] maxtab, mintab = peakdet(series,.3) plot(series) scatter(array(maxtab)[:,0], array(maxtab)[:,1], color='blue') scatter(array(mintab)[:,0], array(mintab)[:,1], color='red') show()

### wyojustin commented Jul 9, 2011

 Thanks for sharing your script. This saves me a huge headache! Justin

### dusita commented Aug 10, 2011

 Thanks a lot for sharing this script. It works well for me too. This script causes me register to github as I would like to say thank you.

### sixtenbe commented Aug 30, 2011

 Great script. Gave me a few good ideas about how to find peaks, but I think you should note that it if your signal starts with a minima you will systematically miss that one. Just remove the need for the: 'If lookoformax' statement and replace with some other logic and you can find the first peak, even if it's a minima. Doing this also means that you almost always will have a false hit on the first or one of the first samples and will need to pop that peak from the results. Changed this logic, added another test criteria to find the peaks and I also made a peak finder routine which divided the signal into different sections by looking for zero crossings and finding the peak in each section.

### endolith commented Sep 9, 2011

 That's great. I'll probably use that version in the future. For the "repeatable sinusoidal signals with some amount of RMS noise tolerable" version, why would you want the peak of the signal+noise? I would think you would want the peak of the signal itself?

### sixtenbe commented Sep 9, 2011

 What I mean by "some amount of RMS noise tolerable'" is that the function also works when you have a noisy signal, I'm not adding any noisy if that's what you think. The sensitivity of the function can be set with the window variable. For a very noisy signal with maybe a thousand or a few thousand samples per period a window of about 20 should be enough, but I set it quite high to get a good margin and it doesn't effect the final result anyway, as long as it can find the zero-crossings properly. Should probably change the description to "repeatable signals, where some noise is tolerated" or something to that end. As for performing interpolation I'm a little split about that, theoretically I think it might be a good practice, but experience hasn't always shown this. A colleague had a labView program for analysing waveforms, where he adapted a sinewave to every peak to find the actual peak smoothing away noise, but when we investigated the peaks we observed that it consistently choose a value lower than the actual peak and offset in time. So the solution we use is to simply find the highest measured peak without trying to do anything smart. Although at another company they have fewer sample per period than we, but they perform an fft of the signal and then padds the fft-array with zero values and calculates the ifft of it to interpolate new values. This way they can have a longer aperture time for each sample increasing their accuracy for each sample and should theoretically retain resolution in time. I would a define a realistic triangle as the following and by realistic I mean that you might actually encounter it on an actual power grid in a house or a factory: ```triangle = lambda T: [1000*sqrt(2)*(sin(t)+0.05*sin(t*3-pi) + 0.05*sin(t*5)+0.02*sin(t*7-pi)+0.01*sin(t*9)) for t in T] #or if you want a fixed frequency of the waveshape and not radians: triangle = lambda T, Hz: [1000*sqrt(2)*(sin(Hz*pi*t)+0.05*sin(Hz*pi*t*3-pi) + 0.05*sin(Hz*pi*t*5)+0.02*sin(Hz*pi*t*7-pi)+0.01*sin(Hz*pi*t*9)) for t in T]```

### endolith commented Sep 9, 2011

 Yeah, that's basically what I was thinking with fitting a sine wave to the noisy sine wave and finding the peak of that instead. And yes, the parabolic interpolation is for sharp narrow peaks in a spectrum, not for finding peaks of sine waves.

### sixtenbe commented Sep 9, 2011

 Oh I notice that I managed to remove the context for my definition of a triangle. The point of the triangle is that a triangle and a sine wave, with some noise can be a good way of testing any function for fitting or interpolating a peak. As for fitting sine waves, as I said I don't think it's worthwhile to fit any sine waves to the peak or interpolating it. For a pure sine it would also be good to compare a RMS calculation of the waveform with the Vpp/sqrt(8) where Vpp = difference between positive and negative peak. This is a good test to see if a function can find peaks for a pure sine wave. The triangle is useful when performing an optical inspection of the peak finding function. For just finding the maximum respecitve the minimum value as done in my peakdetect_zero_crossing function I've verified that on a real pure sine wave it will give results of no worse than 60 ppm, with 1000 samples per period and calculated over 5 periods.

### v3nz3n commented Oct 25, 2013

 thanks for the script. to make it work I had to add asarray and pylab imports: ``````from numpy import NaN, Inf, arange, isscalar, array, asarray from pylab import * `````` as well as actually show the matplotlib graph by adding (in the main function): ``````plot(series) show() ``````

### Zeokat commented Mar 5, 2014

 This piece of code saved Zeokat lots of hours. Thanks.

### endolith commented Mar 18, 2014

 @v3nz3n sorry, I'm not good at listing imports because I just did `from pylab import *`, but I'm trying to be better now. Spyder has pyflakes built in, which helps.

### Alymantara commented Aug 27, 2014

 Works perfectly. Thanks for sharing!!!

### San08 commented Sep 2, 2014

 I have a negative signal with a peak pointing downwards. When I use this code maxtab gives me the start of the peak and mintab gives me the maximum position of the peak. I would like to get the start and end of the peak. Is there a way to find?

### iuribeferrari commented Jan 15, 2015

 Thanks for sharing! A great help.

### MartinSamouiller commented Sep 1, 2015

 Thanks you !

### bancek commented Oct 11, 2015

 Thank you.

### nmningmei commented Mar 21, 2016

 Thank you! Great help for search target waves in EEG data.

### Ruthygg commented Mar 22, 2016

 Thanks a lot for this, do you have this in R?

### golobor commented Jul 1, 2016

 thanks for the code!

### algerianmaster commented Aug 3, 2016

 please how to use this with PyOpenCl ?

### mpiannucci commented Sep 20, 2016

 Here is a `go` version I converted https://github.com/mpiannucci/peakdetect

### Rishie13 commented Dec 3, 2016

 Is there a way I can preprocess csv files and show output for accelerometer?

### atikalidyad commented May 23, 2017

 can i implement it in GNURadio? Thank you very much

### JorgeHormaza commented Dec 4, 2017

 Mate, I just want to say Thank you so much because you save my life with this script.

### PZiso commented Mar 29, 2018

 Great coding ! Thanks a lot!

### LuisAli22 commented Jun 20, 2018

 Thank you so much

### 15738897318 commented Oct 7, 2018

 有关峰值点检测的研究促进这一领域的快速发展，我希望有志之士团结合作共同推进，同时我正在做峰值检测的学术课题，希望与有意者交流探讨，联系方式：m15738897318@gail.com,谢谢

### primehaxor commented Nov 17, 2018

 Thank you for this snippet, it works perfectly in my implementation.

### cgi1 commented Feb 4, 2020

 Thanks man!

### TravisH0301 commented Mar 10, 2021

 Thanks for sharing your code! It works very well for my project.