Algoritmo de Goertzel

El algoritmo de Goertzel es una técnica de procesamiento digital de señales (DSP) para la evaluación eficiente de los términos individuales de la transformada discreta de Fourier (DFT). Resulta útil en ciertas aplicaciones prácticas, como el reconocimiento de tonos de señalización multifrecuencia de doble tono (DTMF) producidos por los botones del teclado de un teléfono analógico tradicional . El algoritmo fue descrito por primera vez por Gerald Goertzel en 1958. ^{[ 1 ]}

Al igual que la DFT, el algoritmo de Goertzel analiza un componente de frecuencia seleccionable de una señal discreta . ^{[ 2 ] [ 3 ] [ 4 ]} A diferencia de los cálculos directos de DFT, el algoritmo de Goertzel aplica un único coeficiente de valor real en cada iteración, utilizando aritmética de valor real para secuencias de entrada de valor real. Para cubrir un espectro completo (excepto cuando se utiliza para un flujo continuo de datos donde los coeficientes se reutilizan para cálculos subsiguientes, lo que tiene una complejidad computacional equivalente a la DFT deslizante ), el algoritmo de Goertzel tiene un orden de complejidad mayor que los algoritmos de transformada rápida de Fourier (FFT), pero para calcular un pequeño número de componentes de frecuencia seleccionados, es numéricamente más eficiente. La estructura simple del algoritmo de Goertzel lo hace muy adecuado para procesadores pequeños y aplicaciones embebidas.

El algoritmo de Goertzel también puede utilizarse "a la inversa" como una función de síntesis sinusoidal, que requiere solo 1 multiplicación y 1 resta por muestra generada. ^{[ 5 ]}

El cálculo principal en el algoritmo de Goertzel tiene la forma de un filtro digital , y por esta razón el algoritmo a menudo se denomina filtro de Goertzel . El filtro opera sobre una secuencia de entrada.${\displaystyle x[n]}$en una cascada de dos etapas con un parámetro${\displaystyle \omega _{0}}$, dando la frecuencia a analizar, normalizada a radianes por muestra.

La primera etapa calcula una secuencia intermedia,${\displaystyle s[n]}$:

La segunda etapa aplica el siguiente filtro a${\displaystyle s[n]}$, produciendo una secuencia de salida${\displaystyle y[n]}$:

Se puede observar que la primera etapa del filtro es un filtro IIR de segundo orden con una estructura de forma directa . Esta estructura particular tiene la propiedad de que sus variables de estado internas son iguales a los valores de salida anteriores de esa etapa. Valores de entrada${\displaystyle x[n]}$para${\displaystyle n<0}$Se presume que todos son iguales a 0. Para establecer el estado inicial del filtro de modo que la evaluación pueda comenzar en la muestra${\displaystyle x[0]}$A los estados del filtro se les asignan valores iniciales.${\displaystyle s[-2]=s[-1]=0}$Para evitar riesgos de aliasing , frecuencia${\displaystyle \omega _{0}}$a menudo se restringe al rango de 0 a π (véase el teorema de muestreo de Nyquist-Shannon ); usar un valor fuera de este rango no carece de sentido, pero es equivalente a usar una frecuencia con alias dentro de este rango, ya que la función exponencial es periódica con un período de 2π en${\displaystyle \omega _{0}}$.

Se puede observar que el filtro de segunda etapa es un filtro FIR , ya que sus cálculos no utilizan ninguna de sus salidas anteriores.

Los métodos de la transformada Z se pueden aplicar para estudiar las propiedades de la cascada de filtros. La transformada Z de la primera etapa del filtro dada en la ecuación (1) es

La transformada Z de la segunda etapa del filtro dada en la ecuación (2) es

La función de transferencia combinada de la cascada de las dos etapas de filtro es entonces

Esto se puede transformar de nuevo en una secuencia equivalente en el dominio del tiempo, y los términos se desenrollan hasta el primer término de entrada en el índice${\displaystyle n=0}$:

Se puede observar que los polos de la transformada Z del filtro se encuentran en${\displaystyle e^{+j\omega _ {0}}}$y${\displaystyle e^{-j\omega _ {0}}}$, en un círculo de radio unitario centrado en el origen del plano de la transformada Z compleja. Esta propiedad indica que el proceso de filtrado es marginalmente estable y vulnerable a la acumulación de errores numéricos cuando se calcula utilizando aritmética de baja precisión y secuencias de entrada largas. ^{[ 6 ]} Christian Reinsch propuso una versión numéricamente estable . ^{[ 7 ]}

Para el caso importante del cálculo de un término DFT, se aplican las siguientes restricciones especiales.

• El filtrado finaliza en el índice${\displaystyle n=N}$, dónde${\displaystyle N}$es el número de términos en la secuencia de entrada de la DFT.
• Las frecuencias elegidas para el análisis de Goertzel están restringidas a la forma especial.

• El número de índice${\displaystyle k}$indicando que el "intervalo de frecuencia" de la DFT se selecciona del conjunto de números de índice

Haciendo estas sustituciones en la ecuación (6) y observando que el término${\displaystyle e^{+j2\pi k}=1}$La ecuación (6) toma entonces la siguiente forma:

Podemos observar que el lado derecho de la ecuación (9) es extremadamente similar a la fórmula definitoria para el término DFT.${\displaystyle X[k]}$, el término DFT para número de índice${\displaystyle k}$, pero no exactamente lo mismo. La suma mostrada en la ecuación (9) requiere${\displaystyle N+1}$términos de entrada, pero solo${\displaystyle N}$Los términos de entrada están disponibles al evaluar una DFT. Un recurso sencillo pero poco elegante es extender la secuencia de entrada.${\displaystyle x[n]}$con un valor artificial más${\displaystyle x[N]=0}$. ^{[ 8 ]} Podemos ver en la ecuación (9) que el efecto matemático sobre el resultado final es el mismo que eliminar el término${\displaystyle x[N]}$a partir de la suma, entregando así el valor DFT deseado.

Sin embargo, existe un enfoque más elegante que evita el paso de filtro adicional. De la ecuación (1), podemos observar que cuando el término de entrada extendido${\displaystyle x[N]=0}$se utiliza en el paso final,

Por lo tanto, el algoritmo se puede completar de la siguiente manera:

• finalizar el filtro IIR después de procesar el término de entrada${\displaystyle x[N-1]}$,
• aplicar la ecuación (10) para construir${\displaystyle s[N]}$a partir de los resultados anteriores${\displaystyle s[N-1]}$y${\displaystyle s[N-2]}$,
• aplicar la ecuación (2) con el valor calculado${\displaystyle s[N]}$valor y con${\displaystyle s[N-1]}$producido por el cálculo directo final del filtro.

Las dos últimas operaciones matemáticas se simplifican combinándolas algebraicamente:

Tenga en cuenta que detener las actualizaciones del filtro en el término${\displaystyle N-1}$y aplicar inmediatamente la ecuación (2) en lugar de la ecuación (11) omite las actualizaciones finales del estado del filtro, lo que produce un resultado con una fase incorrecta. ^{[ 9 ]}

La estructura de filtrado particular elegida para el algoritmo de Goertzel es la clave de sus cálculos DFT eficientes. Podemos observar que solo un valor de salida${\displaystyle y[N]}$se utiliza para calcular la DFT, por lo que se omiten los cálculos para todos los demás términos de salida. Dado que el filtro FIR no se calcula, los cálculos de la etapa IIR${\displaystyle s[0],s[1]}$, etc. pueden descartarse inmediatamente después de actualizar el estado interno de la primera etapa.

Esto parece generar una paradoja: para completar el algoritmo, la etapa del filtro FIR debe evaluarse una vez utilizando las dos últimas salidas de la etapa del filtro IIR, mientras que, para optimizar la eficiencia computacional, la iteración del filtro IIR descarta sus valores de salida. Aquí es donde se aplican las propiedades de la estructura de filtro de forma directa. Las dos variables de estado internas del filtro IIR proporcionan los dos últimos valores de la salida del filtro IIR, que son los términos necesarios para evaluar la etapa del filtro FIR.

• Según la teoría de la complejidad computacional , calcular un conjunto de${\displaystyle M}$Términos DFT utilizando${\displaystyle M}$aplicaciones del algoritmo de Goertzel en un conjunto de datos con${\displaystyle N}$valores con un "costo por operación" de${\displaystyle K}$tiene complejidad${\displaystyle O(KNM)}$.

Para calcular un único bin DFT${\displaystyle X(f)}$para una secuencia de entrada compleja de longitud${\displaystyle N}$El algoritmo de Goertzel requiere${\displaystyle 2N}$multiplicaciones y${\displaystyle 4\ N}$sumas/restas dentro del bucle, así como 4 multiplicaciones y 4 sumas/restas finales, para un total de${\displaystyle 2N+4}$multiplicaciones y${\displaystyle 4N+4}$sumas/restas. Esto se repite para cada uno de los${\displaystyle M}$frecuencias.

• Por el contrario, usar una FFT en un conjunto de datos con${\displaystyle N}$Los valores tienen complejidad${\displaystyle O(KN\log _{2}(N))}$.

Esto es más difícil de aplicar directamente porque depende del algoritmo FFT utilizado, pero un ejemplo típico es una FFT de base 2, que requiere${\displaystyle 2\log _{2}(N)}$multiplicaciones y${\displaystyle 3\log _{2}(N)}$sumas/restas por bin DFT , para cada uno de los${\displaystyle N}$contenedores.

En las expresiones de orden de complejidad, cuando el número de términos calculados${\displaystyle M}$es más pequeño que${\displaystyle \log N}$La ventaja del algoritmo de Goertzel es clara. Pero debido a que el código FFT es comparativamente complejo, el factor "costo por unidad de trabajo"${\displaystyle K}$suele ser mayor para una FFT, y la ventaja práctica favorece al algoritmo de Goertzel incluso para${\displaystyle M}$varias veces más grande que${\displaystyle \log _{2}(N)}$.

Como regla general para determinar si un algoritmo FFT de base 2 o un algoritmo de Goertzel es más eficiente, ajuste el número de términos.${\displaystyle N}$en el conjunto de datos hacia arriba hasta la potencia exacta de 2 más cercana, llamando a esto${\displaystyle N_{2}}$y es probable que el algoritmo de Goertzel sea más rápido si

${\displaystyle M\leq {\frac {5N_{2}}{6N}}\log _{2}(N_{2})}$

Las implementaciones de FFT y las plataformas de procesamiento tienen un impacto significativo en el rendimiento relativo. Algunas implementaciones de FFT ^{[ 11 ]} realizan cálculos internos de números complejos para generar coeficientes sobre la marcha, lo que aumenta significativamente su "costo K por unidad de trabajo". Los algoritmos FFT y DFT pueden usar tablas de valores de coeficientes precalculados para una mayor eficiencia numérica, pero esto requiere más accesos a los valores de los coeficientes almacenados en la memoria externa, lo que puede generar una mayor contención de caché que contrarresta parte de la ventaja numérica.

Ambos algoritmos logran una eficiencia aproximadamente dos veces mayor al usar datos de entrada de valor real en lugar de datos de valor complejo. Sin embargo, estas mejoras son inherentes al algoritmo de Goertzel, pero no se alcanzan con la FFT sin utilizar ciertas variantes del algoritmo especializadas en la transformación de datos de valor real .

• Algoritmo FFT de Bluestein (chirp-Z)
• Modulación por desplazamiento de frecuencia (FSK)
• Modulación por desplazamiento de fase (PSK)

• Proakis, JG; Manolakis, DG (1996), Procesamiento digital de señales: principios, algoritmos y aplicaciones , Upper Saddle River, NJ: Prentice Hall, págs. 480–481 , Bibcode : 1996dspp.book.....P 

• Algoritmo de Goertzel en Wayback Machine (archivado el 28/06/2018)
• Un algoritmo DSP para análisis de frecuencia
• El algoritmo de Goertzel de Kevin Banks
• Análisis del algoritmo de Goertzel realizado por Uwe Beis, en el que lo compara con un filtro analógico de paso bajo de Chebyshev de segundo orden.