2015-06-15 103 views
0

我是信号处理的新手。在这里,我想问如何从python中获取来自FFT的FFT系数。这是我的代码的例子:使用Python的FFT系数

from scipy.fftpack import fft 
# Number of samplepoints 
N = 600 
# sample spacing 
T = 1.0/800.0 
x = np.linspace(0.0, N*T, N) 
y = np.sin(50.0 * 2.0*np.pi*x) + 0.5*np.sin(80.0 * 2.0*np.pi*x) 
yf = fft(y) 
xf = np.linspace(0.0, 1.0/(2.0*T), N/2) 
import matplotlib.pyplot as plt 
plt.plot(xf, 2.0/N * np.abs(yf[0:N/2])) 
plt.grid() 
plt.show() 

enter image description here

回答

0

嗯,我真的不知道信号处理要么但也许这个作品:

from scipy.signal import argrelmax 
f = xf[scipy.signal.argrelmax(yf[0:N/2])] 
Af = np.abs(yf[argrelmax(yf[0:N/2])]) 
0

引用@hotpaw,在this相似回答:

“实数和虚数组合在一起时可以表示一个复杂的数组。在频域中的复杂阵列的每一个复杂的元件可被认为是频率系数,并具有大小SQRT(R R +我 I))”。

所以,系数是在复杂的元件这个函数返回的是fft函数,此外,重要的是使用FFT函数的bin的大小(数量),测试一堆值并选择一个更有意义的值通常情况下,它的样本数量是相同的,这与大多数给出的答案一样,并且产生了很好的和合理的结果,如果需要探究的话,这里是我的代码版本:

%matplotlib inline 
import numpy as np 
import matplotlib.pyplot as plt 
import scipy.fftpack 

fig = plt.figure(figsize=[14,4]) 
N = 600   # Number of samplepoints 
Fs = 800.0 
T = 1.0/Fs  # N_samps*T (#samples x sample period) is the sample spacing. 
N_fft = 80  # Number of bins (chooses granularity) 
x = np.linspace(0, N*T, N)  # the interval 
y = np.sin(50.0 * 2.0*np.pi*x) + 0.5*np.sin(80.0 * 2.0*np.pi*x) # the signal 

# removing the mean of the signal 
mean_removed = np.ones_like(y)*np.mean(y) 
y = y - mean_removed 

# Compute the fft. 
yf = scipy.fftpack.fft(y,n=N_fft) 
xf = np.arange(0,Fs,Fs/N_fft) 

##### Plot the fft ##### 
ax = plt.subplot(121) 
pt, = ax.plot(xf,np.abs(yf), lw=2.0, c='b') 
p = plt.Rectangle((Fs/2, 0), Fs/2, ax.get_ylim()[1], facecolor="grey", fill=True, alpha=0.75, hatch="/", zorder=3) 
ax.add_patch(p) 
ax.set_xlim((ax.get_xlim()[0],Fs)) 
ax.set_title('FFT', fontsize= 16, fontweight="bold") 
ax.set_ylabel('FFT magnitude (power)') 
ax.set_xlabel('Frequency (Hz)') 
plt.legend((p,), ('mirrowed',)) 
ax.grid() 

##### Close up on the graph of fft####### 
# This is the same histogram above, but truncated at the max frequence + an offset. 
offset = 1 # just to help the visualization. Nothing important. 
ax2 = fig.add_subplot(122) 
ax2.plot(xf,np.abs(yf), lw=2.0, c='b') 
ax2.set_xticks(xf) 
ax2.set_xlim(-1,int(Fs/6)+offset) 
ax2.set_title('FFT close-up', fontsize= 16, fontweight="bold") 
ax2.set_ylabel('FFT magnitude (power) - log') 
ax2.set_xlabel('Frequency (Hz)') 
ax2.hold(True) 
ax2.grid() 

plt.yscale('log') 

输出: enter image description here