# Single pole IIR filter, fixed point design

We want to do a fixed point implementation of the single pole IIR filter:

$$y[n] = a\cdot x[n] + (1 - a)\cdot y[n-1] \quad ;\qquad 0

What are the main design considerations to think about ? In particular, is the fixed point design more challenging for smaller a or larger a ?

• Related question: dsprelated.com/thread/10839/… – Cedron Dawg Apr 7 at 15:03
• likely $0< a \ll 1$ so turn that difference equation into $$y[n] = y[n-1] + a \cdot (x[n] - y[n-1])$$ and use the full precision of $y[n-1]$ in the stand-alone term to implement your fraction saving. only where $y[n-1]$ is subtracted from $x[n]$ is where you truncate (that is round toward $-\infty$) bits so that both are single precision and where you multiply by $a$. – robert bristow-johnson Apr 10 at 1:36

One thing to consider when implementing an IIR filter, whatever the order, is quantization and limit cycles.

Let me show you with a quick example with your original filter

$$y[n] = a*x[n]+(1-a)*y[n-1]$$

Let a = 0.005 and say that we use 16-bit signed coefficients.

$$a_{fixedpoint} = a * 32768 = 164$$

Let's assume that the input and output are 16-bit integers. They will go from -32768 to 32767. The previous output $$y[n-1]$$ is also a 16-bit integer in this scenario. The implementation would look like this, with temp being a 32-bit integer. All other variables are 16-bit integers.

$$temp = a_{fixedpoint}*x[n] + (32768 - a_{fixedpoint})*y[n-1]$$

$$y[n] = round(temp/32768)$$

Let's try this implementation with a step response

It looks good, right? Perhaps we should zoom in.

We don't reach 32767 even though the gain of the original filter is 1. This is called a limit cycle. For order-1 IIR filters, limit cycles take the form of a DC offset that you can't get rid of. For IIR filters with complex poles, limit cycles take the form of a small oscillation that you can't get rid of. The problem is caused by the fact that $$y[n-1]$$ is stored on 16 bits. If we increase the resolution of $$y[n-1]$$ to 32 bits, it would solve the problem. Another solution would be to use fraction-saving.

If the input would go from 32767 to 0, we would have the same problem, the output would get close to 0 but would not actually reach 0.

Increased resolution

If instead of storing the output on 16 bits, we store it on 24 bits and name this variable acc (increased resolution of y[n]).

$$temp = a_{fixedpoint}*x[n]*256 + (32768 - a_{fixedpoint})*acc$$

$$acc = temp/32768$$

$$y[n] = round(acc/(256))$$

As we can see, with the increased resolution we can reach 32767 even though the coefficients are still on 16 bits. However the price to pay is having intermediate variables with greater resolution, the calculations might be slower.

fraction-saving

Finally there's a clever approach to minimize the added resolution needed for an IIR, especially when the poles are close to the unit circle. It is called fraction-saving. The idea is to memorize the round-off error of the previous sample and apply it to the next calculation to reduce the quantization effects. An additional variable is needed called error, but you don't need to increase the resolution of the previous output. So if x is represented on 16 bits, then y[n] can be stored on 16 bits too. The error variable would be represented with 16 bits too in this example. As Robert Bristow-Johnson pointed out, fraction-saving is a form of noise-shaping

$$temp = a_{fixedpoint}*x[n] + (32768 - a_{fixedpoint})*y[n-1] - error[n-1]$$

$$y[n] = round(temp/32768)$$

$$error[n] = y[n]*32768 - temp$$

• this really isn't a limit cycle as the predicted gain is (1- a) so you are just seeing that (put in z=1 in the transfer function to see the gain). The normalized exponential averager is simply multiplied by the inverse of this. Limit cycles are an issue in fixed point systems as they cause sustained oscillations. – Dan Boschen Apr 7 at 13:58
• The transfer function is a*z/(z - (1-a)) right? If I replace z by 1, I get a DC gain of 1. Perhaps I made a mistake. – Ben Apr 7 at 14:04
• 1-a = 0.995. The gain I see is 32668/32767 = 0.997 – Ben Apr 7 at 14:33
• Sorry! You're right I missed the gain factor that is already in there. Good demo and look forward to seeing your noise shaping solution. – Would be even better if you could show the case of a limit cycle oscillation since such gain differences are usually inconsequential. – Dan Boschen Apr 7 at 15:16
• yea!! FRACTION SAVING. a.k.a. first-order noise shaping with a zero at $z=1$, or "DC", which means that, at DC, there is no error and infinite signal-to-noise ratio. i am a big proponent of fraction saving for fixed-point FIR or IIR filters. – robert bristow-johnson Apr 10 at 1:38

In particular, is the fixed point design more challenging for smaller a or larger a ?

Smaller is worse: the closer $$a$$ gets to $$0$$, the closer the pol movess to the unit circle and the more time domain ringing you have.

What are the main design considerations to think about ?

A first order low-pass IIR filter is relatively benign. The sum absolute sum of the impulse response is unity so you can't clip the output. One thing to watch out for are potentially limit cycles. This can be controlled by the way you round: "round towards zero" is the safest option.

It used to be that multiplication was a lot more expensive than addition or bitwise operations. For integer (and fixed point) implementation, this can be taken exploited like this:

When $$a=1/2$$ an extremely efficient implementation can be made with a single add and a single shift.

$$y[n] = (x[n] + y[n-1]) >> 1$$

Similarly, fractions with a power of two in the denominator can also be done efficiently. For instance, $$a = 7/8$$ can be done as:

$$y[n] = ((x[n]<<3) - x[n] + y[n-1]) >> 3$$

There are many other possibilities as well. I tend to go with the $$a=1/2$$ option a lot.

Obviously, how much headroom you have is a major consideration, so doing something like $$a = 255/256$$ would need 8 bits of headroom, whereas having fixed point fractions doesn't suffer from that, but it costs you a multiplication or two. Note the equation can also be rewritten as:

$$y[n] = y[n-1] + a * ( x[n] - y[n-1] )$$

This saves a multiplication at the cost of a subtraction.

• thanks. I have edited the equation, I just mis-wrote it. I understand your insights, and that a power of 1/2 makes sense. Is there anything else to be said from maybe a more conceptual level : Say, we could choose any value of a (not just powers of 1/2): is smaller "a" better or larger "a" better from fixed point design ? e.g., a --> 0 implies the pole approaches the unit circle, which I guess is moving the system towards unstability. Does that matter ? – voy82 Apr 7 at 0:27
• @voy82 Somebody else will have to pipe in about the poles and ROC analysis, but since you are defining $0<a<1$ (I assume real) this is going to be a smoothing transform and you don't have to worry about divergence. I have never really played around with complex valued coefficients with this. You may find this article of mine interesting as well: dsprelated.com/showarticle/896.php – Cedron Dawg Apr 7 at 0:36
• @voy82 If $a=0$ then $y[n]=y[0]$ for all $n$. If $a=1$ then $y[n]=x[n]$ for all $n$. Therefore, the parameter $a$ represents how "reactive" the transformed signal is to the input signal. Like longer or shorter moving averages. I don't think there are any fixed point considerations with that, but I'm not really a fixed point jockey. I do use this technique, and the Differ referenced in my article a lot. Heavy smoothing at a cheap price. – Cedron Dawg Apr 7 at 0:46
• Fixed-point implementation can suffer from limit cycles. For an order-1 IIR filter, it means that sometimes the output can never fully reach 0 or 1. A value will be stuck in the accumulator. Noise-shaping can fix this. – Ben Apr 7 at 0:58
• @Ben Thanks Ben. Maybe you should elaborate with a separate answer. It seems that is the kind of think the OP is fishing for. – Cedron Dawg Apr 7 at 1:18

Lyons and Yates later did an IEEE Sig Proc article about DC blocking filters where this was one of the topics.

// let's say sizeof(short) = 2 (16 bits) and sizeof(long) = 4 (32 bits)
short x[], y[];
long acc, A, prev_x, prev_y;
double pole;
unsigned long n, num_samples;
pole = 0.9999;
A = (long)(32768.0*(1.0 - pole));
acc = 0;
prev_x = 0;
prev_y = 0;
for (n=0; n<num_samples; n++)
{
acc   -= prev_x;
prev_x = (long)x[n]<<15;
acc   += prev_x;
acc   -= A*prev_y;
prev_y = acc>>15;               // quantization happens here
y[n]   = (short)prev_y;
// acc has y[n] in upper 17 bits and -e[n] in lower 15 bits
}