Signal strength at some fixed frequency

Ico wrote:

> 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

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 :)

David Ashley wrote:
> rivas.ed wrote:
> > 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.

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[i];
    }

    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

> 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[i];
>     }
>
>     float WNk = exp(-2*PI*frequency/samplerate); // this one ignores
> complex stuff
>     return (Skn - WNk*Skn1);
> }

rivas.ed wrote:
> 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.

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

>
> 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
>
> > 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[i];
> >     }
> >
> >     float WNk = exp(-2*PI*frequency/samplerate); // this one ignores
> > complex stuff
> >     return (Skn - WNk*Skn1);
> > }

linnix wrote:
> 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.

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

-- 
David Ashley                http://www.xdr.com/dash
Embedded linux, device drivers, system architecture

David Ashley wrote:
> linnix wrote:
> > 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.

More like 40mesc.

> > It is certainly doable with a decent micro.  I would
> > not spent too much time trying to convert it to integer.
>
> 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.

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

> So I think integer arithmetic is in order...

If you really need it.

>
> -Dave
>
> --
> David Ashley                http://www.xdr.com/dash
> Embedded linux, device drivers, system architecture

linnix wrote:
>>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.
> 
> 
> 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

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


-- 
David Ashley                http://www.xdr.com/dash
Embedded linux, device drivers, system architecture

"linnix" <me@linnix.info-for.us> wrote in message
news:1161310441.101062.223270@h48g2000cwc.googlegroups.com...
>
> David Ashley wrote:
> > linnix wrote:
> > > 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.
>
> More like 40mesc.

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.

> > > It is certainly doable with a decent micro.  I would
> > > not spent too much time trying to convert it to integer.
> >
> > 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.
>
> 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
>
> > So I think integer arithmetic is in order...
>
> If you really need it.

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

Peter Dickerson wrote:
> "linnix" <me@linnix.info-for.us> wrote in message
> news:1161310441.101062.223270@h48g2000cwc.googlegroups.com...
> >
> > David Ashley wrote:
> > > linnix wrote:
> > > > 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.
> >
> > More like 40mesc.
>
> 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.

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?

David Ashley wrote:
> linnix wrote:
> >>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.
> >
> >
> > 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
>
> 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.

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).