Signal strength at some fixed frequency

Can someone point me to a simple algorithm to compute
the signal strength at some frequency? Say I've got 1024
samples, 16 bit signed values, at 44100 hz, and I want
to know how much of 440 hz is present in the signal.

Thanks--
Dave

 

You can do it by correlation -- multiply your signal by sine and cosine
waves at 440Hz, add the results, then find the RMS of those two numbers.

Or you can do it using the Goertzel algorithm -- it's a favorite of CS
types, it does the job, but it has lots of subtle details to it that
most cookbook solutions overlook. OTOH, there are cookbook solutions
out there...

--

Tim Wescott
Wescott Design Services
http://www.wescottdesign.com

Posting from Google? See http://cfaj.freeshell.org/google/

"Applied Control Theory for Embedded Systems" came out in April.
See details at http://www.wescottdesign.com/actfes/actfes.html

 

With such small number of samples, you just get the amount of signal
in the neighbourhood of 440 Hz, perhaps sufficient for some "real time
analyzer" but definitively not for instrument tuning.

Paul

 

Compute the correlation with a 440 Hz sine wave.

 

I tried Goertzel and the results look promising. I want to use
integer math but the Goertzel sample code I got from
Wikipedia uses doubles, and the numbers balloon up. One
tunable was the number of samples used in the calculation.

Thanks--
Dave

 

It may not be that simple. You will see truncation and leakage effect.
Then you start to think about the window, then you need more than one
spectrum line to calculate the overall rms.

If he wants to achieve an estimation with high accuracy, I have not
seen a simple answer yet.

D.S.

 

Tim Wescott has the right idea. It's also called synchronous
demodulation. It simple and works quite well.

gm

 

I recall an article in Circuit Cellar where they implemented the
Goertzel on an 8 bit micro without using floating point math. Since the
number of frequencies they were looking for was limited (as in your
case) they "pre-calculated" some of the operations needed to compute
the magnitude that you're looking for. They also did some
approximations to simplify the math. Search circuit cellars website. It
may work for your needs.

-Ed Rivas

 

I did find a code snippet on CC related to a design contest, I think
it was a complete phone answering machine. It used
Goertzel for the tone recognition, probably for data entry.
There was a small function referred to but it wasn't complete.
It was obviously integer math but seemed iffy. One section
showed ">>15" and another ">> AMP_BITS". But the code
to set up the coefficients and the definition of AMP_BITS wasn't
included.

I figured I'd just play around with it when I next get a chance.
Thanks--

-Dave

 

Which makes me wonder: does anybody have recomendations on an allround
DSP cookbook for embedded programmers, especially the
non-mathematicians? I'm thinking of something describing the dirty
details of sampling theory (Tim Wescott's "what Nyquist didn't say"!),
introduction to FIR/IIR filters and other alternatives, with practical
guidlines on implemententation, also on smaller CPU's without FPU
hardware, and 'tricks' from the more experienced in the field. The
'hackers delight' version about digital signal processing ?

 

I found "Digital Signal Processing - A Practical Approach", by
Ifeachor and Jervis fairly easy to digest. It does have math, but also
code examples, so if you can't follow the math, you can at least copy
the code

 

Here is one form of Goertzel, doesn't seem to hard to integerize it.

float goertzel(float *x, int N, float frequency, int samplerate) {
float Skn, Skn1, Skn2;
Skn = Skn1 = Skn2 = 0;

for (int i=0; i<N; i++) {
Skn2 = Skn1;
Skn1 = Skn;
Skn = 2*cos(2*PI*frequency/samplerate)*Skn1 - Skn2 + x;
}

float WNk = exp(-2*PI*frequency/samplerate); // this one ignores
complex stuff
return (Skn - WNk*Skn1);
}

 

In the snippet below, Pi = 3.14, freq = 440, samplerate = 44100. You
can calculate cos of this number and use it as a constant. I think this
was the goal of the article I read. Not sure if it was the same one as
you are referring to, I don't think it was though.

Did some searching and I don't think they have the particular article
I'm referring to available on the web. The issue is 182 (Digital
Decoding Simplified: Sequential Exact-Frequency Goertzel Algorithm, by
Eric Kiser, p. 22)

good luck, post up if you get something working.

-Ed

 

I just ran the test program on an 50MHz ARM Cortex (LM3S811). It takes
about 10K bytes of code including floating points and complete in
approx 20msec. It is certainly doable with a decent micro. I would
not spent too much time trying to convert it to integer.

Test case:

R=8000Hz F=941Hz N=200

 

Well it's a bit more complicated, I need to determine signal presence of
almost 100 different frequencies, all at the same time, not just one.
So I think integer arithmetic is in order...

-Dave

 

More like 40mesc.

Within what time frame? What Sampling rates and processing bins? For
example, the 50MHz ARM can process a 10Khz signals in 3 to 4 seconds.
Of course, you can hook-up 100 of them to get results in 40 msec

If you really need it.

 

Actually now that I think about it I can probably get by with
looking at just 24 frequencies at a time. But I need to know
about 10 times per second whether each frequency is present
in the signal.

-Dave

 

For touch-tone you want 8 frequencies, updated in real time. For 200 samples
you've got 25 ms for all eight.

That is abysmal, 0.040 sec * 50 MHz / 200 data points = 10000 cycles per
data point! The algorithm requires an add, a subract and a multiply per
sample. Even for software FP thats aweful.

If you have to go to integer, which I think you must in this case, then even
a 32.32 fix-point solution would work well on a ARM. A 16.16 fix-point might
be quicker but you'd need to analyse the rounding errors. 100
cycles/sample/freq, absolute tops.

The problem I see is that if you have 100 frequencies and that is much
smaller then the FFT then you much be using lots of samples.

Peter

 

Yes, you are right. This is not optimized. I will take the cos() out
of the loop and redo it for 24 frequencies, as required by the OP.
Anyone want to post a integer version?

 

It would need approx. 120MHz for floating points and 80MHz for
integers. Can you sample a subset and do the rest if necessary. For
DTMF, you only need to decode 8 for first harmonics and another 8 for
second harmonics (to confirm pure tones).