0.4 The dft as convolution or filtering  (Page 4/4)

It is possible to modify this algorithm to allow entering the data in forward order rather than reverse order. The difference [link] becomes

if [link] becomes

for $y\left(0\right)=0$ . This is the algorithm programmed later.

The second-order goertzel algorithm

One of the reasons the first-order Goertzel algorithm does not improve efficiency is that the constant in the feedback or recursive path iscomplex and, therefore, requires four real multiplications and two real additions. A modification of the scheme to make it second-order removesthe complex multiplications and reduces the number of required multiplications by two.

Define the variable $q\left(m\right)$ so that

This substituted into the right-hand side of [link] gives

Combining [link] and [link] gives the second order difference equation

which together with the output [link] , comprise the second-order Goertzel algorithm where

for initial conditions $q\left(0\right)=q(-1)=0$ .

A similar development starting with [link] gives a second-order algorithm with forward ordered input as

with

and for $q\left(0\right)=q(-1)=0$ .

Note that both difference [link] and [link] are not changed if $z$ is replaced with $z^{-1}$ , only the output [link] and [link] are different. This means that the polynomial $X\left(z\right)$ may be evaluated at a particular $z$ and its inverse $z^{-1}$ from one solution of the difference [link] or [link] using the output equations

and

Clearly, this allows the DFT of a sequence to be calculated with half the arithmetic since the outputs are calculated two at a time. Thesecond-order DE actually produces a solution $q\left(m\right)$ that contains two first-order components. The output equations are, in effect, zeros thatcancel one or the other pole of the second-order solution to give the desired first-order solution. In addition to allowing the calculating oftwo outputs at a time, the second-order DE requires half the number of real multiplications as the first-order form. This is because thecoefficient of the $q(m-2)$ is unity and the coefficient of the $q(m-1)$ is real if $z$ and $z^{-1}$ are complex conjugates of each other which is true for the DFT.

Analysis of arithmetic complexity and timings

Analysis of the various forms of the Goertzel algorithm from their programs gives the following operation count for real multiplications andreal additions assuming real data.

Timings of the algorithms on a PC in milliseconds are given in the following table.

These timings track the floating point operation counts fairly well.

Conclusions

Goertzel's algorithm in its first-order form is not particularly interesting, but the two-at-a-time second-order form is significantlyfaster than a direct DFT. It can also be used for any polynomial evaluation or for the DTFT at unequally spaced values or for evaluating afew DFT terms. A very interesting observation is that the inner-most loop of the Glassman-Ferguson FFT [link] is a first-order Goertzel algorithm even though that FFT is developed in a very different framework.

In addition to floating-point arithmetic counts, the number of trigonometric function evaluations that must be made or the size of a table to storeprecomputed values should be considered. Since the value of the $W^{nk}$ terms in [link] are iteratively calculate in the IIR filter structure, there is round-off error accumulation that should be analyzed in anyapplication.

It may be possible to further improve the efficiency of the second-order Goertzel algorithm for calculating all of the DFT of a number sequence.Perhaps a fourth order DE could calculate four output values at a time and they could be separated by a numerator that would cancel three of thezeros. Perhaps the algorithm could be arranged in stages to give an $Nlog\left(N\right)$ operation count. The current algorithm does not take into account any of the symmetries of the input index. Perhaps some of theideas used in developing the QFT [link] , [link] , [link] could be used here.

The quick fourier transform (qft)

One stage of the QFT can use the symmetries of the sines and cosines to calculate a DFT more efficiently than directly implementing the definition Multidimensional Index Mapping: Equation 1 . Similar to the Goertzel algorithm, the one-stage QFT is a better $N^2$ DFT algorithm for arbitrary lengths. See The Cooley-Tukey Fast Fourier Transform Algorithm: The Quick Fourier Transform, An FFT based on Symmetries .