Page 1 of 2
How to split high-order IIR to multiple biquads?
Posted: Thu Jan 16, 2014 6:36 pm
by KG_is_back
the tittle says it all... let's say I have coefficients of a N-th order IIR filter and I what to split it into cascade of biquads. How do I compute the coeffiecients? Is there any particular algorithm that can be implemented automatically? (like: I enter IIR coefficients and it outputs an array of multiple biquad coefs)
Re: How to split high-order IIR to multiple biquads?
Posted: Fri Jan 17, 2014 4:08 pm
by martinvicanek
KG, you need to locate the poles and zeros of the original nth order filter, then collect the complex conjugate pairs to form biquads. This amounts to a factorization of polynomials. There are ready-to-use routines for that, some are even accessible online. Dont't know about DLLs which you could integrate with FS, though.
Re: How to split high-order IIR to multiple biquads?
Posted: Fri Jan 17, 2014 11:32 pm
by KG_is_back
To combine two filters one just have to multiply their Transfer functions. That means (b_0 + b_1*z^1+b_2*z^2)*(b_0 + b_1*z^1+b_2*z^2) and the same with A coefficients. Maybe the process can be reversed by dividing the transfer function with biquad transfer function so it eliminates two last coefficients... It would be hard to find such function though (and even harder to create program that finds it). Is there any algorithm that might help? maybe iterative one...
Re: How to split high-order IIR to multiple biquads?
Posted: Sat Jan 18, 2014 12:37 pm
by martinvicanek
KG_is_back wrote:Maybe the process can be reversed by dividing the transfer function with biquad transfer function so it eliminates two last coefficients... It would be hard to find such function though (and even harder to create program that finds it). Is there any algorithm that might help? maybe iterative one...
Yes, that is basically what is known as deflation:
Code: Select all
1. Start with a nth order polynomial P(x) = a0 + a1*x + a2*x^2 + ... + an*x^n
2. Find one root x0 [Note: x0 is a root of P if P(x0) = 0. There are n roots in total.]
3. Divide the original polynomial by (x - x0) to get a new polynomial of degree n-1
4. Repeat steps 2.-3. until all roots are found
Find more information and programs here:
http://en.wikipedia.org/wiki/Root-findi ... olynomialshttp://www.akiti.ca/PolyRootRe.htmlhttp://www.hvks.com/index.htmlhttp://apps.nrbook.com/c/index.html go to page 369
Re: How to split high-order IIR to multiple biquads?
Posted: Sat Jan 18, 2014 2:55 pm
by KG_is_back
This is cool! I can actually do this!
To find the root is simple - Newton method. It would be best if it could be done in ruby as it uses DP-math. I'll probably have to start learning ruby after all.
The division of the polynomials can probably also be done very simply as you are alway dividing by (1*x-x0) where x0 is a constant.
The formula (1*x-x0) can be also written as (b0*z^0 + b1*z^-1) where b0=1 and b1=-x0 so it is nominator of the transfer function of first order filter. The same algorithm would be used for A coefficients.
This way N-th order filter will be split into N first-order filters and they can easily be summed to N/2 biquads.
Re: How to split high-order IIR to multiple biquads?
Posted: Mon Aug 18, 2014 6:12 pm
by KG_is_back
Finally after such a long time I have it figured out. I've created module, that receives coefficients of higher order IIR and creates 5 arrays of biquad coefficients (how much is needed). After casual testing it shows pretty good - the difference is under 60dB across the whole spectrum so it's basically identical.
Here is an example schematic: it has 3 parallel chains:
First one has two biquads of which response you control.
second one is a 4th order filter combined from those original biquads.
third one has two biquads that are reconstructed from the 4th order filter coefficients. Note that these 2 biquads may have different set of coefficients then the original one, because there are multiple ways you can combine 4poles and zeroes into two 2pole-2zeroes filter (biquad). However the filter responses are (almost) identical, so it works fine.
Re: How to split high-order IIR to multiple biquads?
Posted: Mon Aug 18, 2014 6:38 pm
by MyCo
Interesting, but the schematic is probably the worst I've seen in a while
To calculate the difference, I think you have to divide not subtract. You are operating in linear scale (not log)
Re: How to split high-order IIR to multiple biquads?
Posted: Mon Aug 18, 2014 7:29 pm
by KG_is_back
MyCo wrote:Interesting, but the schematic is probably the worst I've seen in a while
To calculate the difference, I think you have to divide not subtract. You are operating in linear scale (not log)
I've set it up in like 5 minutes, just to find out if it's actually working (and it does, so hureyyyyy
).
I'm pretty sure that difference=subtraction not division. Simply because some of the values might be zero (and probably are, near the end of the IR) and dividing by zero gives NaN. I'm working out the difference in time domain, not frequency domain. Still when you subtract the original and reconstructed you get something that peaks around 10^-6 which is basically a rounding error.
Re: How to split high-order IIR to multiple biquads?
Posted: Mon Aug 18, 2014 7:44 pm
by martinvicanek
Respect, that's quite some serious stuff there, KG!
Although I haven't tested it yet (too much spaghetti for me
) Do you think it could be done with a 14th order allpole IIR at hop(1024)? I am thinking of an LPC application...
Re: How to split high-order IIR to multiple biquads?
Posted: Mon Aug 18, 2014 8:00 pm
by KG_is_back
martinvicanek wrote:Respect, that's quite some serious stuff there, KG! Although I haven't tested it yet (too much spaghetti for me ) Do you think it could be done with a 14th order allpole IIR at hop(1024)? I am thinking of an LPC application...
Yes, it should work for 14th order allpole IIR. However, I do not know what you mean with that hop(1024). The thing is done in ruby via event method. But it could possibly be moded to work in sync with frames.
By the way the A and B coefficients are calculated completely separately, so you may look into the module and delete a bunch of stuff so you can input A and B separately. Just remember for the individual A and B sections you have to input ALL coefficients - biquads are usually normalized to a0=1, but here you have to input that too - it is not done automatically. In case your schematic is not normalized to a0 do so or replace the A calculator with B calculator and then normalize the output biquads.