**2. Preliminaries**

In this chapter, we frequently use notations in control systems. For readers who are not familiar to these, we here recall basic notations and facts of control systems used throughout the chapter. We also show MATLAB codes for better understanding.

On the other hand, suppose that we have a transfer function at first:

The first command defines the variable *z* of *Z*-transform with sampling time 1, and the second

Min-Max Design of FIR Digital Filters by Semidefinite Programming 195

*<sup>G</sup>*(*z*) = *<sup>z</sup>*<sup>2</sup> <sup>+</sup> <sup>2</sup>*<sup>z</sup>* <sup>+</sup> <sup>1</sup>

To convert this to state-space matrices *A*, *B*, *C*, and *D*, use the command ss as follows:

*z*<sup>2</sup> + 0.5*z* + 1

.

>> z=tf('z',1);

>> ss(Gz)

a =

b =

c =

d =

>> Gz=(z^2+2\*z+1)/(z^2+0.5\*z+1);

x1 x2

x1 -0.5 -1 x2 1 0

x1 x2 y1 1.5 0

Sampling time (seconds): 1

These outputs shows that the state-space matrices are given by

−0.5 −1 1 0 , *B* = 1 0 

*<sup>G</sup>*(*z*) = *<sup>C</sup>*(*zI* <sup>−</sup> *<sup>A</sup>*)−1*<sup>B</sup>* <sup>+</sup> *<sup>D</sup>* <sup>=</sup>:

with sampling time 1. Note that the state-space representation in **Example 1** is *minimal* in that the state-space model describes the same input/output behavior with the minimum number of states. Such a system

We then introduce a useful notation, called *packed notation* (Vidyasagar, 1988), describing the

, *C* =

1.5 0

 *A B C D* (*z*).

, *D* = 1,

*A* = 

is called *minimal realization* (Rugh, 1996).

transfer function from state-space matrices as

Discrete-time model.

u1 x1 1 x2 0

u1 y1 1

command defines the following transfer function:

Let us begin with a linear system G represented in the following *state-space equation* or *state-space representation* (Rugh, 1996):

$$\mathcal{G}: \begin{cases} \, \, \mathbf{x}[k+1] = A\mathbf{x}[k] + Bu[k], \\ \, \, \, y[k] = \mathbf{C}\mathbf{x}[k] + Du[k], \, \, k = 0, 1, 2, \dots \end{cases} \tag{1}$$

The nonnegative number *<sup>k</sup>* denotes the time index. The vector *<sup>x</sup>*[*k*] <sup>∈</sup> **<sup>R</sup>***<sup>n</sup>* is called the state vector, *u*[*k*] ∈ **R** is the input and *y*[*k*] ∈ **R** is the output of the system G. The matrices *A* ∈ **<sup>R</sup>***n*×*n*, *<sup>B</sup>* <sup>∈</sup> **<sup>R</sup>***n*×1, *<sup>C</sup>* <sup>∈</sup> **<sup>R</sup>**1×*n*, and *<sup>D</sup>* <sup>∈</sup> **<sup>R</sup>** are assumed to be static, that is, independent of the time index *k*. Then the *transfer function G*(*z*) of the system G is defined by

$$\mathcal{G}(z) := \mathcal{C}(zI - A)^{-1}B + D, \quad z \in \mathbb{C}.$$

The transfer function *G*(*z*) is a rational function of *z* of the form

$$G(z) = \frac{b\_n z^n + b\_{n-1} z^{n-1} + \dots + b\_1 z + b\_0}{z^n + a\_{n-1} z^{n-1} + \dots + a\_1 z + a\_0}.$$

Note that *<sup>G</sup>*(*z*) is the *<sup>Z</sup>*-transform of the impulse response {*g*[*k*]}<sup>∞</sup> *<sup>k</sup>*=<sup>0</sup> of the system G with the initial state *x*[0] = 0, that is,

$$G(z) = \sum\_{k=0}^{\infty} g[k]z^{-k} = D + \sum\_{k=1}^{\infty} CA^{k-1}Bz^{-k}.$$

To convert a state-space equation to its transfer function, one can use the above equations or the MATLAB command tf. On the other hand, to convert a transfer function to a state-space equation, one can use realization theory which provides a method to derive the state space matrices from a given transfer function (Rugh, 1996). An easy way to obtain the matrices is to use MATLAB or Scilab with the command ss.

**Example 1.** *We here show an example of MATLAB commands. First, we define state-space matrices:*

```
>A=[0,1;-1,-2]; B=[0;1]; C=[1,1]; D=0;
>G=ss(A,B,C,D,1);
```
*This defines a state-space (ss) representation of* G *with the state-space matrices*

$$A = \begin{bmatrix} 0 & 1 \\ -1 & -2 \end{bmatrix} \text{ , } B = \begin{bmatrix} 0 \\ 1 \end{bmatrix} \text{ , } \mathbf{C} = \begin{bmatrix} 1 \ 1 \end{bmatrix} \text{ , } D = 0.$$

*The last argument* 1 *of* ss *sets the sampling time to be 1.*

To obtain the transfer function *<sup>G</sup>*(*z*) = *<sup>C</sup>*(*zI* <sup>−</sup> *<sup>A</sup>*)−1*<sup>B</sup>* <sup>+</sup> *<sup>D</sup>*, we can use the command tf

```
>> tf(G)
Transfer function:
    z+1
z^2 + 2 z + 1
Sampling time (seconds): 1
```
On the other hand, suppose that we have a transfer function at first:

>> z=tf('z',1); >> Gz=(z^2+2\*z+1)/(z^2+0.5\*z+1);

2 Digital Signal Processing

Let us begin with a linear system G represented in the following *state-space equation* or

The nonnegative number *<sup>k</sup>* denotes the time index. The vector *<sup>x</sup>*[*k*] <sup>∈</sup> **<sup>R</sup>***<sup>n</sup>* is called the state vector, *u*[*k*] ∈ **R** is the input and *y*[*k*] ∈ **R** is the output of the system G. The matrices *A* ∈ **<sup>R</sup>***n*×*n*, *<sup>B</sup>* <sup>∈</sup> **<sup>R</sup>***n*×1, *<sup>C</sup>* <sup>∈</sup> **<sup>R</sup>**1×*n*, and *<sup>D</sup>* <sup>∈</sup> **<sup>R</sup>** are assumed to be static, that is, independent of the

*<sup>G</sup>*(*z*) :<sup>=</sup> *<sup>C</sup>*(*zI* <sup>−</sup> *<sup>A</sup>*)−1*<sup>B</sup>* <sup>+</sup> *<sup>D</sup>*, *<sup>z</sup>* <sup>∈</sup> **<sup>C</sup>**.

*<sup>G</sup>*(*z*) = *bnz<sup>n</sup>* <sup>+</sup> *bn*−1*zn*−<sup>1</sup> <sup>+</sup> ··· <sup>+</sup> *<sup>b</sup>*1*<sup>z</sup>* <sup>+</sup> *<sup>b</sup>*<sup>0</sup>

*g*[*k*]*z*−*<sup>k</sup>* = *D* +

To convert a state-space equation to its transfer function, one can use the above equations or the MATLAB command tf. On the other hand, to convert a transfer function to a state-space equation, one can use realization theory which provides a method to derive the state space matrices from a given transfer function (Rugh, 1996). An easy way to obtain the matrices is to

**Example 1.** *We here show an example of MATLAB commands. First, we define state-space matrices:*

To obtain the transfer function *<sup>G</sup>*(*z*) = *<sup>C</sup>*(*zI* <sup>−</sup> *<sup>A</sup>*)−1*<sup>B</sup>* <sup>+</sup> *<sup>D</sup>*, we can use the command tf

*<sup>z</sup><sup>n</sup>* <sup>+</sup> *an*−1*zn*−<sup>1</sup> <sup>+</sup> ··· <sup>+</sup> *<sup>a</sup>*1*<sup>z</sup>* <sup>+</sup> *<sup>a</sup>*<sup>0</sup>

∞ ∑ *k*=1

, *C* =

1 1

, *D* = 0.

*CAk*−1*Bz*−*<sup>k</sup>*

*<sup>y</sup>*[*k*] = *Cx*[*k*] + *Du*[*k*], *<sup>k</sup>* <sup>=</sup> 0, 1, 2, . . . . (1)

.

.

*<sup>k</sup>*=<sup>0</sup> of the system G with the

*x*[*k* + 1] = *Ax*[*k*] + *Bu*[*k*],

time index *k*. Then the *transfer function G*(*z*) of the system G is defined by

The transfer function *G*(*z*) is a rational function of *z* of the form

Note that *<sup>G</sup>*(*z*) is the *<sup>Z</sup>*-transform of the impulse response {*g*[*k*]}<sup>∞</sup>

∞ ∑ *k*=0

*This defines a state-space (ss) representation of* G *with the state-space matrices*

 , *B* = 0 1 

 0 1 −1 −2

*G*(*z*) =

use MATLAB or Scilab with the command ss.

>A=[0,1;-1,-2]; B=[0;1]; C=[1,1]; D=0;

*A* =

*The last argument* 1 *of* ss *sets the sampling time to be 1.*

*state-space representation* (Rugh, 1996):

initial state *x*[0] = 0, that is,

>G=ss(A,B,C,D,1);

Transfer function: z+1 ------------ z^2 + 2 z + 1

Sampling time (seconds): 1

>> tf(G)

G :

The first command defines the variable *z* of *Z*-transform with sampling time 1, and the second command defines the following transfer function:

$$G(z) = \frac{z^2 + 2z + 1}{z^2 + 0.5z + 1}.$$

To convert this to state-space matrices *A*, *B*, *C*, and *D*, use the command ss as follows:

```
>> ss(Gz)
a =
        x1 x2
  x1 -0.5 -1
  x2 1 0
b =
      u1
  x1 1
  x2 0
c =
       x1 x2
  y1 1.5 0
d =
      u1
  y1 1
Sampling time (seconds): 1
Discrete-time model.
```
These outputs shows that the state-space matrices are given by

$$A = \begin{bmatrix} -0.5 \ -1 \\ 1 \ \ 0 \end{bmatrix}, B = \begin{bmatrix} 1 \\ 0 \end{bmatrix}, \text{ C} = \begin{bmatrix} 1.5 \ 0 \end{bmatrix}, \text{ D} = 1.$$

with sampling time 1.

Note that the state-space representation in **Example 1** is *minimal* in that the state-space model describes the same input/output behavior with the minimum number of states. Such a system is called *minimal realization* (Rugh, 1996).

We then introduce a useful notation, called *packed notation* (Vidyasagar, 1988), describing the transfer function from state-space matrices as

$$G(z) = \mathbb{C}(zI - A)^{-1}B + D =: \left[\frac{A \mid B}{\mathbb{C} \mid D}\right](z).$$

*H*<sup>∞</sup> optimization is thus minimization of the maximum value of a transfer function. This leads to robustness against uncertainty in the frequency domain. Moreover, it is known that the *H*<sup>∞</sup>

Min-Max Design of FIR Digital Filters by Semidefinite Programming 197

*<sup>u</sup>*<sup>∈</sup><sup>2</sup> *u* 0


G*u*<sup>2</sup> *u*<sup>2</sup> ,

norm of a transfer function *<sup>G</sup>*(*z*) is equivalent to the <sup>2</sup>*-induced norm* of <sup>G</sup>, that is,

*u*<sup>2</sup> :=

where *<sup>u</sup>*<sup>2</sup> is the <sup>2</sup> norm of *<sup>u</sup>*:

**3.** *H*<sup>∞</sup> **Design problems of FIR digital filters**

using the *H*<sup>∞</sup> norm defined in the previous section.

(Yamamoto et al., 2003). The problem is formulated as follows:

(*P* − *Q*)*W*<sup>∞</sup> = max

The procedure to solve this problem is shown in Section 4.

**3.1 FIR approximation of IIR filters**

*response) filter Q*(*z*) *which minimizes*

*where W is a given stable weighting function.*

signals.

*G*<sup>∞</sup> = G := sup

 ∞ ∑ *n*=0

The *H*<sup>∞</sup> optimization is minimization of the system gain when the worst case input is applied. This fact implies that the *H*<sup>∞</sup> optimization leads to robustness against uncertainty in input

In this section, we consider two fundamental problems in signal processing: filter approximation and inverse filtering. The problems are formulated as *H*<sup>∞</sup> optimization by

The first problem we consider is *approximation*. In signal processing, there are a number of design methods for IIR (infinite impulse response) filters, e.g., Butterworth, Chebyshev, Elliptic, and so on (Oppenheim & Schafer, 2009). In general, to achieve a given characteristic, IIR filters require fewer memory elements, i.e., *z*−1, than FIR (finite impulse response) filters. However, IIR filters may have a problem of instability since they have feedbacks in their circuits, and hence, we prefer an FIR filter to an IIR one in implementation. For this reason, we employ FIR approximation of a given IIR filter. This problem has been widely studied (Oppenheim & Schafer, 2009). Many of them are formulated by *H*<sup>2</sup> optimization; they aim at minimizing the average error between a given IIR filter and the FIR filter to be designed. This optimal filter works well *averagely*, but in the worst case, the filter may lead a large error. To guarantee the worst case performance, *H*<sup>∞</sup> optimization is applied to this problem

**Problem 1** (FIR approximation of IIR filters)**.** *Given an IIR filter P*(*z*)*, find an FIR (finite impulse*

*<sup>P</sup>*(*ejω*) <sup>−</sup> *<sup>Q</sup>*(*ejω*)

*W*(*ejω*) ,

 

*ω*∈[0,*π*]

Fig. 1. The *<sup>H</sup>*<sup>∞</sup> norm *<sup>G</sup>*<sup>∞</sup> is the maximum gain of the frequency response *<sup>G</sup>*(*ejω*).

By the packed notation, the following formulae are often used in this chapter:

$$
\begin{aligned}
\left[\frac{A\_1 \mid B\_1}{C\_1 \mid D\_1}\right] \times \left[\frac{A\_2 \mid B\_2}{C\_2 \mid D\_2}\right] &= \left[\frac{A\_2 \mid 0}{D\_1 C\_2 \mid A\_1} \frac{B\_2}{B\_1 D\_2}\right] \\
\left[\frac{A\_1 \mid B\_1}{C\_1 \mid D\_1}\right] \pm \left[\frac{A\_2 \mid B\_2}{C\_2 \mid D\_2}\right] &= \left[\frac{A\_1 \mid 0}{0 \mid A\_2} \frac{B\_1}{\pm B\_2}\right] \\
\left[\frac{A\_1 \mid B\_1}{C\_1 \mid C\_2}\right] \pm \left[\frac{B\_2}{C\_1 \mid C\_2}\right] &= \left[\frac{A\_1 \mid B\_1}{C\_1 \mid C\_2}\right] \\
\end{aligned}
$$

Next, we define *stability* of linear systems. The state-space system G in (1) is said to be *stable* if the eigenvalues *λ*1,..., *λ<sup>n</sup>* of the matrix *A* lie in the open unit circle **D** = {*z* ∈ **C** : |*z*| < 1}. Assume that the transfer function *G*(*z*) is irreducible. Then G is stable if and only if the poles of the transfer function *G*(*z*) lie in **D**. To compute the eigenvalues of *A* in MATLAB, use the command eig(A), and for the poles of *G*(*z*) use pole(Gz). The *H*<sup>∞</sup> *norm* of a stable transfer function *G*(*z*) is defined by

$$\|G\|\_{\infty} := \max\_{\omega \in [0,\pi]} \left| G(e^{j\omega}) \right|.$$

This is the maximum gain of the frequency response *<sup>G</sup>*(*ejω*) of <sup>G</sup> as shown in Fig. 1. The MATLAB code to compute the *H*<sup>∞</sup> norm of a transfer function is given as follows:

```
>> z=tf('z',1);
>> Gz=(z-1)/(z^2-0.5*z);
>> norm(Gz,inf)
ans =
    1.3333
```
This result shows that for the stable transfer function

$$G(z) = \frac{z - 1}{z^2 - 0.5z},$$

the *<sup>H</sup>*<sup>∞</sup> norm is given by *<sup>G</sup>*<sup>∞</sup> <sup>≈</sup> 1.3333.

*H*<sup>∞</sup> optimization is thus minimization of the maximum value of a transfer function. This leads to robustness against uncertainty in the frequency domain. Moreover, it is known that the *H*<sup>∞</sup> norm of a transfer function *<sup>G</sup>*(*z*) is equivalent to the <sup>2</sup>*-induced norm* of <sup>G</sup>, that is,

$$\|\|G\|\|\_{\infty} = \|\|\mathcal{G}\|\| := \sup\_{\substack{\boldsymbol{\mu} \in \ell^2 \\ \boldsymbol{\mu} \neq \boldsymbol{0}}} \frac{\|\|\mathcal{G}\boldsymbol{\mu}\|\_2}{\|\boldsymbol{\mu}\|\_2} \boldsymbol{\mu}$$

where *<sup>u</sup>*<sup>2</sup> is the <sup>2</sup> norm of *<sup>u</sup>*:

4 Digital Signal Processing

0 *π*

 = 

 = 

Next, we define *stability* of linear systems. The state-space system G in (1) is said to be *stable* if the eigenvalues *λ*1,..., *λ<sup>n</sup>* of the matrix *A* lie in the open unit circle **D** = {*z* ∈ **C** : |*z*| < 1}. Assume that the transfer function *G*(*z*) is irreducible. Then G is stable if and only if the poles of the transfer function *G*(*z*) lie in **D**. To compute the eigenvalues of *A* in MATLAB, use the

*ω*∈[0,*π*]

This is the maximum gain of the frequency response *<sup>G</sup>*(*ejω*) of <sup>G</sup> as shown in Fig. 1. The

*<sup>G</sup>*(*z*) = *<sup>z</sup>* <sup>−</sup> <sup>1</sup>

*<sup>z</sup>*<sup>2</sup> − 0.5*<sup>z</sup>*

,

 *G*(*ejω*) .

Fig. 1. The *<sup>H</sup>*<sup>∞</sup> norm *<sup>G</sup>*<sup>∞</sup> is the maximum gain of the frequency response *<sup>G</sup>*(*ejω*).

By the packed notation, the following formulae are often used in this chapter:

 *A*<sup>2</sup> *B*<sup>2</sup> *C*<sup>2</sup> *D*<sup>2</sup>

 *A*<sup>2</sup> *B*<sup>2</sup> *C*<sup>2</sup> *D*<sup>2</sup>

*G*<sup>∞</sup> := max

MATLAB code to compute the *H*<sup>∞</sup> norm of a transfer function is given as follows:

*ω*

 ,

 .

*A*<sup>2</sup> 0 *B*<sup>2</sup> *B*1*C*<sup>2</sup> *A*<sup>1</sup> *B*1*D*<sup>2</sup> *D*1*C*<sup>2</sup> *C*<sup>1</sup> *D*1*D*<sup>2</sup>

*A*<sup>1</sup> 0 *B*<sup>1</sup> 0 *A*<sup>2</sup> ±*B*<sup>2</sup> *C*<sup>1</sup> *C*<sup>2</sup> *D*<sup>1</sup> ± *D*<sup>2</sup>

*G*<sup>∞</sup>

 *G*(*ejω*) 

 *A*<sup>1</sup> *B*<sup>1</sup> *C*<sup>1</sup> *D*<sup>1</sup>

 *A*<sup>1</sup> *B*<sup>1</sup> *C*<sup>1</sup> *D*<sup>1</sup>  ×

 ±

command eig(A), and for the poles of *G*(*z*) use pole(Gz). The *H*<sup>∞</sup> *norm* of a stable transfer function *G*(*z*) is defined by

This result shows that for the stable transfer function

the *<sup>H</sup>*<sup>∞</sup> norm is given by *<sup>G</sup>*<sup>∞</sup> <sup>≈</sup> 1.3333.

>> z=tf('z',1);

>> norm(Gz,inf)

1.3333

ans =

>> Gz=(z-1)/(z^2-0.5\*z);

$$\|\|u\|\|\_{2} := \left(\sum\_{n=0}^{\infty} \left|\mu[k]\right|^{2}\right)^{1/2}.$$

The *H*<sup>∞</sup> optimization is minimization of the system gain when the worst case input is applied. This fact implies that the *H*<sup>∞</sup> optimization leads to robustness against uncertainty in input signals.
