Suponga lo siguiente:
- La frecuencia fundamental de una señal se ha estimado utilizando FFT y algunos métodos de estimación de frecuencia y se encuentra entre dos centros de depósito.
- La frecuencia de muestreo es fija
- El esfuerzo computacional no es un problema
Conociendo la frecuencia, ¿cuál es la forma más precisa de estimar el valor pico correspondiente de las señales fundamentales?
Una forma podría ser poner a cero la señal de tiempo para aumentar la resolución FFT de modo que el centro del depósito esté más cerca de la frecuencia estimada. En este escenario, un punto del que no estoy seguro es si puedo hacer cero-pad tanto como quiera o si hay algunos inconvenientes al hacerlo. Otro es el centro del contenedor que debería seleccionar después del relleno cero como el que estoy obteniendo el valor máximo (porque uno no puede alcanzar la frecuencia de interés exactamente, incluso después del relleno cero).
Sin embargo, también me pregunto si existe otro método que pueda ofrecer mejores resultados, por ejemplo, un estimador que utiliza los valores máximos de los dos centros de contenedores circundantes para estimar el valor máximo a la frecuencia de interés.
imax
está el pico FFT) le dará resultados precisosRespuestas:
El primer algoritmo que me viene a la mente es el algoritmo de Goertzel . Ese algoritmo generalmente supone que la frecuencia de interés es un múltiplo entero de la frecuencia fundamental. Sin embargo, este documento aplica el algoritmo (generalizado) al caso que le interesa.
Otro problema es que el modelo de señal es incorrecto. Lo usa
2*%pi*(1:siglen)*(Fc/siglen)
. Debería usarse2*%pi*(0:siglen-1)*(Fc/siglen)
para que la fase salga correctamente.También creo que hay un problema con la frecuencia
Fc=21.3
es muy baja. Las señales de baja frecuencia con valores reales tienden a mostrar sesgo cuando se trata de problemas de estimación de fase / frecuencia.También probé una búsqueda de grilla gruesa para la estimación de fase, y da la misma respuesta que el algoritmo de Goertzel.
A continuación se muestra un gráfico que muestra el sesgo en ambas estimaciones (Goertzel: azul, grueso: rojo) para dos frecuencias diferentes:
Fc=21.3
(sólido) yFc=210.3
(discontinua). Como puede ver, el sesgo de la frecuencia más alta es mucho menor.fuente
Si está dispuesto a usar varios contenedores FFT vecinos, no solo 2, la interpolación Sinc en ventana entre los resultados del contenedor complejo puede producir una estimación muy precisa, dependiendo del ancho de la ventana.
La interpolación Sinc con ventana se encuentra comúnmente en muestreadores de audio de alta calidad, por lo que los documentos sobre ese tema tendrán fórmulas de interpolación adecuadas con análisis de errores.
fuente
[1] JL Flanagan y RM Golden, "Phase vocoder", Bell Systems Technical Journal, vol. 45, págs. 1493–1509, 1966.
[2] K. Dressler, "Extracción sinusoidal utilizando una implementación eficiente de una FFT de resolución múltiple", en Proc. 9th Int. Conf. on Digital Audio Effects (DAFx-06), Montreal, Canadá, septiembre de 2006, págs. 247–252.
fuente
Un método es encontrar el máximo y ajustar una parábola al respecto, y luego usar el máximo de la parábola como la estimación de frecuencia y magnitud. Puede leer todo aquí: https://ccrma.stanford.edu/~jos/sasp/Sinusoidal_Peak_Interpolation.html
fuente
Tuve muchas dificultades con este problema exacto hace un par de años.
Publiqué esta pregunta:
/programming/4633203/extracting-precise-frequencies-from-fft-bins-using-phase-change-between-frames
Terminé haciendo los cálculos desde cero y publiqué una respuesta a mi propia pregunta.
Me sorprende que no haya podido encontrar ninguna exposición similar en Internet.
Publicaré la respuesta nuevamente aquí; tenga en cuenta que el código está diseñado para un escenario en el que estoy superponiendo mi ventana FFT por 4x.
π
Este rompecabezas tiene dos llaves para desbloquearlo.
La primera clave es comprender cómo la superposición de la ventana FFT introduce una rotación en la fase bin.
La segunda clave proviene del Gráfico 3.3 y 3.4 aquí (gracias a Stephan Bernsee por el permiso para copiar las fotos aquí).
Gráfico 3.3:
Gráfico 3.4:
Código:
fuente
Este código de Python le dará un resultado muy preciso (lo usé para muchas notas musicales y obtuve errores de menos del 0,01% de semitono) con interpolación parabólica (método utilizado con éxito por McAulay Quatieri, Serra, etc. en armónico + residual técnicas de separación)
fuente
Las frecuencias con las que está tratando (21.3Hz muestreado a 8kHz) son muy bajas. Debido a que estas son señales de valor real, exhibirán un sesgo en la estimación de fase para ** cualquier ** frecuencia.
Esta imagen muestra una gráfica del sesgo (
phase_est - phase_orig
) paraFc = 210.3;
(en rojo) versus el sesgo paraFc = 21.3;
. Como puede ver, el desplazamiento es mucho más significativo para el21.3
caso.Otra opción es reducir su frecuencia de muestreo. La curva verde muestra el sesgo para en
Fs = 800
lugar de8000
.fuente
goertzel = atan(imag(y(2)),real(y(2)))*180/%pi + 90;
. :-) Excavará un poco más. Mira este espacio.p
conp2 = (abs(y(3)) - abs(y(1)))/(2*(2*abs(y(2)) - abs(y(3)) - abs(y(1)))); phase2 = y0 - 0.25*(ym1-yp1)*p2;
, obtendrá MUCHO mejores respuestas, incluso paraFc=210
. No estoy del todo seguro de que la versión actual dep
le dé algo sensato. La fórmula de interpolación es para la interpolación de la AMPLITUD de una parábola, perop
está interpolando la fase que es simplemente ... extraña.p = (yp1 - ym1)/(2*(2*y0 - yp1 - ym1))
) será incorrecta algunas veces si está usando FASES en lugar de amplitudes. Esto se debe a que las fases pueden saltar alrededor del límite de +/- 180 grados. Todo lo que se necesita para arreglarlo para la fase es cambiar esa línea a mip2
cálculo anterior.