Documentation

This is machine translation

Mouseover text to see original. Click the button below to return to the English version of the page.

goertzel

Discrete Fourier transform with second-order Goertzel algorithm

Syntax

```dft_data = goertzel(data) dft_data = goertzel(data,freq_indices) dft_data = goertzel(data,freq_indices,dim) ```

Description

`dft_data = goertzel(data)` returns the discrete Fourier transform (DFT) of the input data, `data`, using a second-order Goertzel algorithm. If `data` is a matrix, `goertzel` computes the DFT of each column separately. You can specify data types as double or single precision.

`dft_data = goertzel(data,freq_indices)` returns the DFT for the frequency indices `freq_indices`. The values of `freq_indices` can be any whole number or fraction.

`dft_data = goertzel(data,freq_indices,dim)` computes the DFT of the matrix `data` along the dimension `dim`.

Examples

collapse all

Estimate the frequencies of the tone generated by pressing the 1 button on a telephone keypad.

The number 1 produces a tone with frequencies 697 and 1209 Hz. Generate 205 samples of the tone with a sample rate of 8 kHz.

```Fs = 8000; N = 205; lo = sin(2*pi*697*(0:N-1)/Fs); hi = sin(2*pi*1209*(0:N-1)/Fs); data = lo + hi;```

Use the Goertzel algorithm to compute the DFT of the tone. Choose the indices corresponding to the frequencies used to generate the numbers 0 through 9.

```f = [697 770 852 941 1209 1336 1477]; freq_indices = round(f/Fs*N) + 1; dft_data = goertzel(data,freq_indices);```

Plot the DFT magnitude.

```stem(f,abs(dft_data)) ax = gca; ax.XTick = f; xlabel('Frequency (Hz)') title('DFT Magnitude')``` Generate a noisy cosine with frequency components at 1.24 kHz, 1.26 kHz, and 10 kHz. Specify a sample rate of 32 kHz. Reset the random number generator for reproducible results.

```rng default Fs = 32e3; t = 0:1/Fs:2.96; x = cos(2*pi*t*10e3) + cos(2*pi*t*1.24e3) + cos(2*pi*t*1.26e3) ... + randn(size(t));```

Generate the frequency vector. Use the Goertzel algorithm to compute the DFT. Restrict the range of frequencies to between 1.2 and 1.3 kHz.

```N = (length(x)+1)/2; f = (Fs/2)/N*(0:N-1); indxs = find(f>1.2e3 & f<1.3e3); X = goertzel(x,indxs);```

Plot the mean squared spectrum expressed in decibels.

```plot(f(indxs)/1e3,mag2db(abs(X)/length(X))) title('Mean Squared Spectrum') xlabel('Frequency (kHz)') ylabel('Power (dB)') grid``` Algorithms

The Goertzel algorithm implements the DFT as a recursive difference equation. To establish this difference equation, express the DFT as the convolution of an N-point input, x(n), with the impulse response$h\left(n\right)={W}_{N}^{-kn}u\left(n\right)$, where ${W}_{N}^{-kn}={e}^{-j2\pi k/N}$ and u(n) is the unit step sequence.

The Z-transform of the impulse response is

`$H\left(z\right)=\frac{1-{W}_{N}^{k}{z}^{-1}}{1-2\mathrm{cos}\left(2\pi k/N\right)\text{\hspace{0.17em}}{z}^{-1}+{z}^{-2}}.$`

The direct form II implementation is: Alternatives

You can also compute the DFT with:

• `fft`: less efficient than the Goertzel algorithm when you only need the DFT at a few frequencies.

• `czt`, the chirp Z-transform. `czt` calculates the Z-transform of an input on a circular or spiral contour and includes the DFT as a special case.

References

 Proakis, John G., and Dimitris G. Manolakis. Digital Signal Processing: Principles, Algorithms, and Applications. 3rd Edition. Upper Saddle River, NJ: Prentice Hall, 1996, pp. 480–481.

 Burrus, C. Sidney, and Thomas W. Parks. DFT/FFT and Convolution Algorithms: Theory and Implementation. New York: John Wiley & Sons, 1985.