Iowa Hills Software Digital and Analog Filters
IIR Filter Design
Home
IIR Filter Design Using the Bilinear Transform
This page derives the equations used to generate IIR filter coefficients for
low pass, high pass, band pass, and notch filters.
This algorithm
generates second order IIR filter coefficients, not Nth order coefficients.
Direct computation of second order coefficients is preferred for three reasons.
First, we need 2nd order coefficients to implement biquad sections, the preferred, and only practical way to implement an IIR filter.
Second, this approach allows us to use s domain properties, such as the damping coefficient zeta, and the zero locations on the imaginary axis, to properly order the biquad sections. Proper ordering, and proper grouping of the poles and zeros, is critical for numerical reasons.
Third, since we don't start with an Nth order polynomial, we don't need a root solver to get the biquad coefficients we need. This is quite nice for most applications.
This method does however, generate 4th order terms for band pass and notch filters, but 4th order roots are easy to obtain.
The source code for this algorithm is given in part here. The entire algorithm, which includes the necessary 4th order root solver is given in our Code Kit for IIR and FIR Filters.
The Bilinear Transform Design Equations
We start with H(s), redefine the variable s using the bilinear transform, and do the algebra. Our goal is to form H(z) as given above and pick off the various coefficients for use in the IIR filter.
Low Pass Filter Equations
These 2 equations are the desired results, where we have set a0 = 1. Now
all that is needed is to plug in the A’s, B’s .. F’s from the original H(s).
The
value of T is determined by Omega, the desired cutoff frequency for the IIR
filter, which is in terms of Nyquist. For example, if the sampling frequency Fs
is 20 kHz, and the desired corner frequency is 2.0 kHz, then
Ω = 0.2. (2kHz / (Fs/2) = 0.2)
H(z) reduces to the following for all pole filters.
When doing low pass and high pass filters, we get a 2nd order section in H(z)
for each 2nd order section in H(s). You will see however that band pass and notch filters generate a forth order section in H(z) for each
2nd order low pass term of H(s).
If working with odd order polynomials, which have 1 pole sections, do not make
the mistake of assuming that the 2nd order H(z) reduces to the 1st order H(z)
because A = D = 0, it does not. If A and D are set to zero in the 2nd order
H(z), the filter will work correctly, but you will have implemented a single
pole in a 2nd order section. Use the 1 pole H(z) equation for the odd pole.
Pole Ordering
An important aspect to filter design, whether digital or analog, is to
place the 2nd order sections in the correct order.
In the s domain, a 2nd order equation is routinely written as:
H(s) = s^{2 }+ 2ζω s + ω^{2} where ζ is referred to as the
damping factor.
As an example, here are the 2nd order sections for a 6 pole Butterworth. In this
case, ω = 1 for all three.
s^{2} + 0.5176380902 s + 1.0 ζ = 0.258
s^{2} + 1.4142135624 s + 1.0 ζ = 0.707
s^{2} + 1.9318516526 s + 1.0 ζ = 0.965
In this plot, we show the frequency response for these sections and their pole
locations on the s plane.
If building an opamp filter, we order the sections by zeta,
placing the section with ζ = 0.965 first as a way to minimize the voltage swing
of the succeeding sections. This is simply good design practice because it
reduces the chances for the first stage to saturate.
Similarly, when designing an IIR filter to be used on a fixed point processor,
we want to keep the math values minimized so we can prevent register overflow.
To do this, we order the sections by zeta, using the section with the largest
zeta first.
If the pole count is odd, the odd pole resides on the real axis, so its zeta is
essentially 1.0, and should be used first.
Zero Ordering
The Inverse Chebyshev and Elliptic also have s plane zeros. From a
mathematical perspective, the poles and zeros can be grouped together in the 2nd
order sections in any manner desired, but we want to arrange the zeros in a way
that minimizes the peak math values. This can be done by ordering the zeros so
that the inner most zero is used first.
This is a pole zero plot for a 6 pole Inverse Chebyshev.
Here is an example of how incorrect pole zero ordering can affect the peak math values generated by the filter.
Peak Math Values for a 10 Pole Inv Chebyshev Low Pass | |
Correct Ordering | 127.0 |
Zeros Reversed | 198.8 |
Poles Reversed | 513.0 |
Both Reversed | 760.0 |
The following gives the derivations for high pass, band pass, and notch filters.
High Pass Filter Equations
In the s domain, to convert a low pass to a high pass, s becomes 1/s.
Band Pass Filter Equations
You will notice that in this derivation, we set ω =
1. We do this because H(s) is normalized
to ω = 1 and we need to maintain that. Also, because of the frequency warping
occurring here, the definition for Q used here (Q = ω / BW)
works very poorly in the z domain, so we correct this Q calculation with the
code given at the bottom of the page.
Rather than show the expanded view of this equation, we simply give the 10 coefficients. As noted above, this 4th order H(z) would almost certainly need to be factored if implemented in fixed point.
Notch Filter Equations
Q for Bandpass and Notch Filters
Because of the frequency warping that occurs when converting the s domain
polynomials to the z domain, the usual definition for the Q of a bandpass or
notch filter must be modified. Here is the c code we use to correct Q. We
derived this by plotting the bandwidth error as a function of OmegaC and doing a
curve fit.
Without this correction, the resulting filter bandwidth is a strong function of
the filter's center frequency. The filter will be wider than desired at low
center frequencies, and much narrower than desired at high center frequencies.
The bandwidth, BW, like OmegaC, is in terms of Nyquist.
Q = 1.0 + OmegaC; //
OmegaC max = 1.0
if(Q > 1.95)Q = 1.95; // Q must be < 2
for the next line.
Q = 0.8 * tan(Q * Pi / 4); // This is the correction factor.
Q = OmegaC / BW / Q; // This is the corrected Q.
All Pass Filters
In order to design an IIR all pass filter, one starts by adding zeroes to an all
pole transfer function in the usual way, as shown here. Then use the
coefficients A - F in the IIR Low Pass equations given above. Although any all
pole prototype can be used, the Gauss polynomial yields the best results because
of its linear phase. Adjust Omega, in the calculation for T, to achieve the
desired delay time.
Source Code
This is the source code in C for
the equations above.
This Code Kit contains the entire algorithm, including the 4th order root solver required for the band pass and notch filters..
The required 2nd order coefficients for H(s), A, B, C, D, E, and F can be obtained from this Filter Polynomial Program.
For an implementation of the algorithm given here, see our free IIR Filter Designer.
6 Pole Butterworth Low Pass Filter Example
Copyright 2013 Iowa Hills Software