Phân tích tần số Scipy/Numpy FFT

Question

Tôi đang tìm cách xoay trục tần số trong một fft (được thực hiện thông qua scipy.fftpack.fftfreq) thành tần số trong Hertz, chứ không phải là thùng hoặc thùng phân đoạn.

Tôi cố gắng để mã dưới đây để kiểm tra các FFT:

t = scipy.linspace(0,120,4000) 
acc = lambda t: 10*scipy.sin(2*pi*2.0*t) + 5*scipy.sin(2*pi*8.0*t) + 2*scipy.random.random(len(t)) 

signal = acc(t) 

FFT = abs(scipy.fft(signal)) 
FFT = scipy.fftpack.fftshift(FFT) 
freqs = scipy.fftpack.fftfreq(signal.size) 

pylab.plot(freqs,FFT,'x') 
pylab.show() 

Tỷ lệ lấy mẫu nên 4000 mẫu/120 giây = 33,34 mẫu/giây.

Tín hiệu có tín hiệu 2.0 Hz, tín hiệu 8.0 Hz và một số nhiễu ngẫu nhiên.

Tôi lấy FFT, lấy tần số và vẽ nó. Những con số này khá vô nghĩa. Nếu tôi nhân tần số bằng 33,34 (tần số lấy mẫu), thì tôi nhận được các đỉnh vào khoảng 8 Hz và 15 Hz, điều này có vẻ sai (tần số phải là một hệ số 4 tách rời nhau, không phải 2!).

Bất kỳ suy nghĩ nào về những gì tôi đang làm sai ở đây?

Answer 1

Tôi nghĩ rằng bạn không cần phải làm fftshift(), và bạn có thể vượt qua thời gian lấy mẫu để fftfreq():

import scipy 
import scipy.fftpack 
import pylab 
from scipy import pi 
t = scipy.linspace(0,120,4000) 
acc = lambda t: 10*scipy.sin(2*pi*2.0*t) + 5*scipy.sin(2*pi*8.0*t) + 2*scipy.random.random(len(t)) 

signal = acc(t) 

FFT = abs(scipy.fft(signal)) 
freqs = scipy.fftpack.fftfreq(signal.size, t[1]-t[0]) 

pylab.subplot(211) 
pylab.plot(t, signal) 
pylab.subplot(212) 
pylab.plot(freqs,20*scipy.log10(FFT),'x') 
pylab.show() 

từ đồ thị, bạn có thể thấy có hai đỉnh tại 2Hz và 8Hz.

Answer 2

Chiều rộng tần số của mỗi thùng là (sample_freq/num_bins).

Một vấn đề cơ bản nữa là tỷ lệ mẫu của bạn không đủ cho các tín hiệu bạn quan tâm. Tốc độ lấy mẫu của bạn là 8,3 Hz; bạn cần ít nhất 16Hz để chụp được âm đầu vào 8Hz.

---

_{1. Đối với tất cả các chuyên gia DSP; Tôi biết rằng nó thực sự là BW có liên quan, không phải tần số tối đa. Nhưng tôi giả định OP không muốn thực hiện thu thập dữ liệu chưa được lấy mẫu.}

Answer 3

Phương trình của bạn bị rối loạn.

fs = 33.33 
df1 = 2*pi * (2.0/fs) 
df2 = 2*pi * (5.0/fs) 
x = [10*sin(n*df1) + 5*sin(n*df2) + 2*random.random() for n in range(4000)] 

này cung cấp cho bạn 4000 mẫu chức năng của bạn, lấy mẫu tại 33,33 Hz, đại diện cho 120 giây của dữ liệu.

Bây giờ, hãy lấy FFT của bạn. Bin 0 sẽ giữ kết quả DC. Bin 1 sẽ là 33,33, bin 2 sẽ là 66,66, v.v.

Chỉnh sửa: Tôi quên đề cập rằng, vì tốc độ lấy mẫu của bạn là 33,33 Hz, tần suất tối đa có thể biểu diễn sẽ là fs/2 hoặc 16.665 Hz.

Answer 4

scipy.fftpack.fftfreq(n, d) cung cấp cho bạn các tần số trực tiếp. Nếu bạn đặt d=1/33.34, điều này sẽ cho bạn biết tần suất bằng Hz cho mỗi điểm của fft.