PDF letter size

## IIR digital ﬁlter design for low pass ﬁlter based on impulse invariance and bilinear transformation methods using butterworth analog ﬁlter

May 5, 2010   Compiled on January 30, 2024 at 6:45am

### Contents

A low pass digital ﬁlter (IIR) is designed based on given speciﬁcations. Butterworth analog ﬁlter $$H(s)$$ is designed ﬁrst, then it is converted to digital ﬁlter $$H(z)$$

Analytical procedure is illustrated below and simpliﬁed to allow one to more easily program the algorithm. Four diﬀerent numerical examples are used to illustrate the procedure.

### 1 Introduction

This report derives a symbolic procedure to design a low pass IIR digital ﬁlter from an analog Butterworth ﬁlter using 2 methods: impulse invariance and bilinear transformation. Two numerical examples are used to illustrate using the symbolic procedure.

There are a total of 13 steps. A Mathematica program with GUI is written to enable one to use this design for diﬀerent parameters. A Matlab script is written as well.

#### 1.1 Filter speciﬁcations

Filter speciﬁcations are 5 parameters. The frequency speciﬁcations use analog frequencies, while the attenuation for the passband and the stopband are given in db units.

 $$F_{s}$$ The sampling frequency in Hz $$f_{c}$$ The passband cutoﬀ frequency in Hz $$f_{s}$$ The stopband corner frequency in Hz $$\delta _{p}$$ The attenuation in db at $$f_{c}$$ $$\delta _{s}$$ The attenuation in db at $$f_{s}$$

This diagram below illustrates these speciﬁcations

The frequency speciﬁcations are in hz, and they are ﬁrst converted to digital frequencies $$\omega$$ where $$0\leq \left \vert \omega \right \vert \leq \pi$$ before using the attenuation speciﬁcations. The sampling frequency $$F_{s}$$ is used to do this conversion where $$F_{s}$$ corresponds to $$2\pi$$ on the digital frequency scale.

### 2 The design steps

1. Speciﬁcations are converted from analog to digital frequencies.
2. Based on design method (impulse invariance of bilinear) the attenuation criteria is applied to determine $$\Omega _{c}$$ and $$N$$ (the ﬁlter order).
3. Using $$\Omega _{c}$$ and $$N$$ the locations of the poles of the butterworth analog ﬁlter $$H(s)$$ are found.
4. $$H(z)$$ is determined from $$H(s)$$. The method of doing this depends if impulse invariance or bilinear is used. This step is much simpler for the bilinear method as it does not require performing partial fractions decomposition on $$H(s).$$

The analytical design procedure is given below.

#### 2.1 Impulse invariance method

Analog speciﬁcations is converted to digital speciﬁcations: $$\frac {F_{s}}{2\pi }=\frac {f_{p}}{\omega _{p}}$$, hence $$\omega _{p}=2\pi \frac {f_{p}}{F_{s}}$$ and $$\omega _{s}=2\pi \frac {f_{s}}{F_{s}}$$. Next, the criteria relative to the digital normalized scale is found\begin {align*} 20\log \left \vert H\left ( e^{j\omega _{p}}\right ) \right \vert & \geq \delta _{p}\\ 20\log \left \vert H\left ( e^{j\omega _{s}}\right ) \right \vert & \leq \delta _{s} \end {align*}

Therefore\begin {align} \left \vert H\left ( e^{j\omega _{p}}\right ) \right \vert & \geq 10^{\frac {\delta _{p}}{20}}T\tag {A}\\ \left \vert H\left ( e^{j\omega _{s}}\right ) \right \vert & \leq 10^{\frac {\delta _{s}}{20}}T\tag {B} \end {align}

Notice that $$T$$ scaling factor was added above. The Butterworth analog ﬁlter squared magnitude Fourier transform is given by$\left \vert H_{a}\left ( j\Omega \right ) \right \vert ^{2}=\frac {T^{2}}{1+\left ( \frac {j\Omega }{j\Omega _{c}}\right ) ^{2N}}$ Equations (A) and (B) above are now written in terms of the analog Butterworth amplitude frequency response giving\begin {align*} \frac {T^{2}}{1+\left ( \frac {\Omega _{p}}{\Omega _{c}}\right ) ^{2N}} & \geq \left ( 10^{\frac {\delta _{p}}{20}}\right ) ^{2}T^{2}\\ \frac {T^{2}}{1+\left ( \frac {\Omega _{s}}{\Omega _{c}}\right ) ^{2N}} & \leq \left ( 10^{\frac {\delta _{s}}{20}}\right ) ^{2}T^{2} \end {align*}

Hence\begin {align*} \frac {1}{1+\left ( \frac {\Omega _{p}}{\Omega _{c}}\right ) ^{2N}} & \geq \ 10^{\frac {\delta _{p}}{10}}\\ \frac {1}{1+\left ( \frac {\Omega _{s}}{\Omega _{c}}\right ) ^{2N}} & \leq 10^{\frac {\delta _{s}}{10}} \end {align*}

For impulse invariance, $$\Omega _{p}=\frac {\omega _{p}}{T}\,\$$ and $$\Omega _{s}=\frac {\omega _{s}}{T}$$.  Therefore\begin {align} 1+\left ( \frac {\omega _{p}/T}{\Omega _{c}}\right ) ^{2N} & \leq 10^{-\frac {\delta _{p}}{10}}\tag {1}\\ 1+\left ( \frac {\omega _{s}/T}{\Omega _{c}}\right ) ^{2N} & \geq 10^{-\frac {\delta _{s}}{10}}\tag {2} \end {align}

Changing inequalities to equalities and simplifying gives\begin {align*} \left ( \frac {\omega _{p}/T}{\Omega _{c}}\right ) ^{2N} & =10^{-\frac {\delta _{p}}{10}}-1\\ \left ( \frac {\omega _{s}/T}{\Omega _{c}}\right ) ^{2N} & =10^{-\frac {\delta _{s}}{10}}-1 \end {align*}

Dividing the above two equations results in\begin {align*} \left ( \frac {\omega _{p}}{\omega _{s}}\right ) ^{2N} & =\frac {10^{-\frac {\delta _{p}}{10}}-1}{10^{-\frac {\delta _{s}}{10}}-1}\\ 2N\left [ \log \left ( \omega _{p}\right ) -\log \left ( \omega _{s}\right ) \right ] & =\log \left ( 10^{-\frac {\delta _{p}}{10}}-1\right ) -\log \left ( 10^{-\frac {\delta _{s}}{10}}-1\right ) \\ N & =\frac {1}{2}\frac {\log \left [ 10^{-\frac {\delta _{p}}{10}}-1\right ] -\log \left [ 10^{-\frac {\delta _{s}}{10}}-1\right ] }{\log \left ( \omega _{p}\right ) -\log \left ( \omega _{s}\right ) } \end {align*}

The above is later rounded to the nearest integer using the Ceiling function $$N=$$ $$\left \lceil N\right \rceil$$.

For impulse invariance method, using equation (1) above to solve for $$\Omega _{c}$$ gives\begin {align*} \left ( \frac {\omega _{p}/T}{\Omega _{c}}\right ) ^{2N} & =10^{-\frac {\delta _{p}}{10}}-1\\ 2N\left ( \log _{10}\frac {\omega _{p}/T}{\Omega _{c}}\right ) & =\log _{10}\left ( 10^{-\frac {\delta _{p}}{10}}-1\right ) \\ \log _{10}\frac {\omega _{p}/T}{\Omega _{c}} & =\frac {1}{2N}\log _{10}\left ( 10^{-\frac {\delta _{p}}{10}}-1\right ) \\ \frac {\omega _{p}/T}{\Omega _{c}} & =10^{\left ( \frac {1}{2N}\log _{10}\left ( 10^{-\frac {\delta _{p}}{10}}-1\right ) \right ) }\\ \Omega _{c} & =\frac {\omega _{p}/T}{10^{\left ( \frac {1}{2N}\log _{10}\left ( 10^{-\frac {\delta _{p}}{10}}-1\right ) \right ) }} \end {align*}

Now that $$N$$ and $$\Omega _{c}$$ are found, next the poles of $$H\left ( s\right )$$ are determined. Since the butterworth magnitude square of the transfer function is $\left \vert H_{a}\left ( s\right ) \right \vert ^{2}=\frac {T^{2}}{1+\left ( \frac {s}{j\Omega _{c}}\right ) ^{2N}}$ Then $$H\left ( s\right )$$ poles are found by setting the denominator of the above to zero which gives\begin {align*} 1+\left ( \frac {s}{j\Omega _{c}}\right ) ^{2N} & =0\\ \left ( \frac {s}{j\Omega _{c}}\right ) ^{2N} & =-1\\ & =e^{j\left ( \pi +2\pi k\right ) }\ \ \ \ \ k=0,1,2,\cdots 2N-1\\ \frac {s}{j\Omega _{c}} & =e^{j\left ( \frac {\pi +2\pi k}{2N}\right ) }\\ s & =j\Omega _{c}\ e^{j\left ( \frac {\pi +2\pi k}{2N}\right ) }\\ & =\Omega _{c}\ e^{j\frac {\pi }{2}}e^{j\left ( \frac {\pi +2\pi k}{2N}\right ) }\\ & =\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+2k+N\right ) }{2N}\right ) } \end {align*}

Only the LHS poles needs to be found, and these are located at $$i=0\cdots N-1$$. This is because these are the stable poles. Hence the $$i^{th}$$ pole is $s_{i}=\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) }$ For an, using $$i=0$$, $$N=6$$ gives\begin {align*} s_{0} & =\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+N\right ) }{2N}\right ) }\\ & =\Omega _{c}\ e^{j\left ( \frac {\pi 7}{12}\right ) } \end {align*}

Now the analog ﬁlter is generated based on the above selected poles, which for impulse invariance gives\begin {equation} H_{a}\left ( s\right ) =\frac {T\ K}{{\prod \limits _{i=0}^{N-1}}\left ( s-s_{i}\right ) }\tag {3} \end {equation} $$K$$ is found by solving $$H_{a}\left ( 0\right ) =T$$ which gives$k={\prod \limits _{i=0}^{N-1}}\left ( -s_{i}\right )$

The poles are written in non-polar form and substituted into (3) which gives

$s_{i}=\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) }=\Omega _{c}\left ( \cos \frac {\pi \left ( 1+2i+N\right ) }{2N}+j\sin \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) \qquad i=0\cdots N-1$ Therefore\begin {equation} H_{a}\left ( s\right ) =\frac {T\ K}{{\prod \limits _{i=0}^{N-1}}\left ( s-\Omega _{c}\left ( \cos \frac {\pi \left ( 1+2i+N\right ) }{2N}+j\sin \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) \right ) }\tag {4} \end {equation} Where $\Omega _{c}=\frac {\omega _{p}/T}{10^{\left ( \frac {1}{2N}\log _{10}\left ( 10^{-\frac {\delta _{p}}{10}}-1\right ) \right ) }}$ And$N=\left \lceil \frac {1}{2}\frac {\log \left [ 10^{-\frac {\delta _{p}}{10}}-1\right ] -\log \left [ 10^{-\frac {\delta _{s}}{10}}-1\right ] }{\log \left ( \omega _{p}\right ) -\log \left ( \omega _{s}\right ) }\right \rceil$ Now that $$H\left ( s\right )$$ is found, it is converted to $$H\left ( z\right ) .$$ We need to make sure that we multiply poles of complex conjugates with each others to make the result simple to see.

Now that $$H_{a}\left ( s\right )$$ is found, the $$A\rightarrow D$$ conversion is done. This means to obtain $$H\left ( z\right )$$ from the above $$H\left ( s\right )$$. When using impulse invariance, partial fraction decomposition is performed on (4) above in order to write $$H\left ( s\right )$$ as$H\left ( s\right ) ={\sum \limits _{i=0}^{N-1}}\frac {A_{i}}{s-s_{i}}$ For example, to obtain $$A_{j}$$ the result is$A_{j}=\lim _{s\rightarrow s_{j}}H_{a}\left ( s\right ) =\frac {T\ k}{{\prod \limits _{\substack {i=0\\i\neq j}}^{N-1}}\left ( s-s_{i}\right ) }$ Once the $$A^{\prime }s$$ are found, then $$H\left ( z\right )$$ becomes$H\left ( z\right ) ={\sum \limits _{i=0}^{N-1}}\frac {A_{i}}{1-\exp \left ( s_{i}\ T\right ) z^{-1}}$ This completes the design. The above form of $$H\left ( z\right )$$ could also be converted to rational expression as $$H\left ( z\right ) =\frac {N\left ( z\right ) }{D\left ( z\right ) }.$$

#### 2.2 Bilinear transformation method

First step is to convert analog speciﬁcations to digital speciﬁcations: $$\frac {F_{s}}{2\pi }=\frac {f_{p}}{\omega _{p}}$$, hence $$\omega _{p}=2\pi \frac {f_{p}}{F_{s}}$$ and $$\omega _{s}=2\pi \frac {f_{s}}{F_{s}}$$

Converting the criteria relative to the digital normalized scale gives\begin {align*} 20\log \left \vert H\left ( e^{j\omega _{p}}\right ) \right \vert & \geq \delta _{p}\\ 20\log \left \vert H\left ( e^{j\omega _{s}}\right ) \right \vert & \leq \delta _{s} \end {align*}

Hence \begin {align} \left \vert H\left ( e^{j\omega _{p}}\right ) \right \vert & \geq 10^{\frac {\delta _{p}}{20}}\tag {A}\\ \left \vert H\left ( e^{j\omega _{s}}\right ) \right \vert & \leq 10^{\frac {\delta _{s}}{20}}\tag {B} \end {align}

Butterworth analog ﬁlter squared magnitude Fourier transform is given by$\left \vert H_{a}\left ( j\Omega \right ) \right \vert ^{2}=\frac {1}{1+\left ( \frac {j\Omega }{j\Omega _{c}}\right ) ^{2N}}$ Equations (A) and (B) above are now written in terms of the analog Butterworth amplitude frequency response and become\begin {align*} \frac {1}{1+\left ( \frac {\Omega _{p}}{\Omega _{c}}\right ) ^{2N}} & \geq \left ( 10^{\frac {\delta _{p}}{20}}\right ) ^{2}=\ 10^{\frac {\delta _{p}}{10}}\\ \frac {1}{1+\left ( \frac {\Omega _{s}}{\Omega _{c}}\right ) ^{2N}} & \leq \left ( 10^{\frac {\delta _{s}}{20}}\right ) ^{2}=10^{\frac {\delta _{s}}{10}} \end {align*}

Values for $$\Omega _{p}$$ and $$\Omega _{s}$$ are assigned as follows: $$\Omega _{p}=\frac {2}{T}\tan \left ( \frac {\omega _{p}}{2}\right )$$, $$\Omega _{s}=\frac {2}{T}\tan \left ( \frac {\omega _{s}}{2}\right )$$. Hence the above becomes

\begin {align} 1+\left ( \frac {\frac {2}{T}\tan \left ( \frac {\omega _{p}}{2}\right ) }{\Omega _{c}}\right ) ^{2N} & \leq 10^{\frac {\delta _{p}}{10}}\tag {1}\\ 1+\left ( \frac {\frac {2}{T}\tan \left ( \frac {\omega _{s}}{2}\right ) }{\Omega _{c}}\right ) ^{2N} & \geq 10^{\frac {\delta _{s}}{10}}\tag {2} \end {align}

Changing inequalities to equalities and simplifying gives\begin {align*} \left ( \frac {\frac {2}{T}\tan \left ( \frac {\omega _{p}}{2}\right ) }{\Omega _{c}}\right ) ^{2N} & =10^{\frac {\delta _{p}}{10}}-1\\ \left ( \frac {\frac {2}{T}\tan \left ( \frac {\omega _{s}}{2}\right ) }{\Omega _{c}}\right ) ^{2N} & =10^{\frac {\delta _{s}}{10}}-1 \end {align*}

Dividing the above 2 equations results in\begin {align*} \left ( \frac {\tan \left ( \frac {\omega _{p}}{2}\right ) }{\tan \left ( \frac {\omega _{s}}{2}\right ) }\right ) ^{2N} & =\frac {10^{\frac {\delta _{p}}{10}}-1}{10^{\frac {\delta _{s}}{10}}-1}\\ 2N\left [ \log \left ( \tan \left ( \frac {\omega _{p}}{2}\right ) \right ) -\log \left ( \tan \left ( \frac {\omega _{s}}{2}\right ) \right ) \right ] & =\log \left ( 10^{\frac {\delta _{p}}{10}}-1\right ) -\log \left ( 10^{\frac {\delta _{s}}{10}}-1\right ) \\ N & =\frac {1}{2}\frac {\log \left ( 10^{\frac {\delta _{p}}{10}}-1\right ) -\log \left ( 10^{\frac {\delta _{s}}{10}}-1\right ) }{\log \left ( \tan \left ( \frac {\omega _{p}}{2}\right ) \right ) -\log \left ( \tan \left ( \frac {\omega _{s}}{2}\right ) \right ) } \end {align*}

The above is rounded to the nearest integer using the Ceiling function. i.e. $$N=$$ $$\left \lceil N\right \rceil .$$

For bilinear transformation equation (2) is used to ﬁnd $$\Omega _{c}.$$ Solving for $$\Omega _{c}$$ gives\begin {align*} 1+\left ( \frac {\frac {2}{T}\tan \left ( \frac {\omega _{s}}{2}\right ) }{\Omega _{c}}\right ) ^{2N} & =10^{\frac {\delta _{s}}{10}}\\ 2N\left ( \log _{10}\frac {\frac {2}{T}\tan \left ( \frac {\omega _{s}}{2}\right ) }{\Omega _{c}}\right ) & =\log _{10}\left ( 10^{\frac {\delta _{s}}{10}}-1\right ) \\ \log _{10}\frac {\frac {2}{T}\tan \left ( \frac {\omega _{s}}{2}\right ) }{\Omega _{c}} & =\frac {1}{2N}\log _{10}\left ( 10^{\frac {\delta _{s}}{10}}-1\right ) \\ \frac {\frac {2}{T}\tan \left ( \frac {\omega _{s}}{2}\right ) }{\Omega _{c}} & =10^{\left ( \frac {1}{2N}\log _{10}\left ( 10^{\frac {\delta _{s}}{10}}-1\right ) \right ) }\\ \Omega _{c} & =\frac {\frac {2}{T}\tan \left ( \frac {\omega _{s}}{2}\right ) }{10^{\left ( \frac {1}{2N}\log _{10}\left ( 10^{\frac {\delta _{s}}{10}}-1\right ) \right ) }} \end {align*}

Now that $$N$$ and $$\Omega _{c}$$ are found, the poles of $$H\left ( s\right )$$ can be determined. Since for bilinear the magnitude square of the transfer function is$\left \vert H_{a}\left ( s\right ) \right \vert ^{2}=\frac {1}{1+\left ( \frac {s}{j\Omega _{c}}\right ) ^{2N}}$ Hence $$H\left ( s\right )$$ poles are found by setting the denominator of the above to zero \begin {align*} 1+\left ( \frac {s}{j\Omega _{c}}\right ) ^{2N} & =0\\ \left ( \frac {s}{j\Omega _{c}}\right ) ^{2N} & =-1\\ & =e^{j\left ( \pi +2\pi k\right ) }\qquad k=0,1,2,\cdots 2N-1\\ \frac {s}{j\Omega _{c}} & =e^{j\left ( \frac {\pi +2\pi k}{2N}\right ) }\\ s & =j\Omega _{c}\ e^{j\left ( \frac {\pi +2\pi k}{2N}\right ) }\\ & =\Omega _{c}\ e^{j\frac {\pi }{2}}e^{j\left ( \frac {\pi +2\pi k}{2N}\right ) }\\ & =\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+2k+N\right ) }{2N}\right ) } \end {align*}

Only need the LHS poles are needed, which are located at $$i=0\cdots N-1$$, because these are the stable poles. Hence the $$i^{th}$$ pole is $s_{i}=\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) }$ For example for $$i=0$$, $$N=6$$ the result is$s_{0}=\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+N\right ) }{2N}\right ) }=\Omega _{c}\ e^{j\left ( \frac {\pi 7}{12}\right ) }$ For bilinear $$H(s)$$ is given by\begin {equation} H_{a}\left ( s\right ) =\frac {K}{{\prod \limits _{i=0}^{N-1}}\left ( s-s_{i}\right ) }\tag {3} \end {equation} $$K$$ is found by solving $$H_{a}\left ( 0\right ) =1$$ giving$k={\prod \limits _{i=0}^{N-1}}\left ( -s_{i}\right )$ The same expression results for $$k$$ in both cases. Writing poles in non-polar form and substituting them into (3) gives$s_{i}=\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) }=\Omega _{c}\left ( \cos \frac {\pi \left ( 1+2i+N\right ) }{2N}+j\sin \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) \ \ \ \ \ \ i=0\cdots N-1$ Then\begin {equation} H_{a}\left ( s\right ) =\frac {\ K}{{\prod \limits _{i=0}^{N-1}}\left ( s-\Omega _{c}\left ( \cos \frac {\pi \left ( 1+2i+N\right ) }{2N}+j\sin \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) \right ) }\tag {4} \end {equation} Where $\Omega _{c}=\frac {\frac {2}{T}\tan \left ( \frac {\omega _{s}}{2}\right ) }{10^{\left ( \frac {1}{2N}\log _{10}\left ( 10^{\frac {\delta _{s}}{10}}-1\right ) \right ) }}$ And$N=\left \lceil \frac {1}{2}\frac {\log _{10}\left ( 10^{\frac {\delta _{p}}{10}}-1\right ) -\log _{10}\left ( 10^{\frac {\delta _{s}}{10}}-1\right ) }{\log _{10}\left ( \tan \left ( \frac {\omega _{p}}{2}\right ) \right ) -\log _{10}\left ( \tan \left ( \frac {\omega _{s}}{2}\right ) \right ) }\right \rceil$ Now that $$H\left ( s\right )$$ is found, it is converted to $$H\left ( z\right )$$. After ﬁnding $$H(s)$$ as shown above then $$s$$ is replaced by $$\frac {2}{T}\frac {1-z^{-1}}{1+z^{-1}}$$. This is much simpler than the impulse invariance method. Before doing this substitution, we make sure to multiply poles which are complex conjugate of each others in the denominator of $$H(s)$$. After this, then the above substitution is made.

#### 2.3 Summary of analytical derivation method

Table with the derivation equations is made to follow to design in either bilinear or impulse invariance. Note that the same steps are used in both designs except for step $$5,6,8,13$$. This table make it easier to develop a program for the implementation

.

 step Impulse invariance common equation bilinear $$1$$ $$\omega _{p}=2\pi \frac {f_{p}}{F_{s}}$$ $$2$$ $$\omega _{s}=2\pi \frac {f_{s}}{F_{s}}$$ $$3$$ $$\alpha _{p}=\frac {1}{10^{\frac {\delta _{p}}{10}}}$$ $$4$$ $$\alpha _{s}=\frac {1}{10^{\frac {\delta s}{10}}}$$ $$5$$ $$\Omega _{p}=\frac {\omega _{p}}{T}$$ $$\frac {2}{T}\tan \left ( \frac {\omega _{p}}{2}\right )$$ $$6$$ $$\Omega _{s}=\frac {\omega _{s}}{T}$$ $$\frac {2}{T}\tan \left ( \frac {\omega _{s}}{2}\right )$$ $$7$$ $$N=\left \lceil \frac {1}{2}\frac {\log \left [ \alpha _{p}-1\right ] -\log \left [ \alpha _{s}-1\right ] }{\log \left ( \Omega _{p}\right ) -\log \left ( \Omega _{s}\right ) }\right \rceil$$ $$8$$ $$\Omega _{c}=\frac {\Omega _{p}}{10^{\left ( \frac {1}{2N}\log _{10}\left [ \alpha _{p}-1\right ] \right ) }}$$ $$\Omega _{c}=\frac {\Omega _{s}}{10^{\left ( \frac {1}{2N}\log _{10}\left [ \alpha _{s}-1\right ] \right ) }}$$ $$9$$ poles of H(S)$$\ s_{i}=\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) }\qquad i=0\cdots N-1$$ $$10$$ $$k={\prod \limits _{i=0}^{N-1}}\left ( -s_{i}\right )$$ $$11$$ $$H_{a}\left ( s\right ) =\frac {T\ K}{{\prod \limits _{i=0}^{N-1}}\left ( s-s_{i}\right ) }$$ $$H_{a}\left ( s\right ) =\frac {K}{{\prod \limits _{i=0}^{N-1}}\left ( s-s_{i}\right ) }$$ do partial fractions: $$12$$ $$H_{a}\left ( s\right ) ={\sum \limits _{i=0}^{N-1}}\frac {A_{i}}{s-s_{i}}$$ $$13$$ $$H\left ( z\right ) ={\sum \limits _{i=0}^{N-1}}\frac {A_{i}}{1-\exp \left ( Ts_{i}\right ) z^{-1}}$$ $$H\left ( z\right ) =\left . H_{a}\left ( s\right ) \right \vert _{s=\frac {2}{T}\frac {1-z^{-1}}{1+z^{-1}}}$$

### 3 Numerical design examples

#### 3.1 Example 1

Sampling frequency $$F_{s}=20khz$$, passband frequency $$f_{p}=2khz$$, stopband frequency $$f_{s}=3khz$$, with $$\delta _{p}\geq -1db$$ and $$\delta _{stop}\leq -15db$$

##### 3.1.1 using impulse invariance method (using T=1)

 step Impulse invariance $$1$$ $$\omega _{p}=2\pi \frac {f_{p}}{F_{s}}$$ $$\rightarrow \frac {2\pi \left ( 2000\right ) }{20000}\rightarrow 0.2\pi$$ $$2$$ $$\omega _{s}=2\pi \frac {f_{s}}{F_{s}}\rightarrow \frac {2\pi \left ( 3000\right ) }{20000}\rightarrow 0.3\pi$$ $$3$$ $$\alpha _{p}=\frac {1}{10^{\frac {\delta _{p}}{10}}}\rightarrow \frac {1}{10^{\frac {-1}{10}}}\rightarrow$$ $$1.258\,9$$ $$4$$ $$\alpha _{s}=\frac {1}{10^{\frac {\delta s}{10}}}\rightarrow \frac {1}{10^{\frac {-15}{10}}}\rightarrow$$ $$31.623$$ $$5$$ $$\Omega _{p}=\frac {\omega _{p}}{T}\rightarrow \frac {\omega _{p}}{1}\rightarrow 0.2\pi$$ $$6$$ $$\Omega _{s}=\frac {\omega _{s}}{T}\rightarrow \frac {0.3\pi }{1}\rightarrow 0.3\pi$$ $$7$$ $$N=\left \lceil \frac {1}{2}\frac {\log \left [ \alpha _{p}-1\right ] -\log \left [ \alpha _{s}-1\right ] }{\log \left ( \Omega _{p}\right ) -\log \left ( \Omega _{s}\right ) }\right \rceil \rightarrow \frac {1}{2}\frac {\log _{10}\left ( 1.258\,9-1\right ) -\log _{10}\left ( 31.623-1\right ) }{\log _{10}\left ( 0.2\pi \right ) -\log _{10}\left ( 0.3\pi \right ) }\rightarrow$$ $$\left \lceil 5.885\,9\right \rceil \rightarrow 6$$ $$8$$ $$\Omega _{c}=\frac {\Omega _{p}}{10^{\left ( \frac {1}{2N}\log _{10}\left [ \alpha _{p}-1\right ] \right ) }}\rightarrow \frac {0.2\pi }{10^{\left ( \frac {1}{2\times 6}\log _{10}\left ( 1.258\,9-1\right ) \right ) }}$$ $$\rightarrow 0.703\,21$$ $$9$$ poles of H(S)$$\ s_{i}=\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) }\rightarrow s_{i}=0.703\,21\ e^{j\left ( \frac {\pi \left ( 7+2i\right ) }{12}\right ) }i=0\cdots 5$$ $$s_{0}=-0.182\,+j0.679\,25,s_{1}=-0.497\,24+j0.497\,2,s_{2}=-0.679\,25+j0.182\,$$ $$s_{3}=-0.679\,25-j0.182\,,s_{4}=-0.497\,24-j0.497\,24,s_{5}=$$ $$-0.182\,-j0.679\,25$$ $$10$$ $$k={\prod \limits _{i=0}^{N-1}}\left ( -s_{i}\right ) \rightarrow \left ( 0.182\,-j0.679\,25\right ) \left ( 0.497\,24-j0.497\,24\right ) \left ( 0.679\,25-j0.182\,\right )$$ $$\left ( 0.679\,25+j0.182\right ) \left ( 0.497\,24+j0.497\,24\right ) \left ( 0.182\,+j0.679\,25\right ) \rightarrow 0.120\,92$$ $$11$$ $$H_{a}\left ( s\right ) =\frac {T\ K}{{\prod \limits _{i=0}^{N-1}}\left ( s-s_{i}\right ) }\rightarrow$$ $$\frac {1\times \allowbreak 0.120\,92}{\left ( s+0.182\,-j0.679\,25\right ) \left ( s+0.497\,24-j0.497\,24\right ) \left ( s+0.679\,25-j0.182\right ) \left ( s+0.679\,25+j0.182\right ) \left ( s+0.497\,24+j0.497\,24\right ) \left ( s+0.182\,+j0.679\,25\right ) }$$ $$\rightarrow$$multiply complex conjugates$$\rightarrow \frac {\allowbreak 0.120\,92}{\left ( s^{2}+0.364\,s+0.494\,5\right ) \left ( s^{2}+0.994\,48s+0.494\,50\right ) \left ( s^{2}+1.358\,5s+0.494\,5\right ) }$$ $$\rightarrow \frac {\allowbreak 0.120\,92}{s^{6}+2.717\,0s^{5}+3.691\,0s^{4}+3.178\,9\allowbreak s^{3}+1.825\,2s^{2}+0.664\,38\allowbreak s+0.120\,92}$$ $$12$$ partial fraction $$H_{a}\left ( s\right ) ={\sum \limits _{i=0}^{N-1}}\frac {A_{i}}{s-s_{i}}\rightarrow \frac {0.143\,54+j0.248\,61}{s+0.182\,-j0.679\,25}-\frac {1.071\,4-j1.166\,8\times 10^{-5}}{s+0.497\,24+j0.497\,24}+\frac {0.927\,85-j1.607\,1}{s+0.679\,25-j0.182\,}$$ $$+\frac {0.927\,85+j1.607\,1}{s+0.679\,25+j0.182\,}-\frac {1.071\,4+j1.166\,8\times 10^{-5}}{s+0.497\,24-j0.497\,24}+\frac {0.143\,54-j0.248\,61}{s+0.182\,+j0.679\,25}$$ $$13$$ $$H\left ( z\right ) ={\sum \limits _{i=0}^{N-1}}\frac {A_{i}}{1-\exp \left ( Ts_{i}\right ) z^{-1}}\rightarrow \frac {0.143\,54+j0.248\,61}{1-\exp \left ( -0.182\,+j0.679\,25\right ) z^{-1}}-$$ $$\frac {1.071\,4-j1.166\,8\times 10^{-5}}{1-\exp \left ( -0.497\,24-j0.497\,24\right ) z^{-1}}+\frac {0.927\,85-j1.607\,1}{1-\exp \left ( -0.679\,25+j0.182\right ) \,}$$ $$+\frac {0.927\,85+j1.607\,1}{1-\exp \left ( -0.679\,25-j0.182\right ) \,z^{-1}}-\frac {1.071\,4+j1.166\,8\times 10^{-5}}{1-\exp \left ( -0.497\,24+j0.497\,24\right ) z^{-1}}+\frac {0.143\,54-j0.248\,61}{1-\exp \left ( -0.182\,-j0.679\,25\right ) z^{-1}}$$ $$H\left ( z\right ) =\frac {0.143\,54+j0.248\,61}{1-(0.648\,58+j0.523\,68)z^{-1}}$$ $$-\frac {1.071\,4-j1.166\,8\times 10^{-5}}{1-\left ( 0.534\,55-j0.290\,12\right ) z^{-1}}+\frac {0.927\,85-j1.607\,1}{1-\left ( 0.498\,62+j9.176\,5\times 10^{-2}\right ) z^{-1}\,}$$ $$+\frac {0.927\,85+j1.607\,1}{1-\left ( 0.498\,62-j9.176\,5\times 10^{-2}\right ) \,z^{-1}}-\frac {1.071\,4+j1.166\,8\times 10^{-5}}{1-\left ( 0.534\,55+j0.29012\right ) z^{-1}}$$ $$+\frac {0.143\,54-j0.248\,61}{1-\left ( 0.648\,58-j0.52368\right ) z^{-1}}$$

##### 3.1.2 Bilinear method

Using the above design table, these are the numerical values: $$T=\frac {1}{F_{s}}=\frac {1}{20000}$$

 step Bilinear $$1$$ $$\omega _{p}=2\pi \frac {f_{p}}{F_{s}}$$ $$\rightarrow \frac {2\pi \left ( 2000\right ) }{20000}\rightarrow 0.2\pi$$ $$2$$ $$\omega _{s}=2\pi \frac {f_{s}}{F_{s}}\rightarrow \frac {2\pi \left ( 3000\right ) }{20000}\rightarrow 0.3\pi$$ $$3$$ $$\alpha _{p}=\frac {1}{10^{\frac {\delta _{p}}{10}}}\rightarrow \frac {1}{10^{\frac {-1}{10}}}\rightarrow$$ $$1.258\,9$$ $$4$$ $$\alpha _{s}=\frac {1}{10^{\frac {\delta s}{10}}}\rightarrow \frac {1}{10^{\frac {-15}{10}}}\rightarrow$$ $$31.623$$ $$5$$ $$\Omega _{p}=\frac {2}{T}\tan \left ( \frac {\omega _{p}}{2}\right ) \rightarrow 2\times 20000\tan \left ( \frac {0.2\pi }{2}\right ) \rightarrow 12997$$ $$6$$ $$\Omega _{s}=\frac {2}{T}\tan \left ( \frac {\omega _{s}}{2}\right ) \rightarrow 2\times 20000\tan \left ( \frac {0.3\pi }{2}\right ) \rightarrow 20381$$ $$7$$ $$N=\left \lceil \frac {1}{2}\frac {\log \left [ \alpha _{p}-1\right ] -\log \left [ \alpha _{s}-1\right ] }{\log \left ( \Omega _{p}\right ) -\log \left ( \Omega _{s}\right ) }\right \rceil \rightarrow \frac {1}{2}\frac {\log _{10}\left ( 1.258\,9-1\right ) -\log _{10}\left ( 31.623-1\right ) }{\log _{10}\left ( 12997\right ) -\log _{10}\left ( 20381\right ) }\rightarrow$$ $$\left \lceil 5.304\,8\right \rceil \rightarrow 6$$ $$8$$ $$\Omega _{c}=\frac {\Omega _{s}}{10^{\left ( \frac {1}{2N}\log _{10}\left [ \alpha _{s}-1\right ] \right ) }}\rightarrow \frac {20381}{10^{\left ( \frac {1}{2\times 6}\log _{10}\left ( 31.623-1\right ) \right ) }}\rightarrow$$ $$15325.6$$ $$9$$ poles of H(S)$$\ s_{i}=\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) }\rightarrow s_{i}=15325\ e^{j\left ( \frac {\pi \left ( 7+2i\right ) }{12}\right ) }\qquad i=0\cdots 5$$ $$s_{0}=-3966.4+j14803,s_{1}=-10836+j10836.,s_{2}=-14803+j3966.4$$ $$s_{3}=-14803-j3966.4,s_{4}=-10836-j10836.,s_{5}=-3966.4-j14803$$ $$10$$ $$k={\prod \limits _{i=0}^{N-1}}\left ( -s_{i}\right ) \rightarrow \left ( 3966.4-j14803\right ) \left ( 10836-j10836\right ) \left ( 14803-j3966.4\right )$$ $$\left ( 14803+j3966.4\right ) \left ( 10836+j10836\right ) \left ( 3966.4+j14803\right )$$ $$\rightarrow$$  multiply complex conjugate terms $$\rightarrow \left ( 2.348\,6\times 10^{8}\right ) \left ( 2.348\,4\times 10^{8}\right )$$ $$\left ( 2.348\,6\times 10^{8}\right )$$ $$\rightarrow$$ $$1.295\,4\times 10^{25}$$ $$11$$ $$H_{a}\left ( s\right ) =\frac {K}{{\prod \limits _{i=0}^{N-1}}\left ( s-s_{i}\right ) }\rightarrow$$ $$\frac {1.295\,4\times 10^{25}}{\left ( s+3966.4-j14803\right ) \left ( s+10836-j10836\right ) \left ( s+14803-j3966.4\allowbreak \right ) \left ( s+14803+j3966.4\right ) \left ( s+10836+j10836\right ) \left ( s+3966.4+j14803\right ) }$$ $$\rightarrow$$  multiply complex conjugate $$\frac {1.295\,4\times 10^{25}}{\left ( s^{2}+7932.8s+2.348\,6\times 10^{8}\right ) \left ( s^{2}+21\,672s+234\,837\,792\right ) \left ( s^{2}+29606.s+2.348\,6\times 10^{8}\right ) }$$ $$=\frac {6257.4s-1.355\,8\times 10^{8}}{s^{2}+7932.8s+2.348\,6\times 10^{8}}-\frac {46702.s+5.060\,4\times 10^{8}}{s^{2}+21672.s+2.348\,4\times 10^{8}}+\frac {40445.s+8.765\,4\times 10^{8}}{s^{2}+29606.s+2.348\,6\times 10^{8}}\rightarrow$$ $$\rightarrow \frac {1.295\,4\times 10^{25}}{s^{6}+59211.s^{5}+1.753\,0\times 10^{9}s^{4}+3.290\,2\times 10^{13}s^{3}+4.116\,9\times 10^{17}\allowbreak s^{2}+3.265\,8\times 10^{21}s+1.295\,3\times 10^{25}\allowbreak }$$ $$12$$ $$\allowbreak$$ $$13$$ $$H\left ( z\right ) =\left . H_{a}\left ( s\right ) \right \vert _{s=\frac {2}{T}\frac {1-z^{-1}}{1+z^{-1}}}\rightarrow TODO$$ $$\cdots$$

#### 3.2 Example 2

Sampling frequency $$F_{s}=10khz$$, passband corner frequency $$f_{p}=1khz$$, stopband corner frequency $$f_{s}=2khz$$, with criteria $$\delta _{p}\geq -3db$$ and $$\delta _{stop}\leq -10db$$

##### 3.2.1 using impulse invariance method

$$T=1$$

 step Impulse invariance $$1$$ $$\omega _{p}=2\pi \frac {f_{p}}{F_{s}}$$ $$\rightarrow \frac {2\pi \left ( 1000\right ) }{10000}\rightarrow \allowbreak 0.2\pi$$ $$2$$ $$\omega _{s}=2\pi \frac {f_{s}}{F_{s}}\rightarrow \frac {2\pi \left ( 2000\right ) }{10000}\rightarrow \allowbreak 0.4\pi$$ $$3$$ $$\alpha _{p}=\frac {1}{10^{\frac {\delta _{p}}{10}}}\rightarrow \frac {1}{10^{\frac {-3}{10}}}\rightarrow$$ $$1.995\,3$$ $$4$$ $$\alpha _{s}=\frac {1}{10^{\frac {\delta s}{10}}}\rightarrow \frac {1}{10^{\frac {-10}{10}}}\rightarrow$$ $$10.0$$ $$5$$ $$\Omega _{p}=\frac {\omega _{p}}{T}\rightarrow \frac {\omega _{p}}{1}\rightarrow 0.2\pi$$ $$6$$ $$\Omega _{s}=\frac {\omega _{s}}{T}\rightarrow \frac {\allowbreak \allowbreak 0.4\pi }{1}\rightarrow \allowbreak 0.4\pi$$ $$7$$ $$N=\left \lceil \frac {1}{2}\frac {\log \left [ \alpha _{p}-1\right ] -\log \left [ \alpha _{s}-1\right ] }{\log \left ( \Omega _{p}\right ) -\log \left ( \Omega _{s}\right ) }\right \rceil \rightarrow \frac {1}{2}\frac {\log _{10}\left ( 1.995\,3-1\right ) -\log _{10}\left ( 10.0-1\right ) }{\log _{10}\left ( 0.2\pi \right ) -\log _{10}\left ( 0.4\pi \right ) }\rightarrow$$ $$1.588\,4\rightarrow$$ $$2$$ $$8$$ $$\Omega _{c}=\frac {\Omega _{p}}{10^{\left ( \frac {1}{2N}\log _{10}\left [ \alpha _{p}-1\right ] \right ) }}\rightarrow \frac {0.2\pi }{10^{\left ( \frac {1}{2\times 2}\log _{10}\left ( 1.995\,3-1\right ) \right ) }}\rightarrow$$ $$0.629\,06$$ $$9$$ poles of H(S)$$\ s_{i}=\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) }\rightarrow s_{i}=0.629\,06\ e^{j\left ( \frac {\pi \left ( 3+2i\right ) }{4}\right ) }\ \ i=0\cdots 1$$ $$s_{0}=-0.444\,81+j0.444\,81\allowbreak ,s_{1}=-0.444\,81-j0.444\,81\allowbreak$$ $$10$$ $$k={\prod \limits _{i=0}^{N-1}}\left ( -s_{i}\right ) \rightarrow \left ( 0.444\,81-j0.444\,81\right ) \left ( 0.444\,81+j0.444\,81\allowbreak \right ) \rightarrow \allowbreak 0.395\,71$$ $$11$$ $$H_{a}\left ( s\right ) =\frac {T\ K}{{\prod \limits _{i=0}^{N-1}}\left ( s-s_{i}\right ) }\rightarrow \frac {\allowbreak 0.395\,71}{\left ( s+0.444\,81-j0.444\,81\right ) \left ( s+0.444\,81+j0.444\,81\allowbreak \right ) }\rightarrow \frac {\allowbreak 0.395\,71}{s^{2}+0.889\,62s+0.395\,71}$$ $$12$$ partial fraction $$H_{a}\left ( s\right ) ={\sum \limits _{i=0}^{N-1}}\frac {A_{i}}{s-s_{i}}\rightarrow \allowbreak \frac {0.444\,81j}{s+0.444\,81+0.444\,81j}-\frac {\allowbreak 0.444\,81j}{s+0.444\,81-0.444\,81j}\allowbreak$$ $$\rightarrow \allowbreak \frac {0.245\,35z}{z^{2}-1. 157\,2z+0.410\,81}$$ $$13$$ $$H\left ( z\right ) ={\displaystyle \sum \limits _{i=0}^{N-1}} \frac {A_{i}}{1-\exp \left ( Ts_{i}\right ) z^{-1}}\rightarrow \frac {\allowbreak 0.444\,81i}{1-\exp \left ( -0.444\,81-j0.444\,81\right ) z^{-1}}-\allowbreak \frac {0.444\,81i}{1-\exp \left ( -0.444\,81+j0.444\,81\right ) z^{-1}}\allowbreak$$

##### 3.2.2 Using bilinear

$$T=\frac {1}{10000}$$

 step Impulse invariance $$1$$ $$\omega _{p}=2\pi \frac {f_{p}}{F_{s}}$$ $$\rightarrow \frac {2\pi \left ( 1000\right ) }{10000}\rightarrow \allowbreak 0.2\pi$$ $$2$$ $$\omega _{s}=2\pi \frac {f_{s}}{F_{s}}\rightarrow \frac {2\pi \left ( 2000\right ) }{10000}\rightarrow \allowbreak 0.4\pi$$ $$3$$ $$\alpha _{p}=\frac {1}{10^{\frac {\delta _{p}}{10}}}\rightarrow \frac {1}{10^{\frac {-3}{10}}}\rightarrow$$ $$1.995\,3$$ $$4$$ $$\alpha _{s}=\frac {1}{10^{\frac {\delta s}{10}}}\rightarrow \frac {1}{10^{\frac {-10}{10}}}\rightarrow$$ $$10.0$$ $$5$$ $$\Omega _{p}=\frac {2}{T}\tan \left ( \frac {\omega _{p}}{2}\right ) \rightarrow 2\times 10000\tan \left ( \frac {0.2\pi }{2}\right ) \rightarrow$$ $$6498.4$$ $$6$$ $$\Omega _{s}=\frac {2}{T}\tan \left ( \frac {\omega _{s}}{2}\right ) \rightarrow 2\times 10000\tan \left ( \frac {0.4\pi }{2}\right ) \rightarrow 14531.$$ $$7$$ $$N=\left \lceil \frac {1}{2}\frac {\log \left [ \alpha _{p}-1\right ] -\log \left [ \alpha _{s}-1\right ] }{\log \left ( \Omega _{p}\right ) -\log \left ( \Omega _{s}\right ) }\right \rceil \rightarrow \frac {1}{2}\frac {\log _{10}\left ( 1.995\,3-1\right ) -\log _{10}\left ( 10.0-1\right ) }{\log _{10}\left ( 6498.4\right ) -\log _{10}\left ( 14531.\right ) }\rightarrow$$ $$1.368\,1$$ $$\rightarrow 2$$ $$8$$ $$\Omega _{c}=\frac {\Omega _{s}}{10^{\left ( \frac {1}{2N}\log _{10}\left [ \alpha _{s}-1\right ] \right ) }}\rightarrow \frac {14531}{10^{\left ( \frac {1}{2\times 2}\log _{10}\left ( 10.0-1\right ) \right ) }}\rightarrow$$ $$8389.5$$ $$9$$ poles of H(S)$$\ s_{i}=\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) }\rightarrow s_{i}=8389. 5\ e^{j\left ( \frac {\pi \left ( 3+2i\right ) }{4}\right ) }\ \ i=0\cdots 1$$ $$s_{0}=-5932.3+j5932.3,s_{1}=-5932. 3-j5932.3\allowbreak$$ $$10$$ $$k={\displaystyle \prod \limits _{i=0}^{N-1}} \left ( -s_{i}\right ) \rightarrow \left ( 5932. 3-j5932.3\right ) \left ( 5932.3+j5932. 3\allowbreak \right ) \rightarrow \allowbreak 7.038\,4\times 10^{7}$$ $$11$$ $$H_{a}\left ( s\right ) =\frac {K}{{\prod \limits _{i=0}^{N-1}}\left ( s-s_{i}\right ) }\rightarrow \frac {7.038\,4\times 10^{7}}{\left ( s+5932.3-i5932.3\right ) \left ( s+5932.3+i5932.3\allowbreak \allowbreak \right ) }$$ $$12$$ $$\rightarrow$$ $$\frac {7.038\,4\times 10^{7}}{s^{2}+11865.s+7.038\,4\times 10^{7}}$$ $$13$$ $$H\left ( z\right ) =\left . H_{a}\left ( s\right ) \right \vert _{s=\frac {2}{T}\frac {1-z^{-1}}{1+z^{-1}}}\rightarrow \frac {\allowbreak 0.09945\,9z^{2}+0.198\,92z+0.09945\,9}{z^{2}-0.931\,56z+0.329\,38}$$ $$\allowbreak$$

Mathematica GUI program is written to implement the above. It can be downloaded from here

Also a Matlab script nma_filter.m.txt was written (no GUI) to implement this design.

This script (written on Matlab 2007a). This script does not handle the conversion from H(s) to H(z) well yet, need to work more on this... of course, one can just use Matlab butter() function for this.

Example output of the above Matlab script is matlab_output.txt

### 4 IIR design for minimum order ﬁlter

This is another small note on IIR design for minimum order ﬁlter.

This document describes how to design an IIR digital ﬁlter given a speciﬁcation in which the ﬁlter order is speciﬁed.

Given the following diagram, the speciﬁcations for the design will be

1. Digital ﬁlter order (N)
2. $$f_{pass}$$, the passband corner frequency (or the cutoﬀ frequency) at $$-3db$$. This means $$A_{pass}=-3db$$

#### 4.1 Impulse invariance method

We ﬁrst convert analog speciﬁcations to digital speciﬁcations: $$\frac {F_{s}}{2\pi }=\frac {f_{pass}}{\omega _{p}}$$, hence $$\omega _{p}=2\pi \frac {f_{pass}}{F_{s}}$$

Then the cutoﬀ frequency $$\Omega _{c}=\frac {\omega _{p}}{T}$$ for impulse invariance.

Next ﬁnd the poles of $$H\left ( s\right ) .$$ Since the butterworth magnitude square of the transfer function is $\left \vert H_{a}\left ( s\right ) \right \vert ^{2}=\frac {T^{2}}{1+\left ( \frac {s}{j\Omega _{c}}\right ) ^{2N}}$ Hence $$H\left ( s\right )$$ poles are found by setting the denominator of the above to zero$1+\left ( \frac {s}{j\Omega _{c}}\right ) ^{2N}=0$ And as I did in the earlier document, the poles of $$H\left ( s\right ) \,$$ are found at$s=\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+2k+N\right ) }{2N}\right ) }$ We only need to ﬁnd the LHS poles, which are located at $$i=0\cdots N-1$$, because these are the stable poles. Hence the $$i^{th}$$ pole is $s_{i}=\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) }$ Now we can write the analog ﬁlter generated based on the above selected poles, which is, for impulse invariance\begin {equation} H_{a}\left ( s\right ) =\frac {T\ K}{{\displaystyle \prod \limits _{i=0}^{N-1}} \left ( s-s_{i}\right ) } \tag {3} \end {equation} $$K$$ is found by solving $$H_{a}\left ( 0\right ) =T$$  hence$k={\displaystyle \prod \limits _{i=0}^{N-1}} \left ( -s_{i}\right )$ Now we need to write poles in non-polar form and plug them into (3)$s_{i}=\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) }=\Omega _{c}\left ( \cos \frac {\pi \left ( 1+2i+N\right ) }{2N}+j\sin \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) \qquad i=0\cdots N-1$ Hence, \begin {equation} H_{a}\left ( s\right ) =\frac {T\ K}{{\displaystyle \prod \limits _{i=0}^{N-1}} \left ( s-\Omega _{c}\left ( \cos \frac {\pi \left ( 1+2i+N\right ) }{2N}+j\sin \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) \right ) } \tag {4} \end {equation} Now that we have found $$H\left ( s\right )$$ we need to convert it to $$H\left ( z\right )$$

We need to make sure that we multiply poles of complex conjugates with each others to make the result simple to see.

Now that we have $$H_{a}\left ( s\right )$$, we do the $$A\rightarrow D$$ conversion. I.e. obtain $$H\left ( z\right )$$ from the above $$H\left ( s\right )$$. When using impulse invariance, we need to perform partial fraction decomposition on (4) above in order to write $$H\left ( s\right )$$ in this form$H\left ( s\right ) ={\displaystyle \sum \limits _{i=0}^{N-1}} \frac {A_{i}}{s-s_{i}}$ For example, to obtain $$A_{j}$$, we write$A_{j}=\lim _{s\rightarrow s_{j}}H_{a}\left ( s\right ) =\frac {T\ k}{{\displaystyle \prod \limits _{\substack {i=0\\i\neq j}}^{N-1}} \left ( s-s_{i}\right ) }$ Once we ﬁnd all the $$A^{\prime }s$$, we now write $$H\left ( z\right )$$ as follows$H\left ( z\right ) ={\displaystyle \sum \limits _{i=0}^{N-1}} \frac {A_{i}}{1-\exp \left ( s_{i}\ T\right ) z^{-1}}$ This completes the design. We can try to convert the above form of $$H\left ( z\right )$$ to a rational expression as $$H\left ( z\right ) =\frac {N\left ( z\right ) }{D\left ( z\right ) }$$

#### 4.2 bilinear transformation method

We ﬁrst convert analog speciﬁcations to digital speciﬁcations: $$\frac {F_{s}}{2\pi }=\frac {f_{pass}}{\omega _{p}}$$, hence $$\omega _{p}=2\pi \frac {f_{pass}}{F_{s}}$$

Now $$\Omega _{c}$$ $$=\frac {2}{T}\tan \left ( \frac {\omega _{p}}{2}\right )$$

Now that we have $$N$$ and $$\Omega _{c}$$ we ﬁnd the poles of $$H\left ( s\right )$$. Since for bilinear the magnitude square of the transfer function is$\left \vert H_{a}\left ( s\right ) \right \vert ^{2}=\frac {1}{1+\left ( \frac {s}{j\Omega _{c}}\right ) ^{2N}}$ Hence $$H\left ( s\right )$$ poles are found by setting the denominator of the above to zero $1+\left ( \frac {s}{j\Omega _{c}}\right ) ^{2N}=0$ Which leads to poles at $$s_{i}=\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) }$$. We only need to ﬁnd the LHS poles, which are located at $$i=0\cdots N-1$$, because these are the stable poles. For bilinear, $$H(s)$$ is given by\begin {equation} H_{a}\left ( s\right ) =\frac {K}{{\displaystyle \prod \limits _{i=0}^{N-1}} \left ( s-s_{i}\right ) } \tag {3} \end {equation} $$K$$ is found by solving $$H_{a}\left ( 0\right ) =1$$, hence we obtain$k={\displaystyle \prod \limits _{i=0}^{N-1}} \left ( -s_{i}\right )$ We see that the same expression results for $$k$$ for both cases.

Now we need to write poles in non-polar form and plug them into (3)$s_{i}=\Omega _{c}\ e^{j\left ( \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) }=\Omega _{c}\left ( \cos \frac {\pi \left ( 1+2i+N\right ) }{2N}+j\sin \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) \qquad i=0\cdots N-1$ then\begin {equation} H_{a}\left ( s\right ) =\frac {\ K}{{\displaystyle \prod \limits _{i=0}^{N-1}} \left ( s-\Omega _{c}\left ( \cos \frac {\pi \left ( 1+2i+N\right ) }{2N}+j\sin \frac {\pi \left ( 1+2i+N\right ) }{2N}\right ) \right ) }\tag {4} \end {equation} Now that we have found $$H\left ( s\right )$$ we need to convert it to $$H\left ( z\right )$$. After ﬁnding $$H(s)$$ as shown above, we simply replace $$s$$ by $$\frac {2}{T}\frac {1-z^{-1}}{1+z^{-1}}$$. This is much simpler than the impulse invariance method. Before doing this substitution, make sure to multiply poles which are complex conjugate of each others in the denominator of $$H(s)$$. After this, then do the above substitution

### 5 Appendix

The following is listing of the Matlab script, and a sample run output

#### 5.1 Code listing

%Simple matlab script to design an IIR low pass filter using
%butterworth in either impulse inv. or bilinear method.

%
% EE 420, CSU Fullerton
%
% by Nasser M. Abbasi
% 5/5/201

% Filter SPECIFICATIONS
clear all;
close all;
Fs=10000;
fp=1000;
fs=2000;
dbp=-3;
dbs=-10;
BILINEAR=1;  %set this to 0 to do impulse inv.

%%%%%%%%%%%%%%%%%%%%

fprintf('***** STARTING DESIGN *******\n');
fprintf('Sampling frequency=%f Hz\n',Fs);
fprintf('freq at passband=%f Hz\n',fp);
fprintf('freq at stopband=%f Hz\n',fs);
fprintf('db at passband=%f \n',dbp);
fprintf('db at stopband=%f \n',dbs);
if BILINEAR
fprintf('Doing Bilinear method\n');
else
fprintf('Doing impulse invariance method\n');
end

if BILINEAR
T=1/Fs;
else
T=1;
end

fprintf('T=%f\n',T);

wp=2*pi*fp/Fs;
ws=2*pi*fs/Fs;
alphap=1/(10^(dbp/10));
alphas=1/(10^(dbs/10));

if BILINEAR
gammap=2/T * tan(wp/2);
else
gammap=wp/T;
end

if BILINEAR
gammas=2/T * tan(ws/2);
else
gammas=ws/T;
end

oldn=.5*(log10(alphap-1)-log10(alphas-1))/(log10(gammap)-log10(gammas));
n=ceil(oldn);
fprintf('n=%f, rounded to %d\n',oldn,n);

if BILINEAR
gammaC=gammas/(10^(1/(2*n)*log10(alphas-1)));
else
gammaC=gammap/(10^(1/(2*n)*log10(alphap-1)));
end

fprintf('Gamma C =%f\n',gammaC);

poles_of_hs=zeros(n,1);
for i=0:n-1
poles_of_hs(i+1)=gammaC*exp(sqrt(-1)*(pi*(1+2*i+n)/(2*n)));
end

fprintf('POLES Of H(s)\n');
poles_of_hs

k=prod(-poles_of_hs);
den=poly(poles_of_hs);

fprintf('k=%d\n',k);

fprintf('H(s)=\n',k);

if BILINEAR
hs=tf(k,den)
else
hs=tf(T*k,den)
end

[r,p,k]=residue(k,den);
hzp=zeros(n,1);
for i=1:n
hzp(i)=exp(p(i)*T);
end

fprintf('POLES Of H(z)\n');
hzp

[B,A]=residue(r,hzp);
fprintf('H(z)=\n',k);
hz=tf(T*B',A')


#### 5.2 Sample run output

***** STARTING DESIGN *******
Sampling frequency=20000.000000 Hz
freq at passband=2000.000000 Hz
freq at stopband=3000.000000 Hz
db at passband=-1.000000
db at stopband=-15.000000
Doing impulse invariance method
T=1.000000
n=5.885783, rounded to 6
Gamma C =0.703205
POLES Of H(s)

poles_of_hs =

-0.182002858631059 +     0.679243915533887i
-0.497241056902828 +     0.497241056902828i
-0.679243915533887 +     0.182002858631059i
-0.679243915533887 -     0.182002858631059i
-0.497241056902828 -     0.497241056902828i
-0.182002858631059 -     0.679243915533887i

k=1.209183e-001
H(s)=
Warning: Transfer function has complex coefficients.
> In tf.tf at 246
In nma_filter at 92

Transfer function:

(0.1209-1.041e-016i)
----------------------------------------------------------------------------
s^6 + (2.717-4.441e-016i) s^5 + (3.691-1.998e-015i) s^4 + (3.179

-2.22e-015i) s^3 + (1.825-1.554e-015i) s^2 + (0.6644-4.996e-016i) s

+ (0.1209-1.041e-016i)

POLES Of H(z)

hzp =

0.534553736986506 +     0.290115961427623i
0.534553736986506 -     0.290115961427623i
0.648579932539211 +     0.523670977796743i
0.648579932539211 -     0.523670977796743i
0.498626135868541 -    0.0917668874413562i
0.498626135868543 +    0.0917668874413561i

H(z)=
Warning: Transfer function has complex coefficients.
> In tf.tf at 246
In nma_filter at 107

Transfer function:

-(0.9938-0.8577i) s^4 - (0.22-0.703i) s^3 + (1.377+1.095i) s^2 - (0.6574

+2.854i) s - (0.9146-1.114i)

--------------------------------------------------------------------------
-(0.6979+1.411i) s^4 + (0.6009-0.8998i) s^3 - (0.03618-0.9428i) s^2

+ (0.1901+0.7751i) s - (0.6018+0.246i)

=============================================================

***** STARTING DESIGN *******
Sampling frequency=20000.000000 Hz
freq at passband=2000.000000 Hz
freq at stopband=3000.000000 Hz
db at passband=-1.000000
db at stopband=-15.000000
Doing Bilinear method
T=0.000050
n=5.304446, rounded to 6
Gamma C =15324.588619
POLES Of H(s)

poles_of_hs =

-3966.29539304108 +      14802.4159246557i
-10836.1205316146 +      10836.1205316146i
-14802.4159246557 +      3966.29539304109i
-14802.4159246557 -      3966.29539304109i
-10836.1205316146 -      10836.1205316146i
-3966.29539304108 -      14802.4159246557i

k=1.295188e+025
H(s)=
Warning: Transfer function has complex coefficients.
> In tf.tf at 246
In nma_filter at 90

Transfer function:

(1.295e025-1.288e010i)
--------------------------------------------------------------------------
s^6 + (5.921e004-7.276e-012i) s^5 + (1.753e009-3.576e-007i) s^4

+ (3.29e013-0.01172i) s^3 + (4.117e017-224i) s^2 + (3.265e021

-2.49e006i) s + (1.295e025-1.288e010i)

POLES Of H(z)

hzp =

0.498385402712059 +     0.299971817994199i
0.49838540271206 -     0.299971817994199i
0.467705977269569 -    0.0939883945094121i
0.605559876472429 +     0.553064536118524i
0.467705977269569 +    0.0939883945094127i
0.60555987647243 -     0.553064536118524i

H(z)=
Warning: Transfer function has complex coefficients.
> In tf.tf at 246
In nma_filter at 107

Transfer function:

-(0.5053-1.947i) s^4 + (0.02086-1.348i) s^3 - (0.9771+0.04836i) s^2

+ (0.01411+0.02526i) s - (0.3816-0.3931i)

--------------------------------------------------------------------------
(0.1988-1.53i) s^4 + (0.6374+0.9373i) s^3 - (1.064+0.2778i) s^2

- (0.7076-0.6148i) s + (0.4667-0.6282i)

=================================================================

Sampling frequency=10000.000000 Hz
freq at passband=1000.000000 Hz
freq at stopband=2000.000000 Hz
db at passband=-3.000000
db at stopband=-10.000000
Doing impulse invariance method
T=1.000000
n=1.588388, rounded to 2
Gamma C =0.629065
POLES Of H(s)

poles_of_hs =

-0.444816082052343 +     0.444816082052343i
-0.444816082052343 -     0.444816082052343i

k=3.957227e-001
H(s)=
Warning: Transfer function has complex coefficients.
> In tf.tf at 246
In nma_filter at 92

Transfer function:
(0.3957-8.327e-017i)
---------------------------------------------------
s^2 + (0.8896-5.551e-017i) s + (0.3957-8.327e-017i)

POLES Of H(z)

hzp =

0.578571949761589 +     0.275792192443573i
0.578571949761589 -     0.275792192443573i

H(z)=
Warning: Transfer function has complex coefficients.
> In tf.tf at 246
In nma_filter at 107

Transfer function:
(-6.9e-016+1.253i)

======================================

Sampling frequency=10000.000000 Hz
freq at passband=1000.000000 Hz
freq at stopband=2000.000000 Hz
db at passband=-3.000000
db at stopband=-10.000000
Doing Bilinear method
T=0.000100
n=1.368163, rounded to 2
Gamma C =8389.390482
POLES Of H(s)

poles_of_hs =

-5932.19490014964 +      5932.19490014964i
-5932.19490014964 -      5932.19490014964i

k=7.038187e+007
H(s)=
Warning: Transfer function has complex coefficients.
> In tf.tf at 246
In nma_filter at 90

Transfer function:
(7.038e007-2.235e-008i)
---------------------------------------------------------
s^2 + (1.186e004-9.095e-013i) s + (7.038e007-2.235e-008i)

POLES Of H(z)

hzp =

0.458140439181504 -     0.308891358305335i
0.458140439181504 +     0.308891358305335i

H(z)=
Warning: Transfer function has complex coefficients.
> In tf.tf at 246
In nma_filter at 107

Transfer function:
(2.058e-016-1.78i)

EDU>>


### 6 References

1. Alan Oppenheim, Ronald Schafer, Digital Signal Processing,  Prentice-Hall, inc. 1975, Chapter 5.
2. Mostafa Shiva, Electrical engineering department, California state university, Fullerton, Lecture notes, handout H.
3. John Proakis, Dimitris Manolakis, digital signal processing, 3rd edition, section 8.3.