gabe
gabe

Reputation: 193

Finding only the "prominent" local maxima of a 1d array

I have a couple of data sets with clusters of peaks that look like the following: Data sets with clusters of peaks You can see that the main features here a clusters of peaks, each cluster having three peaks. I would like to find the x values of those local peaks, but I am running into a few problems. My current code is as follows:

import numpy as np
import matplotlib.pyplot as plt
from scipy import loadtxt, optimize
from scipy.signal import argrelmax

def rounddown(x):
    return int(np.floor(x / 10.0)) * 10

pixel, value = loadtxt('voltage152_4.txt', unpack=True, skiprows=0)

ax = plt.axes()
ax.plot(pixel, value, '-')
ax.axis([0, np.max(pixel), np.min(value), np.max(value) + 1])

maxTemp = argrelmax(value, order=5)
maxes = []
for maxi in maxTemp[0]:
    if value[maxi] > 40:
        maxes.append(maxi)

ax.plot(maxes, value[maxes], 'ro')

plt.yticks(np.arange(rounddown(value.min()), value.max(), 10))
plt.savefig("spectrum1.pdf")
plt.show()

Which works relatively well, but still isn't perfect. Some peaks labeled: Some peaks labeled The main problem here is that my signal isn't smooth, so a few things that aren't actually my relevant peaks are getting picked up. You can see this in the stray maxima about halfway down a cluster, as well as peaks that have two maxima where in reality it should be one. You can see near the center of the plot there are some high frequency maxima. I was picking those up so I added in the loop only considering values above a certain point.

I am afraid that smoothing the curve will actually make me loose some of the clustered peaks that I want, as in some of my other datasets there are even closer together. Maybe my fears are unfounded, though, and I am just misunderstanding how smoothing works. Any help would be appreciated.

Does anyone have a solution on how to pick out only "prominent" peaks? That is, only those peaks that are quick large compared to the others?

Upvotes: 3

Views: 6081

Answers (2)

Richard
Richard

Reputation: 3394

I'd also recommend scipy.signal.find_peaks for what you're looking for. The other, older, scipy alternate find_peaks_cwt is quite complicated to use.

It will basically do what you're looking for in a single line. Apart from the prominence parameter that lagru mentioned, for your data either the threshold or height parameters might also do what you need.

height = 40 would filter to get all the peaks you like.

Prominence is a bit hard to wrap your head around for exactly what it does sometimes.

Upvotes: 0

lagru
lagru

Reputation: 231

Starting with SciPy version 1.1.0 you may also use the function scipy.signal.find_peaks which allows you to select detected peaks based on their topographic prominence. This function is often easier to use than find_peaks_cwt. You'll have to play around a little bit to find the optimal lower bound to pass as a value to prominence but e.g. find_peaks(..., prominence=5) will ignore the unwanted peaks in your example. This should bring you reasonably close to your goal. If that's not enough you might do your own peak selection based upon peak properties like the left_/right_bases which are optionally returned.

Upvotes: 2

Related Questions