python - fourier - Scipy/Numpy Análisis de frecuencia FFT
fft python example (4)
Estoy buscando cómo convertir el eje de frecuencia en fft (tomado a través de scipy.fftpack.fftfreq) en una frecuencia en hercios, en lugar de bins o bins fraccionales.
Traté de codificar a continuación para probar la 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()
La velocidad de muestreo debe ser 4000 muestras / 120 segundos = 33.34 muestras / seg.
La señal tiene una señal de 2.0 Hz, una señal de 8.0 Hz y algo de ruido aleatorio.
Tomo la FFT, tomo las frecuencias y lo trazo. Los números son bastante absurdos. Si multiplico las frecuencias por 33.34 (la frecuencia de muestreo), entonces obtengo picos de alrededor de 8 Hz y 15 Hz, lo que parece estar mal (¡además, las frecuencias deben ser un factor de 4 aparte, no 2!).
¿Alguna idea sobre lo que estoy haciendo mal aquí?
Creo que no necesita hacer fftshift (), y puede pasar el período de muestreo a 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()
del gráfico se puede ver que hay dos picos a 2Hz y 8Hz.
El ancho de frecuencia de cada contenedor es (sampling_freq / num_bins).
Un problema más fundamental es que su frecuencia de muestreo no es suficiente para sus señales de interés. Su frecuencia de muestreo es 8.3 Hz; necesita al menos 16Hz para capturar un tono de entrada de 8Hz. 1
1. A todos los expertos de DSP; Soy consciente de que en realidad es BW lo que es relevante, no la frecuencia máxima. Pero supongo que OP no quiere hacer una adquisición de datos submuestreada.
Tu ecuación está en mal estado.
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)]
Esto le da 4000 muestras de su función, muestreadas a 33.33 Hz, que representan 120 segundos de datos.
Ahora toma tu FFT. El contenedor 0 mantendrá el resultado de DC. Bin 1 será 33.33, bin 2 será 66.66, etc.
Editar: Olvidé mencionar que, dado que la frecuencia de muestreo es 33.33 Hz, la frecuencia máxima que se puede representar será fs / 2 o 16.665 Hz.
scipy.fftpack.fftfreq (n, d) te da las frecuencias directamente. Si configura d=1/33.34
, esto le indicará la frecuencia en Hz para cada punto de fft.