DSPRelated.com
Forums

Some Goertzel questions

Started by phuture_project June 23, 2004
Hi everyone,

I've got some questions for you. Please consider i'm not a specialist
and so that some questions may appear very simple for you.

My aim is to determine the magnitude of an analog incoming signal at 9
kHz. This incoming signal is made of several frequencies, let's say
from 0 to 30 kHz. First i digitize this signal and then i treat it
thanks to a DSP (the TMS320C5402). The treatment deals with the
determination of this magnitude. I've been told to use the Goertzel
algorithm rather than a classsic FFT since i only need to know the
magnitude at one given frequency. I studied this algorithm and tried
to use it. Now i've got some questions :

Q1 : To make some tests i send a sine at 9 kHz into the ADC input.
This ADC is an 8 bit one.

- What is the magnitude of this sine, knowing that the sine goes from
0.5 V to 2V? Can we say the magnitude is the square root of(Re^2+Im^2)
? If yes, what is the value of Re and Im?


Q2 : I implemented the Goertzel algorithm but the results look
strange.

- The final result i get, the magnitude squared, varies quite widely.
I thought it must remain almost the same?!
Is Goertzel algorithm really efficient? I told one of my professor
about it and he doubts it can give me good results. Is he right or
wrong?
Is my way of making test good (sending a sine wave at 9 kHz)? If no
what should i'd better do?


Q3 : If i use a classic FFT, how can i get the magnitude at 9 kHz? I
guess that if Fs (sampling frequency) is 140 kHz for example and if i
chose N=256, the FFT will compute the real and imaginary part of the
signal for each frequency between 0 and Fs/2, with a space between
these frequencies of Fs/2N=273.438 Hz. Am i wrong? If i'm right, i
guess the results are saved on a table et so that it's possible to
find in which position the results of the about 9 kHz is. And if i
succeed in finding the real and imaginary parts i can determine the
magnitude. Is this right or wrong?
I looked at the functions given by TI in dsplib. There is one called
RFFT256.asm. I thought i could use it, but i can't see how i can find
the results of the FFT. Are they put in the same table in which the
input samples were saved before the FFT?

As you can see my post is very long. That's why I do apreciate your
help.

Thanks in advance.
A Goertzel filter effectively is an IIR implementation FFT in which
information about uninteresting frequencies is not calculated. The
bandwidth can be adjusted by varying the effective number of samples
that affect the result.

In order to measure any quantity, you need to decide how it is to be
defined in your case. In the case of magnitude, what you have is a
number. The way that number relates to the voltage you want to measure
depends on both hardware and software. The sensitivity of the ADC and
gain before it is the hardware part. What your program does with the
numbers that come from it is the software part. These can be lumped
together into a single calibration step of relating the final numbers to
a known voltage that causes it. A sine that goes from 0.5 V to 2 V has a
DC offset of 1.25 V. That is the amount which, if subtracted from it,
would make it symmetrical about zero. The amplitude (a more specific
term than magnitude*) is half the peak-to-peak value, or .75 V.

The Goertzel algorithm works. You may be having one of two problems with
it that I can think of (and many more that I haven't). It's not unusual
for subtle coding errors to persist in newly written code.
# There can be occasional overflows or other problems that make the
result depend on the particulars of the data stream.
# Your difficulty can also arise from the algorithm working too well. It
is, after all, a selective filter, and slight frequency shifts can
strongly affect its output. Unless the signal that interests you has a
very stable and known frequency, you need to pay attention to the
effective bandwidth of the filter. Whether other filter structures may
be more suitable depends on the specifications of your problem and the 
nature of your signal.

Jerry
__________________________________________
* In 'A*sin(wt)', 'A' is the amplitude.
-- 
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������
phuture_project wrote:

> Hi everyone, > > I've got some questions for you. ...
phuture_project wrote:

> Hi everyone, > > I've got some questions for you. Please consider i'm not a specialist > and so that some questions may appear very simple for you. > > My aim is to determine the magnitude of an analog incoming signal at 9 > kHz. This incoming signal is made of several frequencies, let's say > from 0 to 30 kHz. First i digitize this signal and then i treat it > thanks to a DSP (the TMS320C5402). The treatment deals with the > determination of this magnitude. I've been told to use the Goertzel > algorithm rather than a classsic FFT since i only need to know the > magnitude at one given frequency. I studied this algorithm and tried > to use it. Now i've got some questions : > > Q1 : To make some tests i send a sine at 9 kHz into the ADC input. > This ADC is an 8 bit one. > > - What is the magnitude of this sine, knowing that the sine goes from > 0.5 V to 2V? Can we say the magnitude is the square root of(Re^2+Im^2) > ? If yes, what is the value of Re and Im? >
This depends on the ADC reference voltage. Find out what numeric sequence comes out of the ADC. You should be able to fit that to a sine wave and see the amplitude. The value of Re and Im will add up to the amplitude, but unless you know the phase of your signal you can't say just what each will be (that's why you spend all the extra computing power figuring out the square root of the sum).
> > Q2 : I implemented the Goertzel algorithm but the results look > strange. > > - The final result i get, the magnitude squared, varies quite widely. > I thought it must remain almost the same?! > Is Goertzel algorithm really efficient? I told one of my professor > about it and he doubts it can give me good results. Is he right or > wrong? > Is my way of making test good (sending a sine wave at 9 kHz)? If no > what should i'd better do? >
Yes, for what you're doing the Goertzel algorithm is really efficient, but that doesn't mean that it'll give you good results in your environment. You're doing many fewer computations for the same result, which is good, but it will have limitations which may be giving you grief. The algorithm is an infinite-impulse response filter with poles right on the unit circle, so it's response will grow linearly with an input signal of exactly the tuned frequency. If you don't know what the heck I just said pick up a good book on signal processing and read it. Some limitations that you will see are: 1. If your signal frequency isn't exactly equal to your filter frequency then you may not have any response at the end of your sample. The actual filter response will be equal to it's impulse response convolved by the input window -- you're probably just starting to sample then stopping, which will give you a rectangular window, but all the FFT windowing stuff you know applies. 2. Like any resonant FIR filter this algorithm will be subject to coefficient and data precision issues. Make sure that you have enough accuracy (and remember that a 32-bit floating point only has 24 bits of mantissa). 3. Extraction of the real and imaginary parts is a bit obscure and also subject to roundoff -- are you doing this right? If you use the form mag^2 = Q_1^2 + Q_2^2 - Q_1*Q_2*coeff then you'll _really_ be subject to roundoff errors. I suspect that you're either having roundoff problems or you're slightly off frequency -- how exact is your 9kHz, and how long are you collecting data? Your effective bandwidth will be approximately equal to the inverse of your filter duration, but to reject 30Hz (for instance) you'll need to have a long window -- if you're collecting for just 111us your 30Hz content may screw you up, if you collect for 33ms then your filter will only be 30Hz wide. Your test is good for starters -- it's thoroughly grounded in reality, so it tests everything. Once you get it working (any algorithm) with just a 9kHz signal you'll want to throw in other stuff to see how well your algorithm works.
> > Q3 : If i use a classic FFT, how can i get the magnitude at 9 kHz? I > guess that if Fs (sampling frequency) is 140 kHz for example and if i > chose N=256, the FFT will compute the real and imaginary part of the > signal for each frequency between 0 and Fs/2, with a space between > these frequencies of Fs/2N=273.438 Hz. Am i wrong? If i'm right, i > guess the results are saved on a table et so that it's possible to > find in which position the results of the about 9 kHz is. And if i > succeed in finding the real and imaginary parts i can determine the > magnitude. Is this right or wrong? > I looked at the functions given by TI in dsplib. There is one called > RFFT256.asm. I thought i could use it, but i can't see how i can find > the results of the FFT. Are they put in the same table in which the > input samples were saved before the FFT? > > As you can see my post is very long. That's why I do apreciate your > help. > > Thanks in advance.
To get the FFT response you figure the bins that correspond to your signal and look at their amplitude. -- Tim Wescott Wescott Design Services http://www.wescottdesign.com
"Tim Wescott" <tim@wescottnospamdesign.com> wrote in message
news:10dj9ct16ps2jf6@corp.supernews.com...
> phuture_project wrote: > > 2. Like any resonant FIR filter this algorithm will be subject to > coefficient and data precision issues. Make sure that you have enough > accuracy (and remember that a 32-bit floating point only has 24 bits of > mantissa).
Just a clarification/correction, 32-bit IEEE floating-point has 25 bits of mantissa due to the "hidden bit". (This assumes you count the sign bit as part of the mantissa, a very reasonable assumption since this is routinely done with fixed-point.)
Jon Harris wrote:

> "Tim Wescott" <tim@wescottnospamdesign.com> wrote in message > news:10dj9ct16ps2jf6@corp.supernews.com... > >>phuture_project wrote: >> >>2. Like any resonant FIR filter this algorithm will be subject to >>coefficient and data precision issues. Make sure that you have enough >>accuracy (and remember that a 32-bit floating point only has 24 bits of >>mantissa). > > > Just a clarification/correction, 32-bit IEEE floating-point has 25 bits of > mantissa due to the "hidden bit". (This assumes you count the sign bit as part > of the mantissa, a very reasonable assumption since this is routinely done with > fixed-point.) > >
Jeez -- and I thought _I_ was a picky so-and-so. Good point. I'll remember it. It doesn't (much) change the fact that 32-bit floating point can be disappointingly imprecise. -- Tim Wescott Wescott Design Services http://www.wescottdesign.com
"Tim Wescott" <tim@wescottnospamdesign.com> wrote in message
news:10djhl69rs3hlbb@corp.supernews.com...
> Jon Harris wrote: > > > "Tim Wescott" <tim@wescottnospamdesign.com> wrote in message > > news:10dj9ct16ps2jf6@corp.supernews.com... > > > >>2. Like any resonant FIR filter this algorithm will be subject to > >>coefficient and data precision issues. Make sure that you have enough > >>accuracy (and remember that a 32-bit floating point only has 24 bits of > >>mantissa). > > > > > > Just a clarification/correction, 32-bit IEEE floating-point has 25 bits of > > mantissa due to the "hidden bit". (This assumes you count the sign bit as
part
> > of the mantissa, a very reasonable assumption since this is routinely done
with
> > fixed-point.) > > > Jeez -- and I thought _I_ was a picky so-and-so. > > Good point. I'll remember it. It doesn't (much) change the fact that > 32-bit floating point can be disappointingly imprecise.
This misconception happens quite often, so I guess I've made it my personal mission to correct it. :-) There is some debate in the audio community about the relative merits of fixed- vs. floating-point processing, and the fixed-point proponents always seem to short-change floating-point by one bit, which irks me! Also, it means that IEEE floating point _always_ beats 24-bit fixed point by at least 1 bit! At any rate, I've found 32-bit floating point implemented with 40-bit accumulators (SHARC) to be quite adequate for the audio tasks I've been using it for. I would imagine that it would be more than enough for a simple Goertzel algorithm on audio data (which is probably only has 16 or at most 20 useful bits to begin with).
Tim Wescott wrote:

   ...

> The value of Re and Im will add up to the amplitude ...
I'm picky too. Do you want that to stand? ... Jerry -- Engineering is the art of making what you want from things you can get. &#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;
Jerry Avins wrote:

> Tim Wescott wrote: > > ... > >> The value of Re and Im will add up to the amplitude ... > > > I'm picky too. Do you want that to stand? > > ... > > Jerry
No. How about "You can find the amplitude from Re and Im". Better? -- Tim Wescott Wescott Design Services http://www.wescottdesign.com
Jon Harris wrote:

> "Tim Wescott" <tim@wescottnospamdesign.com> wrote in message > news:10djhl69rs3hlbb@corp.supernews.com... > >>Jon Harris wrote: >> >> >>>"Tim Wescott" <tim@wescottnospamdesign.com> wrote in message >>>news:10dj9ct16ps2jf6@corp.supernews.com... >>> >>> >>>>2. Like any resonant FIR filter this algorithm will be subject to >>>>coefficient and data precision issues. Make sure that you have enough >>>>accuracy (and remember that a 32-bit floating point only has 24 bits of >>>>mantissa). >>> >>> >>>Just a clarification/correction, 32-bit IEEE floating-point has 25 bits of >>>mantissa due to the "hidden bit". (This assumes you count the sign bit as > > part > >>>of the mantissa, a very reasonable assumption since this is routinely done > > with > >>>fixed-point.) >>> >> >>Jeez -- and I thought _I_ was a picky so-and-so. >> >>Good point. I'll remember it. It doesn't (much) change the fact that >>32-bit floating point can be disappointingly imprecise. > > > This misconception happens quite often, so I guess I've made it my personal > mission to correct it. :-) There is some debate in the audio community about > the relative merits of fixed- vs. floating-point processing, and the fixed-point > proponents always seem to short-change floating-point by one bit, which irks me! > Also, it means that IEEE floating point _always_ beats 24-bit fixed point by at > least 1 bit! > > At any rate, I've found 32-bit floating point implemented with 40-bit > accumulators (SHARC) to be quite adequate for the audio tasks I've been using it > for. I would imagine that it would be more than enough for a simple Goertzel > algorithm on audio data (which is probably only has 16 or at most 20 useful bits > to begin with). > >
Probably because you're doing all FIR filters. IIR filters are very sensitive to truncation, particularly if they are casually written (such as the Goertzel examples that I've found that use the controller canonical form). Thinking 24 (ok, ok -- 25) bits is enough is a common pitfall, and one I've dug myself out of numerous times. -- Tim Wescott Wescott Design Services http://www.wescottdesign.com
Tim Wescott wrote:

> Jerry Avins wrote: > >> Tim Wescott wrote: >> >> ... >> >>> The value of Re and Im will add up to the amplitude ... >> >> >> >> I'm picky too. Do you want that to stand? >> >> ... >> >> Jerry > > No. > > How about "You can find the amplitude from Re and Im". Better?
Sure. Jerry -- Engineering is the art of making what you want from things you can get. &#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;