Close

An Integer Version of Goertzel Algorithm

A project log for Spectrum Analyser Code

Code for a spectrum analyzer.

agpcooperagp.cooper 10/23/2016 at 10:252 Comments

ATTiny85 Migration

Migrating the test code to ATTiny85 requires some optimisation to suit the limitations of the processor program, RAM and execution speed.

The ATTiny85 has;

SRAM will be the main limitation.

*** Big mistake ***

The 125 kHz specification is the ADC clock speed, not the conversion speed.

It takes 25 cycles for the first conversion and then 13 cycles per conversion.

This is 9.6 kHz!

All is not lost, I can speed (http://www.microsmart.co.za/technical/2014/03/01/advanced-arduino-adc/) clocked his code on an AVR at an average of 20.5 us (48.8 kHz).

Not enough for a 40 kHz ultrasonic tone converter but enough for an audio spectrum analyzer.

Integer Goertzel

Can I avoid floating point maths?

My first attempt a while back failed.

I did not really try to work it out then.

This time I set up variable shift fixed point maths.

Using long (32 bit) failed for low frequencies (i.e. 100 Hz).

Using long long (64 bit) worked for shifts between 15 and 20:

Integer Goertzel Numerical Stability
Frequency 100 Hz
Data Type Magnitude Error
Single 105.0406 0.001%
Double 105.0415 0.000%
Long long fixed point Shift
13 92.559 11.883%
14 100.272 4.541%
15 104.327 0.680%
16 104.329 0.678%
17 104.329 0.678%
18 104.846 0.186%
19 104.846 0.186%
20 104.976 0.062%
21 121.582 -15.747%
22 108.954 -3.725%

Single precision (float) work okay as well.

As long long takes as long as float to process there is no advantage to using integers for this algorithm!

If you are wondering what variable shift fixed point code looks like, here is the critical bit of code:

    coeff=(long long)(2*cos(2*M_PI*freq/SampleFreq)*(1<<shift));
    s_prev=0;
    s_prev2=0;
    for (i=0;i<N;i++) {
      // Goertzel
      s=(window[i]*samples[i]>>shift)+(s_prev*coeff>>shift)-s_prev2;
      s_prev2=s_prev;
      s_prev=s;
    }

Yes, all you do is scale everything up before use, using "*(1<<shift)" and scale any multiplications down using ">>shift".

AlanX

Discussions

agp.cooper wrote 10/23/2016 at 15:43 point

Hi Eric,

You are right with regard to long long and float but the difference is not great for add and multiply but double for divide.

The fixed point maths is straight forward:

int shift=8;

long long a=123*(1<<shift);

long long b=345*(1<<shift);

                 b=a+b;

long long c=(a*b>>shift);

                printf("%lld\n",(c>>shift));

etc.

Regards AlanX

  Are you sure? yes | no

Eric Hertz wrote 10/23/2016 at 15:06 point

Long-long taking longer than float... Seemed surprising 'till I realized this was a Tiny, not a Mega (so doesn't have a mult instruction). 

I wonder how much speed-improvement would be gained by somehow inlining the multiplication procedure, rather than having to jump to a mult-function (and likely multiple lower-bit mult-function-calls?)...
might require assembly :/


I definitely want to compare your integer-goertzel code to mine... but I'll have to wrap my head back into
this stuff first... Might be some time, I'mma link this page in that project's header-file :)

Not sure I understand what you mean by "shifts between..." are you talking about the value of the variable "shift"? Wouldn't having a lower-number of shifts be the equivalent of storing more precision? Or is there a point where even your long-long is overflowing? (is it possible that >> is higher-precidence than *? I always put these things in parens)

  Are you sure? yes | no