How to design digital filter from transfer function

Question

I have two transfer functions and would like to implement an IIR with coefficients.

What is the process to get from my complex transfer function to IIR coefficients?

\begin{align} F_W &=\frac{1+j\displaystyle \frac{\omega r}{c}-\frac 13\left(\frac{\omega r}{c}\right)^2}{1+\displaystyle \frac 13 j\frac{\omega r}{c}}\\ F_{XYZ} &=\sqrt{6}\frac{1+\displaystyle\frac 13j \frac{\omega r}{c}-\frac 13\left(\frac{\omega r}{c}\right)^2}{1+\displaystyle \frac 13 j\frac{\omega r}{c}}\\ \end{align}

These are the Gerzon/Craven transfer functions for Ambisonic B-Format audio

Answer 1

Following OP comment requesting more details to Maximilian Matthé's answer:

Just like mentioned, the key to your question is the Bilinear Transform, where you substitute $s$ (or $j\omega$ in your case) with $\frac{2}{T}\frac{z-1}{z+1}$

This transform will have for effect to convert your transfer function from the continuous $s$ domain to the discrete $z$ domain.

I assume that you are fimiliar with the $z$ transform and that you know that a delay in the time domain is equivalent to multiplying by $z^{-1}$ in the z domain. $Z\{f[x-n]\} = Z\{f[x]\}z^{-n}$

The implementation of a discrete filter (FIR or IIR) will happen in the time domain. Your software will implement something like

$$y[n] = b_{0}x[n]+b_{1}x[n-1]+b_{2}x[n-2]...+b_{N}x[n-N]-a_{1}y[n-1] - a_{2}y[n-2]...-a_{N}y[n-N]$$

In my previous equations :

• $x$ is the input
• $y$ is the output
• $n$ is the index of the actual value in your buffer. $n-1$ being the previous value and so on.

If you bring that equation in the $z$ domain, you will have an equivalent like that.

$$Y = b_{0}X + b_{1}Xz^{-1}+b_{2}Xz^{-2}...+b_{N}Xz^{-N} -a_{1}Yz^{-1}-a_{2}Yz^{-2}...-a_{N}Yz^{-N}$$

After some algebra, you'll get the following form, which correspond to a transfer function (ratio of output on input)

$$\frac{Y}{X} = \frac{b_{0}+b_{1}z^{-1}+b_{2}z^{-2}...+b_{N}z^{-N}}{1+a_{1}z^{-1}+a_{2}z^{-2}...+a_{N}z^{-N}}$$

This transfer function is what Maximilian Matthé answer refers to. In it's condensed form :

$$\frac{Y}{X}=\frac{\sum_{i=0}^{N}b_{i}z^{-i}}{1+\sum_{i=1}^{N}a_{i}z^{-i}}$$

Designing a discrete filter correspond to determining the coefficients of that transfer function. e.g. Finding the values of $a$ and $b$.

The value of your coefficients will change according to your sampling frequency, and you can confirm that considering that the Bilinear Transform introduce the sampling time $T$ in your equation.

You might want to consider the fact that the Bilinear transform is an approximation. If you look at Wikipedia's page, you will find this :

To compensate this approximation, you can substitute $\omega$ with

$$\omega_d = \frac{2}{T}\tan^{-1}(\frac{2\omega}{T})$$

This substitution is called frequency warping and will have a slight effect on the frequency response of your filter, but will most likely be negligible.

Finally, if you have access to matlab and the signal processing toolbox, you might want to have a look to bilinear() which will do all of this in a single command.

Answer 2

A common way is the application of the Bilineartransform.

In your equations you first substitute $j\omega=s$ and then

$$ s=\frac{2}{T}\frac{z-1}{z+1}, $$

where $1/T$ is your sampling frequency. Then you have a digital filter which you can directly read the coefficients from to implement it.

So, for the first equation, let's reformulate it as

$$ F_W(j\omega) = \frac{1+j\omega\tau+\frac{1}{3}(j\omega)^2\tau^2}{1+\frac{1}{3}j\omega\tau} $$

where I have set $r/c=\tau$. Now, we substitute $j\omega=s$ to get

$$ F_W(s)=\frac{1+s\tau+\frac{1}{3}s^2\tau^2}{1+\frac{1}{3}s\tau} $$

You now have to substitute $s=\frac{2}{T}\frac{z-1}{z+1}$ to get an expression like

$$ F_W(z) = \frac{\sum_{i=0}^N b_i z^{-i}}{\sum_{i=0}^{N} a_i z^{-i}} $$

(I'm currently to lazy to calculate $a_i$ and $b_i$ explicitely, that's your task or ask Mathematical et.al.) from which you can directly read the filter coefficients for the implementation. Have a look here for understanding what $F_W(z)$ means and how to implement it (or you use a ready-made library, such as Matlab's filter or scipy.signal lfilter, where you can directly put b and a as parameters. Also, have a look at the example in the wikipedia page on Bilineartransform.

Note that the resulting filter will be IIR (because your analog filter is also IIR). You can approximate this with a FIR filter (if necessary) by calculating the impulse response and then truncate the response at some point. The resulting impulse response directly gives you the FIR filter coefficients.