DSPRelated.com
Forums

Simulating higher order IIR filters using cascaded low order sections

Started by Unknown November 27, 2012
I'm trying to make a tool that implements various types of IIR filters up to a certain maximum filter order (8 in my case). I've been learning about the bilinear transform and how to use it to turn an analog transfer function into a something that can be implemented as a digital filter.

I already have working implementations for 1-5th order butterworth low pass filters, but I've noticed some issues with filter stability for filter orders above 3 when the cutoff frequency is a small with respect to the sample rate. Going to double precision helps a little bit but there are still filter instabilities.

I hear that it should be possible to implement higher order filters as cascades of low order filters to improve stability, but I could use some tips for how to do that. I know that cascading two n-order butterworth filters produces a 2n-order linkwitz-riley filter, but how can I determine the coefficients for the cascaded filters so that when combined they produce a butterworth response (or other canonical type of response)?

The math involved in this stuff is really new to me (and maybe a little out of my league) so please explain things in as simple of terms as you can!
On Mon, 26 Nov 2012 20:14:24 -0800, carl.schissler wrote:

> I'm trying to make a tool that implements various types of IIR filters > up to a certain maximum filter order (8 in my case). I've been learning > about the bilinear transform and how to use it to turn an analog > transfer function into a something that can be implemented as a digital > filter. > > I already have working implementations for 1-5th order butterworth low > pass filters, but I've noticed some issues with filter stability for > filter orders above 3 when the cutoff frequency is a small with respect > to the sample rate. Going to double precision helps a little bit but > there are still filter instabilities. > > I hear that it should be possible to implement higher order filters as > cascades of low order filters to improve stability, but I could use some > tips for how to do that. I know that cascading two n-order butterworth > filters produces a 2n-order linkwitz-riley filter, but how can I > determine the coefficients for the cascaded filters so that when > combined they produce a butterworth response (or other canonical type of > response)? > > The math involved in this stuff is really new to me (and maybe a little > out of my league) so please explain things in as simple of terms as you > can!
You have rediscovered a basic principle. Polynomial roots get more and more sensitive to coefficient values as the polynomial degree goes up; with filters, this translates to higher-order direct-form filters having pole positions that get more and more sensitive to precision as the filter order goes up. So you want to realize higher-order filters as cascades. This isn't just a "should be possible" thing -- this is a "just do all the time" thing. How to do? Assuming all lowpass: Find a reference that gives you the poles for an n-order Butterworth filter. For n even they'll all be complex pole pairs; for n odd they'll be one real pole and a bunch of complex pole pairs. For each pole pair, make up a 2nd-order filter. For the one real pole (if there), make up a 1st-order filter. Prewarp, and use the bilinear transform to get into the z domain. Realize your filters. -- My liberal friends think I'm a conservative kook. My conservative friends think I'm a liberal kook. Why am I not happy that they have found common ground? Tim Wescott, Communications, Control, Circuits & Software http://www.wescottdesign.com
On Tuesday, November 27, 2012 1:02:05 AM UTC-5, Tim Wescott wrote:
> On Mon, 26 Nov 2012 20:14:24 -0800, carl.schissler wrote: > > > > > I'm trying to make a tool that implements various types of IIR filters > > > up to a certain maximum filter order (8 in my case). I've been learning > > > about the bilinear transform and how to use it to turn an analog > > > transfer function into a something that can be implemented as a digital > > > filter. > > > > > > I already have working implementations for 1-5th order butterworth low > > > pass filters, but I've noticed some issues with filter stability for > > > filter orders above 3 when the cutoff frequency is a small with respect > > > to the sample rate. Going to double precision helps a little bit but > > > there are still filter instabilities. > > > > > > I hear that it should be possible to implement higher order filters as > > > cascades of low order filters to improve stability, but I could use some > > > tips for how to do that. I know that cascading two n-order butterworth > > > filters produces a 2n-order linkwitz-riley filter, but how can I > > > determine the coefficients for the cascaded filters so that when > > > combined they produce a butterworth response (or other canonical type of > > > response)? > > > > > > The math involved in this stuff is really new to me (and maybe a little > > > out of my league) so please explain things in as simple of terms as you > > > can! > > > > You have rediscovered a basic principle. Polynomial roots get more and > > more sensitive to coefficient values as the polynomial degree goes up; > > with filters, this translates to higher-order direct-form filters having > > pole positions that get more and more sensitive to precision as the > > filter order goes up. > > > > So you want to realize higher-order filters as cascades. This isn't just > > a "should be possible" thing -- this is a "just do all the time" thing. > > > > How to do? Assuming all lowpass: > > > > Find a reference that gives you the poles for an n-order Butterworth > > filter. For n even they'll all be complex pole pairs; for n odd they'll > > be one real pole and a bunch of complex pole pairs. > > > > For each pole pair, make up a 2nd-order filter. For the one real pole > > (if there), make up a 1st-order filter. > > > > Prewarp, and use the bilinear transform to get into the z domain. > > > > Realize your filters. > > > > -- > > My liberal friends think I'm a conservative kook. > > My conservative friends think I'm a liberal kook. > > Why am I not happy that they have found common ground? > > > > Tim Wescott, Communications, Control, Circuits & Software > > http://www.wescottdesign.com
Thanks, that helps I think. So, tell me if I am correct in my thinking here before I spend some time implementing it: For a 4th order butterworth filter, the factored approximate transfer function is: H(s) = 1 / ((s^2 + 0.7654*s + 1)(s^2 + 1.8478*s + 1)) So, is it true that I can separate the filter into its two 2nd order components: H1(s) = 1 / (s^2 + 0.7654*s + 1) H2(s) = 1 / (s^2 + 1.8478*s + 1) then apply the bilinear transform and determine their filter coefficients, then apply each 2nd filter in series to the signal to produce the same output as a single 4th order filter?
On Mon, 26 Nov 2012 22:51:24 -0800, Carl Schissler wrote:

> On Tuesday, November 27, 2012 1:02:05 AM UTC-5, Tim Wescott wrote: >> On Mon, 26 Nov 2012 20:14:24 -0800, carl.schissler wrote: >> >> >> >> > I'm trying to make a tool that implements various types of IIR >> > filters >> >> > up to a certain maximum filter order (8 in my case). I've been >> > learning >> >> > about the bilinear transform and how to use it to turn an analog >> >> > transfer function into a something that can be implemented as a >> > digital >> >> > filter. >> >> >> > >> > I already have working implementations for 1-5th order butterworth >> > low >> >> > pass filters, but I've noticed some issues with filter stability for >> >> > filter orders above 3 when the cutoff frequency is a small with >> > respect >> >> > to the sample rate. Going to double precision helps a little bit but >> >> > there are still filter instabilities. >> >> >> > >> > I hear that it should be possible to implement higher order filters >> > as >> >> > cascades of low order filters to improve stability, but I could use >> > some >> >> > tips for how to do that. I know that cascading two n-order >> > butterworth >> >> > filters produces a 2n-order linkwitz-riley filter, but how can I >> >> > determine the coefficients for the cascaded filters so that when >> >> > combined they produce a butterworth response (or other canonical type >> > of >> >> > response)? >> >> >> > >> > The math involved in this stuff is really new to me (and maybe a >> > little >> >> > out of my league) so please explain things in as simple of terms as >> > you >> >> > can! >> >> >> >> You have rediscovered a basic principle. Polynomial roots get more and >> >> more sensitive to coefficient values as the polynomial degree goes up; >> >> with filters, this translates to higher-order direct-form filters >> having >> >> pole positions that get more and more sensitive to precision as the >> >> filter order goes up. >> >> >> >> So you want to realize higher-order filters as cascades. This isn't >> just >> >> a "should be possible" thing -- this is a "just do all the time" thing. >> >> >> >> How to do? Assuming all lowpass: >> >> >> >> Find a reference that gives you the poles for an n-order Butterworth >> >> filter. For n even they'll all be complex pole pairs; for n odd >> they'll >> >> be one real pole and a bunch of complex pole pairs. >> >> >> >> For each pole pair, make up a 2nd-order filter. For the one real pole >> >> (if there), make up a 1st-order filter. >> >> >> >> Prewarp, and use the bilinear transform to get into the z domain. >> >> >> >> Realize your filters. >> >> >> >> -- >> >> My liberal friends think I'm a conservative kook. >> >> My conservative friends think I'm a liberal kook. >> >> Why am I not happy that they have found common ground? >> >> >> >> Tim Wescott, Communications, Control, Circuits & Software >> >> http://www.wescottdesign.com > > Thanks, that helps I think. So, tell me if I am correct in my thinking > here before I spend some time implementing it: > > For a 4th order butterworth filter, the factored approximate transfer > function is: > > H(s) = 1 / ((s^2 + 0.7654*s + 1)(s^2 + 1.8478*s + 1)) > > So, is it true that I can separate the filter into its two 2nd order > components: > > H1(s) = 1 / (s^2 + 0.7654*s + 1) > H2(s) = 1 / (s^2 + 1.8478*s + 1) > > then apply the bilinear transform and determine their filter > coefficients, then apply each 2nd filter in series to the signal to > produce the same output as a single 4th order filter?
You got it! The only caveat is that you'll run afoul of the same polynomial coefficient sensitivity issue in your math package as you do when implementing your filters. So you do _not_ want to use a method that numerically generates a great big polynomial, then factor it: you want to use a method that algorithmically generates the poles, then select pairs to make your polynomials. These algorithms exist for both Butterworth and Chebychev polynomials. The Butterworth is easy: they're just evenly distributed on a half circle centered at s = 0. I can't remember about the Chebychev, so you'll have to look for yourself. -- Tim Wescott Control system and signal processing consulting www.wescottdesign.com
On Tuesday, November 27, 2012 11:10:42 AM UTC-5, Tim Wescott wrote:
> On Mon, 26 Nov 2012 22:51:24 -0800, Carl Schissler wrote: > > > > > On Tuesday, November 27, 2012 1:02:05 AM UTC-5, Tim Wescott wrote: > > >> On Mon, 26 Nov 2012 20:14:24 -0800, carl.schissler wrote: > > >> > > >> > > >> > > >> > I'm trying to make a tool that implements various types of IIR > > >> > filters > > >> > > >> > up to a certain maximum filter order (8 in my case). I've been > > >> > learning > > >> > > >> > about the bilinear transform and how to use it to turn an analog > > >> > > >> > transfer function into a something that can be implemented as a > > >> > digital > > >> > > >> > filter. > > >> > > >> > > >> > > > >> > I already have working implementations for 1-5th order butterworth > > >> > low > > >> > > >> > pass filters, but I've noticed some issues with filter stability for > > >> > > >> > filter orders above 3 when the cutoff frequency is a small with > > >> > respect > > >> > > >> > to the sample rate. Going to double precision helps a little bit but > > >> > > >> > there are still filter instabilities. > > >> > > >> > > >> > > > >> > I hear that it should be possible to implement higher order filters > > >> > as > > >> > > >> > cascades of low order filters to improve stability, but I could use > > >> > some > > >> > > >> > tips for how to do that. I know that cascading two n-order > > >> > butterworth > > >> > > >> > filters produces a 2n-order linkwitz-riley filter, but how can I > > >> > > >> > determine the coefficients for the cascaded filters so that when > > >> > > >> > combined they produce a butterworth response (or other canonical type > > >> > of > > >> > > >> > response)? > > >> > > >> > > >> > > > >> > The math involved in this stuff is really new to me (and maybe a > > >> > little > > >> > > >> > out of my league) so please explain things in as simple of terms as > > >> > you > > >> > > >> > can! > > >> > > >> > > >> > > >> You have rediscovered a basic principle. Polynomial roots get more and > > >> > > >> more sensitive to coefficient values as the polynomial degree goes up; > > >> > > >> with filters, this translates to higher-order direct-form filters > > >> having > > >> > > >> pole positions that get more and more sensitive to precision as the > > >> > > >> filter order goes up. > > >> > > >> > > >> > > >> So you want to realize higher-order filters as cascades. This isn't > > >> just > > >> > > >> a "should be possible" thing -- this is a "just do all the time" thing. > > >> > > >> > > >> > > >> How to do? Assuming all lowpass: > > >> > > >> > > >> > > >> Find a reference that gives you the poles for an n-order Butterworth > > >> > > >> filter. For n even they'll all be complex pole pairs; for n odd > > >> they'll > > >> > > >> be one real pole and a bunch of complex pole pairs. > > >> > > >> > > >> > > >> For each pole pair, make up a 2nd-order filter. For the one real pole > > >> > > >> (if there), make up a 1st-order filter. > > >> > > >> > > >> > > >> Prewarp, and use the bilinear transform to get into the z domain. > > >> > > >> > > >> > > >> Realize your filters. > > >> > > >> > > >> > > >> -- > > >> > > >> My liberal friends think I'm a conservative kook. > > >> > > >> My conservative friends think I'm a liberal kook. > > >> > > >> Why am I not happy that they have found common ground? > > >> > > >> > > >> > > >> Tim Wescott, Communications, Control, Circuits & Software > > >> > > >> http://www.wescottdesign.com > > > > > > Thanks, that helps I think. So, tell me if I am correct in my thinking > > > here before I spend some time implementing it: > > > > > > For a 4th order butterworth filter, the factored approximate transfer > > > function is: > > > > > > H(s) = 1 / ((s^2 + 0.7654*s + 1)(s^2 + 1.8478*s + 1)) > > > > > > So, is it true that I can separate the filter into its two 2nd order > > > components: > > > > > > H1(s) = 1 / (s^2 + 0.7654*s + 1) > > > H2(s) = 1 / (s^2 + 1.8478*s + 1) > > > > > > then apply the bilinear transform and determine their filter > > > coefficients, then apply each 2nd filter in series to the signal to > > > produce the same output as a single 4th order filter? > > > > You got it! > > > > The only caveat is that you'll run afoul of the same polynomial > > coefficient sensitivity issue in your math package as you do when > > implementing your filters. So you do _not_ want to use a method that > > numerically generates a great big polynomial, then factor it: you want to > > use a method that algorithmically generates the poles, then select pairs > > to make your polynomials. > > > > These algorithms exist for both Butterworth and Chebychev polynomials. > > The Butterworth is easy: they're just evenly distributed on a half circle > > centered at s = 0. I can't remember about the Chebychev, so you'll have > > to look for yourself. > > > > -- > > Tim Wescott > > Control system and signal processing consulting > > www.wescottdesign.com
Thanks a lot! I've now got butterworth hi and lowpass filters working for arbitrarily high orders! I'm just using the transfer function definition from wikipedia to get the factored polynomials directly, generating the polynomial coefficients at runtime and then realizing the filters from there. It looks like this method will work well for other filter types too.