I need to convert a fourth-order IIR filter made using two cascaded second-order sections with the transfer function: $$H(z) = \frac{b2_{1}z^{-2} + b1_{1}z^{-1} + b0_{1}}{a2_{1}z^{-2} + a1_{1}z^{-1} + a0_{1}}\cdot \frac{b2_{2}z^{-2} + b1_{2}z^{-1} + b0_{2}}{a2_{2}z^{-2} + a1_{2}z^{-1} + a0_{2}}$$ to a parallel form like this: $$H(z) = \frac{b2_{3}z^{-2} + b1_{3}z^{-1} + b0_{3}}{a2_{3}z^{-2} + a1_{3}z^{-1} + a0_{3}}+\frac{b2_{4}z^{-2} + b1_{4}z^{-1} + b0_{4}}{a2_{4}z^{-2} + a1_{4}z^{-1} + a0_{4}}$$ while preseving $$H(z)$$.

I thought the key was finding the coefficients $$b0_{3}, b1_{3}, b2_{3}, b0_{4}, b1_{4}, b2_{4}, a0_{3}, a1_{3}, a2_{3}, a0_{4}, a1_{4}$$ and $$a2_{4}$$ so I tried solving for them using the following equations:

\begin{align} a0_{1}a0_{2} &= a0_{1}a0_{2}\\ a2_{1}a2_{2} &= a2_{1}a2_{2}\\ a0_{1}a1_{2} + a0_{2}a1_{1} &= a0_{1}a1_{2} + a0_{2}a1_{1}\\ a0_{1}b0_{4} + a0_{2}b0_{3} &= b0_{1}b0_{2}\\ a1_{1}a2_{2} + a1_{2}a2_{1} &= a1_{1}a2_{2} + a1_{2}a2_{1}\\ a2_{1}b2_{4} + a2_{2}b2_{3} &= b2_{1}b2_{2}\\ a0_{1}a2_{2} + a0_{2}a2_{1} + a1_{1}a1_{2} &= a0_{1}a2_{2} + a0_{2}a2_{1} + a1_{1}a1_{2}\\ a0_{1}b1_{4} + a0_{2}b1_{3} + a1_{1}b0_{4} + a1_{2}b0_{3} &= b0_{1}b1_{2} + b0_{2}b1_{1}\\ a1_{1}b2_{4} + a1_{2}b2_{3} + a2_{1}b1_{4} + a2_{2}b1_{3} &= b1_{1}b2_{2} + b1_{2}b2_{1}\\ a0_{1}b2_{4} + a0_{2}b2_{3} + a1_{1}b1_{4} + a1_{2}b1_{3} + a2_{1}b0_{4} + a2_{2}b0_{3} &= b0_{1}b2_{2} + b0_{2}b2_{1} + b1_{1}b1_{2}\\ \end{align} to no avail...

Is there another method to get the parallel biquads' coefficients?

EDIT

The problem is fully solved now as suggested by Hilmar.

The solution will be posted below.

• mathworks.com/help/signal/ref/sos2tf.html – Juha P Jun 28 at 7:36
• @JuhaP That doesn't return the coefficients of the individual biquads (second order sections) for summation as needed. – Nizar Nizar Jun 28 at 8:32
• Yes, I know, it returns the cascaded 4th order filter. Anyway, what type of filter you have there? You sure know that when you cascade to second order sections you use different Q value for each section. Q value calculator: earlevel.com/main/2016/09/29/cascading-filters – Juha P Jun 28 at 11:18
• I have no problem getting the coefficients of each cascaded section or the resulting 4th-order filter. I need to solve this problem symbolically because the original cascaded filters must be calculated on the fly in the target application so they can be any kind of filters, then I have to create a routine to convert the cascaded system to a parallel (summed) one as shown in the second H(z). – Nizar Nizar Jun 28 at 11:40

You need to do partial fraction expansion.

Most of the time, the result will be 3 parallel sections with the same poles as the original filter but only with a real zero in each parallel section. The first section is just a single real scale factor.

$$H(z) = r0_0 + \frac{r0_1 + r1_1 \cdot z^{-1}}{a0_1 + a1_1 \cdot z^{-1} + a2_1 \cdot z^{-2}} + \frac{r0_2 + r1_2 \cdot z^{-1}}{a0_2 + a1_2 \cdot z^{-1} + a2_2 \cdot z^{-2}}$$

Note that the "a" coefficients are the same as in the cascaded form, i.e. you only need to calculate the "r" coefficients.

This is all relatively straight forward if your poles are unique and complex conjugated pairs. Repeated poles and single real poles require special treatment.

The Matlab function residuez() can typically help with this although it can have numerical problems if the poles are too close to the unit circle. If you need to do it by hand take a look at https://ccrma.stanford.edu/~jos/filters/Partial_Fraction_Expansion.html

Since you only have two sections, the following approach should work

1. Write both versions (parallel and serial) as a single fraction
2. Discard the denominators, they are the same
3. Both numerators are a 4th order polynomial in $$z^{-1}$$, i.e. they have 5 coefficients.
4. Each power of $$z^{-1}$$ needs to have the same coefficient in both numerators.
5. That gives you five linear equations with 5 unknowns that you can solve for $$r0_0, r0_1, r1_1, r0_2, r1_2$$
• Repeated poles made $cd=0$. What should be done in that case? BTW I've tried the same with 4 second order sections and it worked well but with the same issue. In the case of more than 2 biquad sections, the term $cd$ becomes the poduct of each $cd$ where the parallel biquad in question is involved. I've made a function to generate it if someone is interested.. – Nizar Nizar Jul 3 at 22:26

There is the solution as announced in the question text: The transfer function $$H(z) = \frac{b2_{1}z^{-2} + b1_{1}z^{-1} + b0_{1}}{a2_{1}z^{-2} + a1_{1}z^{-1} + a0_{1}}.\frac{b2_{2}z^{-2} + b1_{2}z^{-1} + b0_{2}}{a2_{2}z^{-2} + a1_{2}z^{-1} + a0_{2}}$$ was converted to the following parallel form composed of the biquads $$bq1$$ and $$bq2$$ $$H(z) = r0 + \frac{1}{cd}\frac{b1_{bq1}z^{-1} + b0_{bq1}}{a2_{bq1}z^{-2} + a1_{bq1}z^{-1} + a0_{bq1}} + \frac{1}{cd}\frac{b1_{bq2}z^{-1} + b0_{bq2}}{a2_{bq2}z^{-2} + a1_{bq2}z^{-1} + a0_{bq2}}$$ with $$r0 = \frac{b2_{1}b2_{2}}{a2_{1}a2_{2}},\\ cd = (a0_{1}a2_{2}-a0_{2}a2_{1})^2+(-a1_{1}a2_{2}+a1_{2}a2_{1})(a0_{1}a1_{2}-a0_{2}a1_{1}),\\ b0_{bq1} = -a0_{1}^3a2_{2}b2_{1}b2_{2} + ((b2_{1}b2_{2}a0_{2} + (b0_{1}b2_{2} + b0_{2}b2_{1} + b1_{1}b1_{2})a2_{2} - a1_{2}(b1_{1}b2_{2} + b1_{2}b2_{1}))a2_{1} + b2_{1}b2_{2}a1_{1}a1_{2})a0_{1}^2 + (((-b0_{1}b2_{2} - b0_{2}b2_{1} - b1_{1}b1_{2})a0_{2} - b0_{1}b0_{2}a2_{2} + a1_{2}(b0_{1}b1_{2} + b0_{2}b1_{1}))a2_{1}^2 + ((b1_{1}b2_{2} + b1_{2}b2_{1})a0_{2} - a2_{2}(b0_{1}b1_{2} + b0_{2}b1_{1}))a1_{1}a2_{1} - b2_{1}b2_{2}a0_{2}a1_{1}^2)a0_{1} + a2_{1}b0_{1}b0_{2}(a0_{2}a2_{1}^2 + a1_{1}^2a2_{2} - a1_{1}a1_{2}a2_{1}),\\ b1_{bq1} = ((b0_{1}b1_{2} + b0_{2}b1_{1})a0_{2} - a1_{2}b0_{1}b0_{2})a2_{1}^3 + (((-b1_{1}b2_{2} - b1_{2}b2_{1})a0_{2} + a1_{2}b0_{2}b2_{1} + a1_{2}b0_{1}b2_{2} + b1_{2}a1_{2}b1_{1} - a2_{2}(b0_{1}b1_{2} + b0_{2}b1_{1}))a0_{1} - ((b0_{1}b2_{2} + b0_{2}b2_{1} + b1_{1}b1_{2})a0_{2} - b0_{1}b0_{2}a2_{2})a1_{1})a2_{1}^2 + (((-a1_{2}b2_{2} + a2_{2}b1_{2})b2_{1} + b2_{2}a2_{2}b1_{1})a0_{1}^2 - (-2b2_{1}b2_{2}a0_{2} + a1_{2}(b1_{1}b2_{2} + b1_{2}b2_{1}))a1_{1}a0_{1} + a0_{2}a1_{1}^2(b1_{1}b2_{2} + b1_{2}b2_{1}))a2_{1} - b2_{1}b2_{2}a1_{1}(a0_{1}^2a2_{2} - a0_{1}a1_{1}a1_{2} + a0_{2}a1_{1}^2),\\ a0_{bq1} = a0_{1}a2_{1},\\ a1_{bq1} = a1_{1}a2_{1},\\ a2_{bq1} = a2_{1}^2,\\ b0_{bq2} = -a0_{2}^3a2_{1}b2_{1}b2_{2} + ((a0_{1}b2_{1}b2_{2} + (b0_{1}b2_{2} + b0_{2}b2_{1} + b1_{1}b1_{2})a2_{1} - a1_{1}(b1_{1}b2_{2} + b1_{2}b2_{1}))a2_{2} + b2_{1}b2_{2}a1_{1}a1_{2})a0_{2}^2 + (((-b0_{1}b2_{2} - b0_{2}b2_{1} - b1_{1}b1_{2})a0_{1} - a2_{1}b0_{1}b0_{2} + a1_{1}(b0_{1}b1_{2} + b0_{2}b1_{1}))a2_{2}^2 + a1_{2}((b1_{1}b2_{2} + b1_{2}b2_{1})a0_{1} - a2_{1}(b0_{1}b1_{2} + b0_{2}b1_{1}))a2_{2} - a0_{1}b2_{1}b2_{2}a1_{2}^2)a0_{2} + b0_{1}b0_{2}a2_{2}(a0_{1}a2_{2}^2 - a1_{1}a1_{2}a2_{2} + a1_{2}^2a2_{1}),\\ b1_{bq2} = ((b0_{1}b1_{2} + b0_{2}b1_{1})a0_{1} - a1_{1}b0_{1}b0_{2})a2_{2}^3 + (((-b1_{1}b2_{2} - b1_{2}b2_{1})a0_{1} + a1_{1}b0_{2}b2_{1} + a1_{1}b0_{1}b2_{2} + b1_{2}a1_{1}b1_{1} - a2_{1}(b0_{1}b1_{2} + b0_{2}b1_{1}))a0_{2} - ((b0_{1}b2_{2} + b0_{2}b2_{1} + b1_{1}b1_{2})a0_{1} - a2_{1}b0_{1}b0_{2})a1_{2})a2_{2}^2 + (((-a1_{1}b2_{2} + a2_{1}b1_{2})b2_{1} + b2_{2}a2_{1}b1_{1})a0_{2}^2 + 2a1_{2}(a0_{1}b2_{1}b2_{2} - a1_{1}(b1_{1}b2_{2} + b1_{2}b2_{1})/2)a0_{2} + a0_{1}a1_{2}^2(b1_{1}b2_{2} + b1_{2}b2_{1}))a2_{2} - b2_{1}b2_{2}a1_{2}(a0_{1}a1_{2}^2 + a0_{2}^2a2_{1} - a0_{2}a1_{1}a1_{2}),\\ a0_{bq2} = a0_{2}a2_{2},\\ a1_{bq2} = a1_{2}a2_{2}\text{ and}\\ a2_{bq2} = a2_{2}^2\\$$ as schematized below. You can use the original denominator coefficients with a small modification: $$H(z)=r0+ \frac{1}{a2_{1}cd} \frac{b1_{bq1}z^{-1}+b0_{bq1}}{z^{-2}a2_{1}+z^{-1}a1_{1}+a0_{1}}+ \frac{1}{a2_{2}cd} \frac{b1_{bq2}z^{-1}+b0_{bq2}}{z^{-2}a2_{2}+z^{-1}a1_{2}+a0_{2}}$$ I think there is room for simplification concerning numerator coefficients.

Both transfer functions have been tested for equality symbolically and here is a sample plot made using two cascaded peaking filters showing correct transformation:  The test was made using the following original biquad coefficients: \begin{align} b0_{1}&=1.084035017747234\\ b1_{1}&=-0.2664645875629038\\ b2_{1}&=0.7470948934664502\\ a0_{1}&=1\\ a1_{1}&=-0.2664645875629038\\ a2_{1}&=0.8311299112136844\\ b0_{2}&=1.031331836346721\\ b1_{2}&=-1.8154546438670163\\ b2_{2}&=0.9057061969928599\\ a0_{2}&=1\\ a1_{2}&=-1.8154546438670163\\ a2_{2}&=0.9370380333395809,\\ \end{align} giving the following coefficients, after transformation for parallel application: \begin{align} r0&=0.8688343370\\ cd&=1.961691184\\ b0_{bq1}&=0.3012806494\\ b1_{bq1}&=-0.0328276431\\ a0_{bq1}&=0.8311299114\\ a1_{bq1}&=-0.2214666890\\ a2_{bq1}&=0.6907769291\\ b0_{bq2}&=0.118338990\\ b1_{bq2}&=-0.1125429904\\ a0_{bq2}&=0.9370380333\\ a1_{bq2}&=-1.701150049\\ a2_{bq2}&=0.8780402757\\ \end{align}