implementing exponential moving average filter

I'm trying to apply an exponential moving average filter to an analog input.

The formula for an EMA filter is as follows:

value = measurementalpha + previous value(1-alpha)

where alpha is some number between 0 and 1.

Because I'd like to avoid floating value math, I've implemented it as shown below, and it works quite well. It allows me to choose alpha values of integers from 1-10.

//  GAUGE DATA 2 UPDATE FUNCTION  //
int gauge2 () {
  int alpha_2 = 1;                                         // smoothing factor (between 0-10)
  A8_raw = analogRead(A8);                              // read analog values from input_8
  A8_sm = (A8_raw*alpha_2 + A8_sm*(10-alpha_2))/10;    // calculate smoothed value (EMA filter)
  int  angle = map( A8_sm, 0, 1024, 0, STEPS2);             // calculate angle of motor
  return angle;                                             // return angle of motor
}

The limitation in the code above is that I can't fine tune the filter with only ten possible values. I'd prefer to chose an integer from 1-100, but when I change the code to allow that (see below), the A8_sm value begins to swing wildly, if A8_raw gets above about 325.

//  GAUGE DATA 2 UPDATE FUNCTION  //
int gauge2 () {
  int alpha_2 = 10;                                         // smoothing factor (between 0-100)
  A8_raw = analogRead(A8);                              // read analog values from input_8
  A8_sm = (A8_raw*alpha_2 + A8_sm*(100-alpha_2))/100;    // calculate smoothed value (EMA filter)
  int  angle = map( A8_sm, 0, 1024, 0, STEPS2);             // calculate angle of motor
  return angle;                                             // return angle of motor
}

What's going on here? This function is used in a code that also has an interrupt every 200 microseconds. My theory is that the division by 100 is taking so long that the interrupt is causing problems with this function, and somehow, division by 10 doesn't take as long... But I really don't know, I'm still pretty new to arduino and weak in programming.

Maybe consider the limitations of sixteen bit signed arithmetic?

AWOL:
Maybe consider the limitations of sixteen bit signed arithmetic?

Perhaps you know of a good place for a beginner to learn about these limitations?

Consider this line.

  A8_sm = (A8_raw*alpha_2 + A8_sm*(100-alpha_2))/100; 
}

Imagine A8_sm has the value 400.
400 * 90 . . .doesn't fit in a sixteen bit int.
Do the arithmetic in long (32 bit) arithmetic.

The majority of Arduinos use a 16-bit integer with a range of

-32,768 : 32767

EDIT:
Something else to keep in mind is that the results of multiplying two integers occupies the sum of the number of bits used to represent those integers. And since integers come in 8/16/32 and sometimes 64 bit flavors ...

AWOL:
Imagine A8_sm has the value 400.
400 * 90 . . .doesn't fit in a sixteen bit int.

Thank you! This makes perfect sense now!

It is much better to use floating point math with that type of filter.

With integers you have rounding problems that often lead to instability.

jremington:
It is much better to use floating point math with that type of filter.

With integers you have rounding problems that often lead to instability.

Maybe I'm being a little too paranoid about the calculation overhead of floating point math?

neongreen:
Maybe I'm being a little too paranoid about the calculation overhead of floating point math?

I use a lag filter like that in a phase lock loop that is update once a second so there is no problem. How often are you doing the calculation? That is what matters.

Hi,

If you are doing purely integer math then 100 is a poor choice for division. You should use either 64 or 128 as then you can use a shift operation which should be faster than any division.

As to the type of integer, if you are dealing with all positive integers then you may be able to use unsigned int instead of just int. If any of your calculations go over the limit of 64k though then you will have to resort to using unsigned long.

stowite, this function just sits in the main loop with another function that is basically the same. There are also a couple other functions on timers in the main loop like display update(called at 5Hz), a GPS data parser (also called at 5Hz), and then a stepper motor update function (called at 5000Hz). So I'm not sure how often it runs, but it is run anytime the other functions aren't being run.

neongreen:
Maybe I'm being a little too paranoid about the calculation overhead of floating point math?

Whenever I want a really low overhead filter implementation I use alpha = 2**(-k). By which I mean two raised to the power of negative k.

So the basic filter state equation: x = (1-alpha)x + (alpha) u

becomes: x = ( 1-2**(-k) )x + 2**(-k) u

Right now that doesn't look all that simple, but multiply both sides by 2**k and it becomes:

(2k)x = (2k)x - x + u

Now instead of taking "x" as our state variable we just use the multiple instead. Say for example let w=(2**k)x then the state equations become:

w = w - x + u
x = w >> k

Believe it or not, there is your exponential filter without a divide or even a multiplication in sight. :slight_smile:

Just in case you need a concrete example, let k=4, so the state variable is 16 times the smoothed output. Let's call it smoothed16. The state equations for this particular case would be:

smoothed16 = smoothed16 - smoothed + raw_input
smoothed = smoothed16 >> 4

I know that's probably getting a bit "anal" about reducing mpu cycles, but it does work.

Hey stuart0

Many thanks for your solution, works like a charm for my problem.
You're awesome!!!

Thanks

i’ve done leaky integration using fixed point arithmetic, similar to above but provides a little more precision in terms of the filtering factor (your alpha)

a new sample, S, is added to an existing average, S, at some rate, R

A += (S - A) / R

since you would like greater accuracy and avoid floats, S and A can be represented using fixed point values where the binary point, where the bit to its right represents 1/2, can be arbitrary.

For example, S can be shifted left by some number of bits, N. This results in an average which is shifted, multiplied by 2N. An integer result (binary pt at zero) is obtained by shifting the average right N bits. Both S and A need to be an appropriate size integer (i.e. 8, 16, 32, 64) to avoid overflow.

S = analogRead(somePin) << N;
A += (S - A) / R
a = A >> N

the filtering coefficient is R / 2N

the filtering coefficient is R / 2^N

It’s better to write that as " R / 2N " - the XOR operator makes the first form a little confusing.