我想知道如何在 Python 的 FFT 图中检测新的峰值。 假设我有这个简单的情节: 我想自动测量噪声信号中的“相似度”或峰值位置,我尝试使用余弦相似度,但我的真实信号太嘈杂了,即使我在信号中添加了一个新峰值,我继续获得 0.9 的余弦值,因为它只有一个峰值。 这是我的真实信号的一个例子,我也有一个问题,我的信号可以在测量范围内移动,所以我无法获得稳定的频率阵列,它们可以在 +/- 100 Hz 的窗口内: 这是用于第一个 Plot 的代码:
import numpy as np
from pylab import *
import scipy.fftpack
# Number of samplepoints
N = 600
# sample spacing
T = 1.0 / 800.0
x = np.linspace(0.0, N*T, N)
y1 = np.sin(50.0 * 2.0*np.pi*x) + 0.5*np.sin(80.0 * 2.0*np.pi*x)+ 0.7*np.sin(30.0 * 2.0*np.pi*x)+ 0.5*np.sin(10.0 * 2.0*np.pi*x)
y2 = np.sin(50.0 * 2.0*np.pi*x) + 0.5*np.sin(80.0 * 2.0*np.pi*x)+ 0.2*np.sin(60.0 * 2.0*np.pi*x)+ 0.4*np.sin(40.0 * 2.0*np.pi*x)
yf1 = scipy.fftpack.fft(y1)
yf2 = scipy.fftpack.fft(y2)
xf = np.linspace(0.0, 1.0/(2.0*T), N/2)
fig, ax = plt.subplots()
plot(xf, 2.0/N * np.abs(yf1[:N/2]))
plot(xf, 2.0/N * np.abs(yf2[:N/2]))
xlabel('Freq (Hz)',fontsize=16,weight='bold')
ylabel('|Y(freq)|',fontsize=16,weight='bold')
ax = gca()
fontsize = 14
for tick in ax.xaxis.get_major_ticks():
tick.label1.set_fontsize(fontsize)
tick.label1.set_fontweight('bold')
for tick in ax.yaxis.get_major_ticks():
tick.label1.set_fontsize(fontsize)
tick.label1.set_fontweight('bold')
grid(True)
show()
def cosine_similarity(v1,v2):
"compute cosine similarity of v1 to v2: (v1 dot v2)/{||v1||*||v2||)"
sumxx, sumxy, sumyy = 0, 0, 0
for i in range(len(v1)):
x = v1[i]; y = v2[i]
sumxx += x*x
sumyy += y*y
sumxy += x*y
return sumxy/math.sqrt(sumxx*sumyy)
print 'Cosine Similarity', cosine_similarity(2.0/N * np.abs(yf1[:N/2]),2.0/N * np.abs(yf2[:N/2]))
虽然我也设置了一个阈值,但有时实际信号中的峰值可能小于预定义的阈值。 有什么想法吗?
最佳答案
有很多方法可以找到峰值,甚至可以插入它们的子样本位置。 找到峰后,只需检查是否找到新峰即可。
您可以使用 peakutils
包找到山峰。您可以在那里设置峰值之间的阈值和最小距离。
import numpy as np
from pylab import *
import scipy.fftpack
# Number of samplepoints
N = 600
# sample spacing
T = 1.0 / 800.0
x = np.linspace(0.0, N*T, N)
y1 = np.sin(50.0 * 2.0*np.pi*x) + 0.5*np.sin(80.0 * 2.0*np.pi*x)+ 0.7*np.sin(30.0 * 2.0*np.pi*x)+ 0.5*np.sin(10.0 * 2.0*np.pi*x)
y2 = np.sin(50.0 * 2.0*np.pi*x) + 0.5*np.sin(80.0 * 2.0*np.pi*x)+ 0.2*np.sin(60.0 * 2.0*np.pi*x)+ 0.4*np.sin(40.0 * 2.0*np.pi*x)
yf1 = scipy.fftpack.fft(y1)
yf2 = scipy.fftpack.fft(y2)
xf = np.linspace(0.0, 1.0/(2.0*T), N/2)
v1 = 2.0/N * np.abs(yf1[:N/2])
v2 = 2.0/N * np.abs(yf2[:N/2])
# Find peaks
import peakutils
peaks_ind1 = peakutils.indexes(v1, thres=0.2, min_dist=5)
peaks_ind2 = peakutils.indexes(v2, thres=0.2, min_dist=5)
dist_th_for_new_peaks = 3
new_peaks = []
for p in peaks_ind2:
found_new_peak = np.all(np.abs(p - peaks_ind1) > dist_th_for_new_peaks)
if found_new_peak:
new_peaks.append(p)
print("New Peak!! - %d" % p)
fig, ax = plt.subplots()
plot(xf, v1, color='blue')
plot(xf, v2, color='green')
for p in peaks_ind1:
ax.scatter(xf[p], v1[p], s=40, marker='s', color='blue', label='v1')
for p in peaks_ind2:
ax.scatter(xf[p], v2[p], s=40, marker='s', color='green', label='v2')
for p in new_peaks:
ax.scatter(xf[p], v2[p], s=40, marker='s', color='red', label='new peaks')
xlabel('Freq (Hz)',fontsize=16,weight='bold')
ylabel('|Y(freq)|',fontsize=16,weight='bold')
ax = gca()
fontsize = 14
for tick in ax.xaxis.get_major_ticks():
tick.label1.set_fontsize(fontsize)
tick.label1.set_fontweight('bold')
for tick in ax.yaxis.get_major_ticks():
tick.label1.set_fontsize(fontsize)
tick.label1.set_fontweight('bold')
ax.set_xlim([0,400])
ax.set_ylim([0,0.8])
grid(True)
show()
关于python - 检测 FFT 图中的峰值,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/37719613/