filtro de paso bajo y FFT para principiantes con Python

23

Soy nuevo en el procesamiento de señales y especialmente en FFT, por lo tanto, no estoy seguro de si estoy haciendo lo correcto aquí y estoy un poco confundido con el resultado.

Tengo una función real discreta (datos de medición) y quiero configurar un filtro de paso bajo en eso. La herramienta de elección es Python con el paquete numpy. Sigo este procedimiento:

  • calcular el fft de mi función
  • cortar altas frecuencias
  • realizar el inverso fft

Aquí está el código que estoy usando:

import numpy as np
sampling_length = 15.0*60.0 # measured every 15 minutes
Fs = 1.0/sampling_length
ls = range(len(data)) # data contains the function
freq = np.fft.fftfreq(len(data), d = sampling_length)
fft = np.fft.fft(data)
x = freq[:len(data)/2] 
for i in range(len(x)):
if x[i] > 0.005: # cut off all frequencies higher than 0.005
    fft[i] = 0.0
    fft[len(data)/2 + i] = 0.0
inverse = np.fft.ifft(fft)

¿Es este el procedimiento correcto? El resultado inversecontiene valores complejos, eso me confunde.

Hasta B
fuente
1
Cuando estaba aprendiendo FFT, encontré esta publicación de blog muy útil. glowingpython.blogspot.com/2011/08/…
David Poole

Respuestas:

23

Es de esperar que el resultado sea complejo. Quiero señalar un par de cosas:

Está aplicando un filtro de dominio de frecuencia de pared de ladrillo a los datos, intentando poner a cero todas las salidas FFT que corresponden a una frecuencia mayor que 0.005 Hz, luego transformando inversamente para obtener una señal de dominio de tiempo nuevamente. Para que el resultado sea real, la entrada a la FFT inversa debe ser simétrica conjugada . Esto significa que para una longitud- FFT,norte

X[k]=X[norte-k],k=1,2,...,norte2-1(nortemivminorte)

X[k]=X[norte-k],k=1,2,...,norte2(norteorere)
  • Tenga en cuenta que para par, y no son iguales en general, pero ambos son reales. Para impar , debe ser real.X [ 0 ] X [ NnorteX[0 0]X[norte2]norteX[0 0]

Veo que intentaste hacer algo como esto en tu código anterior, pero no es del todo correcto. Si aplica la condición anterior en la señal que pasa a la FFT inversa, entonces debería obtener una señal real.

Mi segundo punto es más filosófico: lo que estás haciendo funcionará, ya que suprimirá el contenido del dominio de frecuencia que no deseas. Sin embargo, esta no suele ser la forma en que se implementaría un filtro de paso bajo en la práctica. Como mencioné antes, lo que está haciendo es esencialmente aplicar un filtro que tiene una respuesta de magnitud de pared de ladrillo (es decir, perfectamente rectangular). La respuesta al impulso de dicho filtro tiene una forma . Dado que la multiplicación en el dominio de la frecuencia es equivalente a (en el caso de utilizar la circunvolución DFT, circular) en el dominio del tiempo, esta operación es equivalente a convolucionar la señal del dominio del tiempo con una función .syonortedo(X)syonortedo

¿Por qué es esto un problema? Recuerde cómo se ve la función en el dominio del tiempo (imagen de abajo tomada prestada descaradamente de Wikipedia):syonortedo

trama de la función sinc

La función tiene un soporte muy amplio en el dominio del tiempo; decae muy lentamente a medida que te alejas a tiempo de su lóbulo principal. Para muchas aplicaciones, esta no es una propiedad deseable; cuando convoluciona una señal con un , los efectos de los lóbulos laterales que decaen lentamente a menudo serán evidentes en la forma de dominio de tiempo de la señal de salida filtrada. Este tipo de efecto a menudo se denomina timbre . Si sabe lo que está haciendo, hay algunos casos en los que este tipo de filtrado podría ser apropiado, pero en el caso general, no es lo que desea.syonortedosyonortedo

Hay medios más prácticos para aplicar filtros de paso bajo, tanto en el dominio del tiempo como en el de la frecuencia. Los filtros de respuesta de impulso finito y de respuesta de impulso infinito se pueden aplicar directamente usando su representación de ecuación de diferencia . O, si su filtro tiene una respuesta de impulso suficientemente larga, a menudo puede obtener beneficios de rendimiento utilizando técnicas de convolución rápida basadas en la FFT (aplicando el filtro multiplicando en el dominio de frecuencia en lugar de convolución en el dominio de tiempo), como la superposición. guardar y superponer-agregar métodos.

Jason R
fuente
Sin embargo, la función sinc es el filtrado ideal, ¿no? Eso es lo que buscan todos los demás filtros, pero no lo logran. Es malo para el procesamiento de imágenes porque las imágenes no tienen antialias primero, por lo que produce un timbre que se ve horrible, pero para el audio u otras señales que fueron antialias filtradas antes del muestreo, ¿no es el mejor filtro que puede obtener?
endolito el
1
Sí, mi resultado no fue conjugado simétrico. Corrija el código, ahora todo funciona bien. ¡Gracias!
Hasta B
3
@endolith: un Sinc es un interpolador ideal para ciertos tipos de interpolación, pero puede estar lejos de ser ideal como filtro para la mayoría de los requisitos comunes de filtro, como la planitud de respuesta de banda de paso, detener el rechazo de banda, etc.
hotpaw2
+1 por la buena explicación sobre "por qué las personas no implementan el filtro como lo hace el PO"
Sibbs Gambling
Tienes que usar una ventana sinc. Si no tiene limitaciones de tiempo, este es el filtro óptimo, mucho mejor que Chebichev.