Forums

Goertzel Algorithm: Help Needed

Started by Sohaib Afzal December 31, 2009
Hi

I am trying to use the goertzel algo to make a real time audio
spectrum analyzer using a PIC controller. I have started out by trying
it on Matlab, but i am getting this problem:

The high frequency peaks in the power graph are much smaller than the
low frequency peaks. For example, the peak at 60Hz is, say 50 units on
the plot, then the peak at 8000Hz is only 0.023 units or so. Windows
Media Player and other players show the peaks are approximately the
same size throughout the sound file (.wav) that i am testing.

Please tell me if this algorithm is appropriate for what i am doing,
and if this difference in the magnitudes of peaks for power is usual
or due to some error in coding.

Thank you for your time.
On Thu, 31 Dec 2009 09:39:02 -0800, Sohaib Afzal wrote:

> Hi > > I am trying to use the goertzel algo to make a real time audio spectrum > analyzer using a PIC controller. I have started out by trying it on > Matlab, but i am getting this problem: > > The high frequency peaks in the power graph are much smaller than the > low frequency peaks. For example, the peak at 60Hz is, say 50 units on > the plot, then the peak at 8000Hz is only 0.023 units or so. Windows > Media Player and other players show the peaks are approximately the same > size throughout the sound file (.wav) that i am testing. > > Please tell me if this algorithm is appropriate for what i am doing, and > if this difference in the magnitudes of peaks for power is usual or due > to some error in coding. > > Thank you for your time.
Well, from the information you've given I can conclusively say that either the algorithm that you've chosen is screwed, or your code is. Hmm. That's a bit ambiguous. Could you post the difference equation that you're trying to execute? Please don't post all your Matlab code, unless it's no more than about ten lines. If you _do_ feel a need to post a Matlab program that's as long as War and Peace, extract the essential bit _first_, then put in the rest for context (and expect the "rest" to get snipped by many respondents). What does Matlab do when you perform an FFT on the data and graph that? -- www.wescottdesign.com
What you are doing is wrong in so many ways.
If your goal is nice looking visualization, just differentiate a signal, 
measure distance between zero crossings and display the periodogram.

Vladimir Vassilevsky
DSP and Mixed Signal Design Consultant
http://www.abvolt.com




Sohaib Afzal wrote:

> Hi > > I am trying to use the goertzel algo to make a real time audio > spectrum analyzer using a PIC controller. I have started out by trying > it on Matlab, but i am getting this problem: > > The high frequency peaks in the power graph are much smaller than the > low frequency peaks. For example, the peak at 60Hz is, say 50 units on > the plot, then the peak at 8000Hz is only 0.023 units or so. Windows > Media Player and other players show the peaks are approximately the > same size throughout the sound file (.wav) that i am testing. > > Please tell me if this algorithm is appropriate for what i am doing, > and if this difference in the magnitudes of peaks for power is usual > or due to some error in coding. > > Thank you for your time.
Sohaib Afzal wrote:
> Hi > > I am trying to use the goertzel algo to make a real time audio > spectrum analyzer using a PIC controller. I have started out by trying > it on Matlab, but i am getting this problem: > > The high frequency peaks in the power graph are much smaller than the > low frequency peaks. For example, the peak at 60Hz is, say 50 units on > the plot, then the peak at 8000Hz is only 0.023 units or so. Windows > Media Player and other players show the peaks are approximately the > same size throughout the sound file (.wav) that i am testing. > > Please tell me if this algorithm is appropriate for what i am doing, > and if this difference in the magnitudes of peaks for power is usual > or due to some error in coding.
The Goertzel algorithm is effectively a stripped Fourier transform that looks at only a few frequencies. It is useful, for example, for detecting DTMF tones. If there is a good way to use it for displaying a broad spectrum, I don't know of it. Jerry -- Engineering is the art of making what you want from things you can get. �����������������������������������������������������������������������
On Thu, 31 Dec 2009 14:41:14 -0500, Jerry Avins wrote:

> Sohaib Afzal wrote: >> Hi >> >> I am trying to use the goertzel algo to make a real time audio spectrum >> analyzer using a PIC controller. I have started out by trying it on >> Matlab, but i am getting this problem: >> >> The high frequency peaks in the power graph are much smaller than the >> low frequency peaks. For example, the peak at 60Hz is, say 50 units on >> the plot, then the peak at 8000Hz is only 0.023 units or so. Windows >> Media Player and other players show the peaks are approximately the >> same size throughout the sound file (.wav) that i am testing. >> >> Please tell me if this algorithm is appropriate for what i am doing, >> and if this difference in the magnitudes of peaks for power is usual or >> due to some error in coding. > > The Goertzel algorithm is effectively a stripped Fourier transform that > looks at only a few frequencies. It is useful, for example, for > detecting DTMF tones. If there is a good way to use it for displaying a > broad spectrum, I don't know of it.
If you were willing to sweep it across the frequency range it could be made to work -- but that's not what the OP says he's doing. Hmm. Need more detail... -- www.wescottdesign.com
I am probably doing everything wrong:

The frequency i am searching for is stored in K. this code runs for K=
120, 260, ...12000
NN is my bin width, x is the signal.

for N=NN:NN:size of x
    for each K
       find coeff, put Q1 and Q2 = 0
       for n=NN consecutive integers, determined using N
            Q0(K) = coeff*Q1(K) - Q2(K) + x(n);
            Q2(K) = Q1(K);
            Q1(K) = Q0(K);
       end;
       pow(K) = Q1(K)*Q1(K) + Q2(K)*Q2(K) - coeff*Q1(K)*Q2(K);

i plot pow. Different plots for different values of N.
Ten lines of code :)

When i did fft, i got large peaks at 2000, 4500, 6000 as well as at
some smaller frequencies. But my code gives only tiny or no peaks
above 1200.

My implementation is probably/definitely wrong, but please also tell
me if it is recommended to use goertzel algo to view the audio
spectrum in real time, or would i be better-off using fft? I think
differentiation etc would be very difficult in the PIC controller.

By "sweeping it across the frequency range" do you mean i need to
apply goertzel to all values of K between 120 and 12000, i.e. also to
121, 122, 123... etc?

Thanks for helping!
Sohaib Afzal wrote:
> I am probably doing everything wrong: > > The frequency i am searching for is stored in K. this code runs for K= > 120, 260, ...12000 > NN is my bin width, x is the signal. > > for N=NN:NN:size of x > for each K > find coeff, put Q1 and Q2 = 0 > for n=NN consecutive integers, determined using N > Q0(K) = coeff*Q1(K) - Q2(K) + x(n); > Q2(K) = Q1(K); > Q1(K) = Q0(K); > end; > pow(K) = Q1(K)*Q1(K) + Q2(K)*Q2(K) - coeff*Q1(K)*Q2(K); > > i plot pow. Different plots for different values of N. > Ten lines of code :) > > When i did fft, i got large peaks at 2000, 4500, 6000 as well as at > some smaller frequencies. But my code gives only tiny or no peaks > above 1200. > > My implementation is probably/definitely wrong, but please also tell > me if it is recommended to use goertzel algo to view the audio > spectrum in real time, or would i be better-off using fft? I think > differentiation etc would be very difficult in the PIC controller. > > By "sweeping it across the frequency range" do you mean i need to > apply goertzel to all values of K between 120 and 12000, i.e. also to > 121, 122, 123... etc?
What do you expect from your code? If you think of an FT as a banl of filters, each tuned to the center frequency of a "bin"m then a Goertzel is a single filter tuned to a single bin. You may be enlightened by http://www.mstarlabs.com/dsp/goertzel/goertzel.html Jerry -- Engineering is the art of making what you want from things you can get. �����������������������������������������������������������������������
Jerry Avins wrote a jumble.

What do you expect from your code? If you think of an FT as a bank of 
filters, each tuned to the center frequency of a "bin", then a Goertzel 
is a single filter tuned to a single bin. You may be enlightened by 
http://www.mstarlabs.com/dsp/goertzel/goertzel.html (Particularly by the 
figure with the green border that shows the response to a wide spectrum.

Jerry
-- 
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������
On Jan 2, 1:07&#2013266080;am, Sohaib Afzal <504...@gmail.com> wrote:
> I am probably doing everything wrong: > > The frequency i am searching for is stored in K. this code runs for K= > 120, 260, ...12000 > NN is my bin width, x is the signal. > > for N=NN:NN:size of x > &#2013266080; &#2013266080; for each K > &#2013266080; &#2013266080; &#2013266080; &#2013266080;find coeff, put Q1 and Q2 = 0 > &#2013266080; &#2013266080; &#2013266080; &#2013266080;for n=NN consecutive integers, determined using N > &#2013266080; &#2013266080; &#2013266080; &#2013266080; &#2013266080; &#2013266080; Q0(K) = coeff*Q1(K) - Q2(K) + x(n); > &#2013266080; &#2013266080; &#2013266080; &#2013266080; &#2013266080; &#2013266080; Q2(K) = Q1(K); > &#2013266080; &#2013266080; &#2013266080; &#2013266080; &#2013266080; &#2013266080; Q1(K) = Q0(K); > &#2013266080; &#2013266080; &#2013266080; &#2013266080;end; > &#2013266080; &#2013266080; &#2013266080; &#2013266080;pow(K) = Q1(K)*Q1(K) + Q2(K)*Q2(K) - coeff*Q1(K)*Q2(K); > > i plot pow. Different plots for different values of N. > Ten lines of code :) > > When i did fft, i got large peaks at 2000, 4500, 6000 as well as at > some smaller frequencies. But my code gives only tiny or no peaks > above 1200. > > My implementation is probably/definitely wrong, but please also tell > me if it is recommended to use goertzel algo to view the audio > spectrum in real time, or would i be better-off using fft? I think > differentiation etc would be very difficult in the PIC controller. > > By "sweeping it across the frequency range" do you mean i need to > apply goertzel to all values of K between 120 and 12000, i.e. also to > 121, 122, 123... etc? > > Thanks for helping!
Sohaib, I am assuming that your are running the Goertzel aganist K for each tone you need, and NN is the number of samples taken from K in consecutive order in NN steps. What I use to give a fairly accurate amplitude measure is (using your notation) y(N:N+NN-1) = sqrt(Q1^2 + Q2^2 - (Q1 * Q2 * coeff))*2/NN For an approximate of power square y, pow(N:N+NN-1) = (Q1^2 + Q2^2 - (Q1 * Q2 * coeff))*4/NN^2. This might work. Also, there's a lot going on where you say for each K find coeff in Matlab lingo, you should have something like k = fix(0.5 + NN*f/Fs); w = (2*pi/NN)*k; cosine = cos(w); coeff = 2 * cosine; CAUTION: in Matlab the use of fix( ) could be important. Maurice Givens
Thanks for your interest.
This is what i am trying to do:
K is an array of the frequencies i am looking for (sorry about the
confusing naming).
K = [120, 260, 660, 990, 1200, 2400, 6000, 10000]; (in Hertz)
I take NN samples from my input signal x, and then for each value of
K, i apply the code i gave above. I use the formula you have given for
magnitude calculation (although i didnt scale with 2/NN)
I do find the coefficients exactly as you have described, (although i
hav named f as K...)

After i get the values for magnitude, y, for NN consecutive samples of
signal x, i plot all the magnitudes y (corresponding to each frequency
K) on the same plot.

Ascii art:
|
|
|
|
|  o
|o |
|| |
|| |
|| |
|| |
|| |
|| | o
|| | | o
|| | | |         o
||_|_|_|_________|
_____________________________o_______________________________o_______
0    1000    2000    3000    4000    5000    6000    7000    8000
9000    10000

The problem is that for low frequencies, y has a very large value, and
for high frequencies, a very low value. I plot this graph for each NN
samples of x, and all the plots have this same problem.

What i expect from my code: plots for each NN consecutive values. The
x axis of the plots has values of frequency, K, and the y axis has
amplitude. The amplitudes should be the same order of magnitude, but
aren't.

When i apply the goertzel function in matlab:
y = goertzel(x, K);
plot(K, abs(y));
then i get the kind of plot i am looking for.

If i succeed in making my own code work, i eventually aim to implement
it in the PIC microcontroller.

....
The sampling frequency is 22050 Hz. I have tried taking various values
of NN. Larger NN, say 900, give an even greater difference between the
values of y corresponding to high and low values of frequency K.
My code works quite well if i generate a signal myself in Matlab (sum
of a few sinusoids of different frequencies). But the problem is
happening when i get the signal from the function wavread (which reads
a .wav file and generates an array containing samples of its
magnitude).

Two questions:
Why does the goertzel function of Matlab not need sampling frequency
as an input?
And why does it give an error when i give it NN samples of x, and a
frequency K greater than NN? Does Goertzel algorithm require the
sample window NN to be larger than the greatest frequency being
searched for?

Thank you all for helping!