Design of IIR Digital Filters with Arbitrary Flatness Using Iterative Quadratic Programming

This paper presents a design method of Chebyshev-type and inverse-Chebyshev-type infinite impulse response (IIR) filters with an approximately linear phase response. In the design of Chebyshev-type filters, the flatness condition in the stopband is preincorporated into a transfer function, and an equiripple characteristic in the passband is achieved by iteratively solving the QP problem using the transfer function. In the design of inverse-Chebyshev-type filters, the flatness condition in the passband is added to the constraint of the QP problem as the linear matrix equality, and an equiripple characteristic in the stopband is realized by iteratively solving the QP problem. To guarantee the stability of the obtained filters, we apply the extended positive realness to the QP problem. As a result, the proposed method can design the filters with more high precision than the conventional methods. The effectiveness of the proposed design method is illustrated with some examples.


Introduction
Digital filters are needed in many applications in signal processing.Infinite impulse response (IIR) filters are useful for high-speed processing, and IIR filters with lower order can be realized that are comparable to the magnitude responses of finite impulse response (FIR) filters.However, its design is more difficult than FIR filters because IIR filters have a rational transfer function and are not always stable.For that reason, the study on IIR filter design is still a hot topic in the area of digital signal processing and many design methods on the IIR filters have been presented [1][2][3][4][5][6][7][8].In this paper, we will treat two types of filters: Chebyshev type which is equiripple in the passband and flat in the stopband, Inverse Chebyshev type which is equiripple in the stopband and flat in the stopband.These filters are needed in many image processing applications to suppress ringing and chessboard distortion and are used in each situation appropriately.Moreover, in many signal processing applications, such as wavelet transform and channel equalization, the complex coefficient digital filters are required to meet some specifications that cannot be achieved by real coefficient filters, such as asymmetric spectral response [9][10][11].
In [12,13], the design methods based on the Remez algorithm were proposed for the Chebyshev-type and inverse-Chebyshev-type IIR filters with an approximately linear phase response.These methods can design the filters with small computational complexity.However, the filters that can be designed using these methods are restricted greatly because of a condition imposed on setting the initial value.Moreover, these methods cannot guarantee the stability of the filter obtained.By using the linear semi-infinite programming and the extended positive realness, the design method of stable inverse Chebyshev-type IIR filters with an approximately linear phase response has been proposed [14].This method needs to set an initial value appropriately and the performance of the filter obtained is dependent on the initial value.
In this paper, a design method based on quadratic programming (QP) is presented for stable IIR filters and for FIR filters with prescribed flatness in the passband or stopband and an approximately linear phase response in the passband.To guarantee the stability of the filter obtained, we apply the extended positive realness to the QP problem.As a result, the proposed method can design the filters that cannot be designed by the conventional methods.This paper is organized as follows: in Section 2, the flatness conditions in the passband and stopband are described.In Section 3, the design algorithm for inverse-Chebyshev-type filters is described.In Section 4, the design algorithm for Chebyshevtype filters is described.To verify the effectiveness of the proposed method, several design examples are given in Section 5. Section 6 is the conclusions of this work.

IIR Digital Filters and Flatness Conditions
The frequency response (  ) of an IIR digital filter is defined as where  and  are the orders of the numerator and denominator, respectively.  =   +   and   =   +   are the filter coefficients and  0 = 1 in general.The desired frequency response   (  ) can be expressed as where   is a desired group delay.Then, the flatness conditions of the magnitude and group delay at  =   in the passband are given as follows: where   is a parameter expressing the flatness in the stopband.

Update of the Weighting Function 𝑊(𝜔).
It has been well known that the filters obtained under weighted least square criterion have a large ripple near the band edges.So in order to realize the equiripple characteristics in the passband or stopband or both, the weighting function used at every iteration is adjusted using the modified Lawson's method [16] and the QP problem is solved to obtain the coefficients.In this paper, the weighting function () in th iteration step is updated as follows: where the envelope function   () is given as the function of straight line formed by joining together all the extremal points of the same frequency band of interest on the error function which is expressed as Using the extremal points ω of   (),   () can be calculated by where ω denotes the th extremal frequency of the error function   ().

Stability Constraint.
To obtain the stable IIR filters, the stability condition based on a positive realness has been applied to many design methods.However, the use of this condition may exclude the candidate for the transfer function with excellent performance because this condition is a sufficient condition to assure the stability and is often too restrictive.In [17], an iterative method for the stability guarantee based on the positive realness was proposed in order to get a better transfer function.In this method, a stability condition is given by where  < 1.If the maximum pole radius  max of the filter obtained using a given  is greater than prescribed maximum allowable pole radius   ,  is updated and redesign is carried out using the updated .The update of  and the redesign are repeated until satisfying | max −   | ≤   , where   is a positive small value.The update procedure of  is described in Section 4.2.
Using the discrete angular frequency   ( = 1, . . ., ), (18) can be expressed as the linear matrix inequality where Thus, the design problem in which the stability constraint and the update of the weighting function were considered becomes a standard QP problem as below: sub. to Γh  ≥ , where This problem can be solved using a powerful QP tool, such as  in MATLAB.

Design of Chebyshev-Type Filters
In order to meet the flatness condition of (22), it is necessary to place   multiple zeros at  =   .Hence, the frequency response Ĥ(  ) can be expressed as where   =   +  .To obtain the equiripple responses in the passband, we consider the following iterative design formula: (24) Thus, the design problem in which the stability constraint and the update of the weighting function were considered becomes a standard QP problem as below: where ĥ = [ 0 , . .

Design Procedure.
The design procedure of the proposed method is summarized as follows.
Step 0. Set the design specifications.
Step 1. Solve the QP problem to obtain the filter coefficient h  .
Step 4. Calculate   = min all  Re{()} from the obtained filter and then set to   =   and   = 1.
Step 6. Solve the QP problem using the updated .
Step 9.If  max <   , set to   = ; if  max >   , set to   = , and then go back to Step 5.

Design Examples
In this section, some design examples are given to illustrate the effectiveness of the proposed method.In all the following examples,  = 10 −7 ,   = 10 −5 , and the initial value of  is −10 2 .The program of the proposed design algorithm was coded by MATLAB and was run on the PC which has Core i5-M520 processor and 2 GB memory.

Example 1.
To compare with the conventional method [12], we design the Chebyshev-type IIR filters which have an equiripple characteristic in passband and a flat characteristic in stopband.The filter specifications are as follows: In the proposed method, we set the maximum allowable pole radius   to 0.98, and the number of grid points is  = 500 and  = 1000.Note that the conventional method [12] cannot guarantee the stability of the obtained filter.The performance of the obtained filters is shown in Table 1 which are obtained by the conventional method [12] based on Remez algorithm.In Table 1,   is the maximum magnitude error in passband,  err is the maximum group delay error in passband, and  max is the maximum pole radius of the obtained filter.From Table 1, it is confirmed that the proposed method needs more computational cost (CPU time and iteration number) than the conventional method, but the characteristic of the obtained filters by the proposed method is about the same as or better than them by the conventional method.Particularly, the filter obtained by the conventional method is unstable when   = 9.The stability of the filters obtained by the method which cannot guarantee the stability depends on a given group delay specification.In fact, when the group delay   is less than 9.5, all filters obtained by conventional method were unstable.In contrast, the proposed method can always obtain the stable filters because the maximum pole radius of the filter obtained can be prespecified.Therefore, the proposed method can design the filters that cannot be designed by the conventional method.

Example 2.
Next we design the inverse Chebyshev-type IIR filters with the real coefficients in order to compare with the conventional method [13].The filter specifications are as follows: In the proposed method, we set the maximum allowable pole radius   to 0.98, and the number of grid points is  = 500 and  = 1000.Note that the conventional method [13] cannot guarantee the stability of the obtained filter.
The numerical performance of the obtained filters is shown in Table 2 which are obtained by the conventional method [13] based on Remez algorithm.In Table 2,   is the minimum stopband attenuation in dB.From Table 2, it is confirmed that the proposed method needs more computational cost than the conventional method but gives slightly better characteristics.Moreover, when the group delay   is less than 7.2, the filters obtained by conventional method were unstable.In contrast, the proposed method can always obtain the stable filters.Therefore, the proposed method can design the filters that cannot be designed by the conventional method.
Next, we designed the filter with   = 9 for comparison with [14].The minimum stopband attenuation of the filter obtained was 41.38 dB in the proposed method and 40.04 dB in [14].It is thought that this difference is because the performance of the filter obtained by the conventional method is dependent on the given initial value.
In this example, we use   = 1.0,  = 500, and  = 1000.The proposed method took an average of 0.78 seconds for convergence.The magnitude response and group delay response of the obtained filter are shown in Figures 1, 2, and  3. From these figures, it is confirmed that the magnitude and group delay responses have a flat characteristic at  = 0.1 and the magnitude response in stopband is equiripple.Table 3 is the numerical performance of the obtained filter.In this example, we set the number of the grit points  to 10,000.The proposed method took an average of 2,794 seconds for convergence.The magnitude response and group delay response of the obtained filter are shown in Figures 4, 5, and 6.From these figures, it is confirmed that the magnitude and group delay responses have a flat characteristic at  = 0 and the magnitude response in stopband is equiripple.The minimum stopband attenuation of the obtained filters was 72.98 (dB) for   = 10, 59.11 (dB) for   = 12, and 47.19 (dB) for   = 14.We confirmed that these filters cannot be designed by the conventional method [13].Hence, the proposed method can also design the high-order filters that cannot be designed by the conventional method.

Conclusion
In this paper, a design method based on quadratic programming has been proposed for approximately linear phase IIR filters and FIR filters with prescribed flatness in passband or stopband.To guarantee the stability of the filter obtained, we applied the extended positive realness to the QP problem.With this method, stable IIR filters and FIR filters that are Chebyshev type and inverse Chebyshev type can be easily designed.In the design examples, we showed that the performance of the filter obtained by the proposed method is about the same as or better than the conventional method, and the filters which cannot be designed by the conventional method can also be designed.

Figure 2 :
Figure 2: Passband magnitude response of the proposed IIR filters in Example 3.

Figure 3 :
Figure 3: Group delay response of the proposed IIR filters in Example 3.

Figure 4 :
Figure 4: Overall magnitude response of the proposed FIR filters in Example 4.

Figure 5 :
Figure 5: Passband magnitude response of the proposed FIR filters in Example 4.

Figure 6 :
Figure 6: Group delay response of the proposed FIR filters in Example 4.

Table 3 :
The performance of the filters in Example 3.