Comprensión de la salida FFT

87

Necesito ayuda para comprender la salida del cálculo DFT / FFT.

Soy un ingeniero de software experimentado y necesito interpretar algunas lecturas del acelerómetro de teléfonos inteligentes, como encontrar las frecuencias principales. Desafortunadamente, dormí durante la mayoría de mis clases de EE en la universidad hace quince años, pero he estado leyendo sobre DFT y FFT durante los últimos días (aparentemente con poco éxito).

Por favor, no responda "ve a tomar una clase de EE". De hecho, planeo hacerlo si mi empleador me paga. :)

Pues aqui esta mi problema:

Capturé una señal a 32 Hz. Aquí hay una muestra de 1 segundo de 32 puntos, que he trazado en Excel.

ingrese la descripción de la imagen aquí

Luego obtuve un código FFT escrito en Java de la Universidad de Columbia (después de seguir las sugerencias en una publicación sobre " FFT confiable y rápido en Java ").

El resultado de este programa es el siguiente. Creo que está ejecutando una FFT en el lugar, por lo que reutiliza el mismo búfer tanto para la entrada como para la salida.

Before: 

Re: [0.887  1.645  2.005  1.069  1.069  0.69  1.046  1.847  0.808  0.617  0.792  1.384  1.782  0.925  0.751  0.858  0.915  1.006  0.985  0.97  1.075  1.183  1.408  1.575  1.556  1.282  1.06  1.061  1.283  1.701  1.101  0.702  ]

Im: [0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  ]

After: 

Re: [37.054  1.774  -1.075  1.451  -0.653  -0.253  -1.686  -3.602  0.226  0.374  -0.194  -0.312  -1.432  0.429  0.709  -0.085  0.0090  -0.085  0.709  0.429  -1.432  -0.312  -0.194  0.374  0.226  -3.602  -1.686  -0.253  -0.653  1.451  -1.075  1.774  ]

Im: [0.0  1.474  -0.238  -2.026  -0.22  -0.24  -5.009  -1.398  0.416  -1.251  -0.708  -0.713  0.851  1.882  0.379  0.021  0.0  -0.021  -0.379  -1.882  -0.851  0.713  0.708  1.251  -0.416  1.398  5.009  0.24  0.22  2.026  0.238  -1.474  ]

Entonces, en este punto, no puedo sacar ni cara ni cruz de la salida. Entiendo los conceptos de DFT, como que la parte real son las amplitudes de las ondas del coseno componentes y la parte imaginaria son las amplitudes de las ondas sinusoidales componentes. También puedo seguir este diagrama del gran libro " La guía del científico e ingeniero para el procesamiento de señales digitales ": ingrese la descripción de la imagen aquí

Entonces mis preguntas específicas son:

  1. A partir de la salida de la FFT, ¿cómo encuentro las "frecuencias más frecuentes"? Esto es parte de mi análisis de los datos de mi acelerómetro. ¿Debo leer las matrices real (coseno) o imaginaria (seno)?

  2. Tengo una entrada de 32 puntos en el dominio del tiempo. ¿No debería la salida de la FFT ser una matriz de 16 elementos para reales y una matriz de 16 elementos para imaginarios? ¿Por qué el programa me da salidas de matriz reales e imaginarias de tamaño 32?

  3. En relación con la pregunta anterior, ¿cómo analizo los índices en las matrices de salida? Dada mi entrada de 32 muestras muestreadas a 32 Hz, tengo entendido que una salida de matriz de 16 elementos debe tener su índice distribuido uniformemente hasta la mitad de la frecuencia de muestreo (de 32 Hz), así que estoy en lo correcto al entender que cada elemento de la matriz representa (32 Hz * 1/2) / 16 = 1 Hz?

  4. ¿Por qué la salida FFT tiene valores negativos? Pensé que los valores representan amplitudes de una sinusoide. Por ejemplo, la salida de Real [3] = -1.075 debería significar una amplitud de -1.075 para una onda cosenoidal de frecuencia 3. ¿Es correcto? ¿Cómo puede ser negativa una amplitud?

java audio signal-processing accelerometer fft stackoverflowuser2010
fuente
¿Qué le gustaría calcular a partir de las lecturas del acelerómetro: velocidad, distancia? El ruido de las lecturas del acelerómetro sigue la distribución gaussiana y no puedo ver cómo el ajuste de una onda sinusoidal solucionaría eso.
Ali
2
la etiqueta java debe eliminarse ya que es más genérica que para un idioma específico
user3791372
Mirando el origen de la Universidad de Columbia, no es nada eficiente. Es una implementación simple y no optimizada de Cooley-Tucky con tablas de búsqueda de mariposas, y la inversión de bits se realiza manualmente en lugar de usar las funciones de biblioteca existentes
Mark Jeronimus
@MarkJeronimus: ¿Puede recomendar una implementación FFT eficiente en Java? Si mal no recuerdo, la razón por la que elegí el código de la Universidad de Columbia fue que la biblioteca FFTW era demasiado compleja para ejecutarla en un teléfono inteligente Android.
stackoverflowuser2010
Encontré algunas implementaciones 'optimizadas' dispersas, pero son básicamente un algoritmo por tamaño N, por lo que si necesita un rango de tamaños, necesita todas esas rutinas. En la práctica, he utilizado principalmente Intel Integrated Performance Primitives (sí, de Java, a través de JNA), pero eso no es gratuito. En casa utilizo básicamente el mismo algoritmo que vinculó, pero escrito desde cero en 2005 usando un libro de texto. Es solo FFT (Transformada rápida de Fourier), nada tan "rápido" para justificar el nombre de "Fast FFT".
Mark Jeronimus

Respuestas:

85

  1. No debe buscar la parte real o imaginativa de un número complejo (que es su matriz real e imaginaria). En su lugar, desea buscar la magnitud de la frecuencia que se define como sqrt (real * real + imag * imag). Este número siempre será positivo. Ahora todo lo que tiene que buscar es el valor máximo (ignore la primera entrada en su matriz. Esa es su compensación de CC y no contiene información dependiente de la frecuencia).

  2. Obtiene 32 salidas reales y 32 imaginarias porque está utilizando una FFT compleja a compleja. Recuerde que ha convertido sus 32 muestras en 64 valores (o 32 valores complejos) extendiéndolos con cero partes imaginarias. Esto da como resultado una salida FFT simétrica donde el resultado de frecuencia ocurre dos veces. Una vez que esté listo para usar en las salidas 0 a N / 2, y una vez reflejado en las salidas N / 2 a N.En su caso, es más fácil simplemente ignorar las salidas N / 2 a N.No las necesita, son sólo un artefacto sobre cómo calcula su FFT.

  3. La frecuencia de la ecuación fft-bin es (bin_id * freq / 2) / (N / 2) donde freq es su frecuencia de muestra (también conocida como 32 Hz, y N es el tamaño de su FFT). En su caso, esto se simplifica a 1 Hz por contenedor. Los contenedores N / 2 a N representan frecuencias negativas (concepto extraño, lo sé). Para su caso, no contienen ninguna información significativa porque son solo un espejo de las primeras frecuencias N / 2.

  4. Tus partes reales e imaginarias de cada contenedor forman un número complejo. Está bien si las partes reales e imaginarias son negativas mientras que la magnitud de la frecuencia en sí es positiva (vea mi respuesta a la pregunta 1). Le sugiero que lea sobre números complejos. Explicar cómo funcionan (y por qué son útiles) excede lo que es posible explicar en una sola pregunta de stackoverflow.

Nota: es posible que también desee leer qué es la autocorrelación y cómo se utiliza para encontrar la frecuencia fundamental de una señal. Tengo la sensación de que esto es lo que realmente quieres.

Nils Pipenbrinck
fuente
1
Gracias. Con respecto a 1: vi esta página de Matlab que muestra un espectro de frecuencia ( mathworks.com/help/techdoc/ref/fft.html ). En esa página hay un gráfico con el título de "Espectro de amplitud de un solo lado de y (t)". ¿Está graficando la magnitud de la frecuencia como sugirió, sqrt (real ^ 2 + img ^ 2)? Respecto a 3: todavía no obtengo el resultado de 2Hz por bin. En mi caso, N = 32 y frecuencia = 32, ¿verdad? Entonces, hay N / 2 = 32/2 = 16 bins, y la frecuencia más alta (Nyquist) es freq / 2 = 32/2 = 16 Hz, lo que da como resultado 16 Hz por 16 bins, dando 1 Hz por bin.
stackoverflowuser2010
1
Sí, el gráfico muestra la magnitud del espectro - | Y (f) |. Las barras de valor absoluto significan magnitud. Ancho del contenedor = frecuencia de muestreo / tamaño de FFT. Su frecuencia de muestreo es de 32 hz, su tamaño de FFT es de 32. ¡Sí, tiene razón sobre el ancho del contenedor!
Matt Montag
Se corrigió la frecuencia del contenedor.
André Chalella
1
Buena respuesta, gracias! Perdón por llegar un poco tarde a la fiesta, pero tal vez podrías responderme ¿qué unidad es la magnitud de la frecuencia (como mencionaste en el punto 1), en general? En mi caso, en una señal de valores de un acelerómetro (es m / s ^ 2). No puedo entenderlo del todo.
Markus Wüstenberg
¡Fascinante! Mis barras de frecuencia de visualización de música salían todas reflejadas de izquierda a derecha; la respuesta # 2 explica esto !! ¡¡Loco!!
Ryan S
11

Ya tiene algunas buenas respuestas, pero solo agregaré que realmente necesita aplicar una función de ventana a sus datos de dominio de tiempo antes de la FFT; de lo contrario, obtendrá artefactos desagradables en su espectro, debido a la fuga espectral .

Paul R
fuente
Aprecio que haya pasado bastante tiempo desde esta respuesta. Sin embargo, ¿podrías explicar a qué tipo de artefactos te refieres?
MattHusz
1
@MattHusz: el término general para el origen de estos artefactos es "fuga espectral". Ahora agregué un enlace a la respuesta que explica esto. Sin embargo, la mejor manera en que puedo describir el efecto es que su espectro se "manchará" debido a la ventana rectangular implícita.
Paul R
6

1) Busque los índices en la matriz real con los valores más altos, además del primero (que es el componente DC). Probablemente necesite una frecuencia de muestreo considerablemente superior a 32 Hz y un tamaño de ventana más grande para obtener resultados significativos.

2) La segunda mitad de ambas matrices es el espejo de la primera mitad. Por ejemplo, observe que el último elemento de la matriz real (1.774) es el mismo que el segundo elemento (1.774), y el último elemento de la matriz imaginaria (1.474) es el negativo del segundo elemento.

3) La frecuencia máxima que puede captar a una frecuencia de muestreo de 32 Hz es 16 Hz ( límite de Nyquist ), por lo que cada paso es de 2 Hz. Como se señaló anteriormente, recuerde que el primer elemento es 0 Hz (es decir, el desplazamiento de CC).

4) Claro, una amplitud negativa tiene mucho sentido. Simplemente significa que la señal está "invertida": una FFT estándar se basa en un coseno, que normalmente tiene un valor = 1 en t = 0, por lo que una señal que tiene un valor = -1 en el momento = 0 tendría una amplitud negativa .

anochecer-inactivo-
fuente
Gracias por la respuesta. (1) ¿Quiere decir que puedo ignorar la matriz imaginaria (seno) y, de ser así, por qué? ¿Seguramente el componente del seno debe ser importante? (2) ¿Por qué ocurre este reflejo? ¿Es solo el resultado del algoritmo FFT? ¿La mayoría de la gente ignora la mitad reflejada? (3) ¿Cómo calculó los pasos de 2Hz? Entiendo el límite de Nyquist de 16 Hz, así que si hay 16 elementos de matriz (no reflejados), ¿cada elemento debe ser de 16 Hz / 16 = 1 Hz cada uno? (4) Para encontrar las frecuencias principales, ¿simplemente tomo el valor absoluto de los valores de amplitud en las matrices de salida?
stackoverflowuser2010
No debe buscar en la matriz real el valor más alto, y no puede ignorar la matriz sinusoidal / I. En su lugar, desea la magnitud del vector complejo compuesto. La duplicación se produce porque la mitad de la entrada (la matriz I) son todos ceros, por lo que el resultado tiene la mitad de grados de libertad. Puede ignorarlo si sus datos son estrictamente reales.
hotpaw2
@duskwuff Muchas gracias: respondiste a una pregunta que iba a publicar, si no hubiera encontrado tu respuesta: cómo interpretar la SEGUNDA parte de la FFT. Quiero modificar los datos y realizar lo inverso y seguí obteniendo solo la mitad de los resultados, porque modifiqué los datos incorrectos en esa parte. Gracias de nuevo.
Martin
(3), el valor de step = 2Hz permanece implícito para mí hasta ahora. Tenemos 16 contenedores, representados por una matriz de longitud = 16. Necesitamos describir todas las frecuencias de 0Hz a 16Hz. Supongo que cada contenedor describe una parte de ese rango, ¿no es así?
Krafter
@krafter Creo que se ha reducido a la mitad porque no se puede deducir una frecuencia de un solo valor (ya que no hay repetición).
JVE999
5

Tenga en cuenta que la "frecuencia más frecuente" se puede dispersar en varios contenedores FFT, incluso con una función de ventana. Por lo tanto, es posible que deba usar una ventana más larga, múltiples ventanas o interpolación para estimar mejor la frecuencia de cualquier pico espectral.

hotpaw2
fuente