2015-06-15 111 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