Signal Processing 91 (2011) 2589–2594
Contents lists available at ScienceDirect
Signal Processing journal homepage: www.elsevier.com/locate/sigpro
Method for FIR filters design with compressed cosine using Chebyshev’s norm Peter Apostolov Institute for Special Technical Equipment – MI, POB 83, Sofia 1799, Bulgaria
a r t i c l e i n f o
abstract
Article history: Received 8 April 2011 Accepted 22 May 2011 Available online 27 May 2011
This paper considers a new method for FIR filters design. The method uses an LN optimality norm. To achieve a better approximating effect, a new modulating function which compresses the oscillations of the cosine is proposed. A parameter sets the gradient of the modulating function, with respect to the oscillations’ compression. The approximating polynomial is carried out using Remez’ exchange algorithm. An optimal polynomial with lowest possible (four) degree, that approximates an ideal filter’s response with high precision is proposed. With the proposed method a FIR filter with arbitrary specifications can be designed. Design examples of FIR filters with a minimization of calculation are performed. The obtained filter’s responses are close to the ideal response. The design examples demonstrate that the proposed approach may be a good alternative in several applications. & 2011 Elsevier B.V. All rights reserved.
Keywords: FIR digital filters Frequency response Polynomial approximation Chebyshev norm
polynomial defines the error function
1. Introduction The filters’ design is a mathematical problem for an approximation of the ideal filter’s response (Fig. 1). This is an ideal function, comprising rectangular shape with two areas: pass band and stop band. The analytical expression of the ideal low pass filter response is ( 1, o 2 ½0, oc DðoÞ ¼ , ð1Þ 0, o 2 ðoc , p where oc is the cut-off frequency. In the FIR filter design, the approximation function is a cosine polynomial m m X X 1 AðoÞ ¼ pffiffiffi a0 þ an cosðnoÞ9 an cosðnoÞ: 2 n¼1 n¼0
ð2Þ
The coefficients of the filter’s impulse response are determined by the coefficients of the polynomial. The difference between the ideal function and the approximating E-mail address:
[email protected] 0165-1684/$ - see front matter & 2011 Elsevier B.V. All rights reserved. doi:10.1016/j.sigpro.2011.05.014
EðoÞ ¼ DðoÞAðoÞ,
o 2 ½0, p:
ð3Þ
When the approximation is optimal, the graph of the approximating polynomial oscillates with equal amplitude close to the ideal function’s graph. The amplitude of the oscillations determines the approximation’s error. The three most popular norms used in FIR design are as follows: 1. Weighted L1 norm. The approximation’s error is Z p WðoÞEðoÞdo: :EðoÞ:1 ¼
ð4Þ
0
2. L2 error-weighted integral least-squares norm is given by Z p 2 2 :EðoÞ:2 ¼ WðoÞEðoÞ do: ð5Þ 0
3. LN error-weighted Chebyshev’s norm :EðoÞ:1 ¼ max WðoÞEðoÞ : o2½o, p
ð6Þ
2590
P. Apostolov / Signal Processing 91 (2011) 2589–2594
Fig. 1. Ideal low pass filter’s response.
Fig. 2. Approximations with 32-degree polynomials using L1, L2 and LN norms.
In the all above cases W(o) is a positive weighting function, used in order to weight certain frequencies. When W(o)¼1, the maximal error in the pass band and stop band is equal. In this case the approximation using LN norm is equiripple. Fig. 2 shows a comparison of 32-degree approximating polynomials using L1, L2 and LN norms. It is seen that in the FIR filters design a suitable trade-off between flatness and the transition bandwidth must be made. In all the criteria, the functions have the oscillations in the pass band and stop band. These oscillations are undesirable. The goal in the design is to obtain a rectangular shape of the ideal filter’s response, that is maximally flat pass band and stop band, and narrowest possible transition band. In L1 and L2 cases the oscillations increase near to the function’s transition band. This is due to the Gibbs’ phenomenon [1]. The filter design using L2 criterion is studied in details in many publications [2–8]. A various methods for
reduction of the Gibbs’ phenomenon have been proposed: with window functions, ‘‘don’t care’’ transition band, change the transition band, etc. Methods for L2 FIR filters design, where the function’s overshot near the transition band is constrained, are presented in Refs. [9–12]. The design is known as ‘‘Constrained Least Square’’ filters. The synthesis is performed with exchange algorithm, having set the maximum value of the overshot. As a result an alignment of the oscillations’ amplitude near the transition band is obtained. This leads to the transition’s band expansion. In Ref. [10] the transition band is either not set as an input parameter in the synthesis, or set only at one of the frequencies that define it. The transition bandwidth is the result of the approximation, and is called ‘‘induced’’ transition band. In Ref. [13] computationally efficient method for L2 FIR filter design is proposed. There are methods for L2 FIR filters design in which the transition band is defined. In Ref. [14] the proposed algorithm is not convergent for arbitrary filter’s specification. In Ref. [15] L2 filters can be designed with arbitrary specification with ‘‘nearly linear phase’’. With approximation using L1 norm [16] filters with flatter response, but with wider transition band are obtained (Fig. 2). The LN norm is the most suitable for the filter design problem. The design is performed with the well known Parks–McClellan’s method [17,18]. This is minimax approximation using the Chebyshev’s norm. With this method, the transition band can be accurately defined and arbitrarily narrow. However this increases the approximation’s error and the amplitude of oscillations (Fig. 2). An important advantage of the method is that the synthesis is performed with Remez’ exchange algorithm that has fast convergence and low computational complexity. It was found that with the same specification (pass band ripple, stop band attenuation and transition band width) with Parks–McClellan’s method, an approximation with polynomial of lowest degree has been done. That means that the LN filters will be implemented with the least number of coefficients. A new polynomial approximation using LN norm producing FIR filters with response close to the ideal will be proposed in the article.
2. Background As noted in all criteria, the approximating polynomial is a sum of cosines. The LN method’s theoretical base is the Alternation theorem [19].
If a function f(o) is defined and continuous in a closed
definition area, it can be approximated by trigonometric polynomial Am(o) by power m, with the basis function cosð:Þ. The polynomial is the unique and best approximation, if the error function E(o) has at least m þ2 extremes in the definition domain. All extremes are alternative and their modules are equal to a positive number e.
P. Apostolov / Signal Processing 91 (2011) 2589–2594
Such approximation can be done if the polynomial is a linear combination of Chebyshev’s polynomials (Fig. 3) AðoÞ ¼
m X n¼0
o , an cos narccos
p
o 2 ½0, p:
ð7Þ
The approximation is optimal, but inefficient, because the function’s graph compresses its oscillations on both ends of the definition domain, while in the transition band they are most sparse. Thus, a high slope of the function in the transition band cannot be obtained. From Eq. (7) it is obvious that the oscillation’s compression is owing to the arccosð:Þ function, which modulates the cosð:Þ function’s argument, as is shown in Fig. 4. It is logical to assume if another modulating function with inverse slope to arccosð:Þ, and higher gradient of the graph is applied, a polynomial with higher oscillations’ compression in the transition band, and greater slope will be obtained. A new function called ‘‘Third basis function’’
in Ref. [20] is proposed
p b yðxÞ ¼ cos m tanh bx þ1 ; 2 2
2591
x 2 ½0,1:
ð8Þ
The modulating function is tanhð:Þ, which comprises parameter b 44. Changing the parameter sets the slope of the modulating function, and the oscillations’ compression, as is shown in Fig. 5. Similar to the Parks–McClellan’s method, the approximating polynomial with Remez’ exchange algorithm [21] is performed. Concerning the specific requirements of the algorithm, the degree of the polynomial m is an even number. For the purposes of the FIR filters design the polynomial has the form
m þ1 X p 2bðooc Þ Am ðoÞ ¼ tanh ak cos ðk1Þ þ1 ; 2 osampl k¼1
o 2 ½0, osampl =2,
ð9Þ
where osampl is the sampling frequency. The goal in the proposed method is that the amplitude of the function’s overshot decreases when the parameter b increases, without changing the width of the transition band, or increasing the degree of the polynomial. In the other methods it is not so. Fig. 6 shows the graph of polynomial with the lowest possible (four) degree for two values of the parameter b. The pass band ripple DP and stop band attenuation DS are determined by the approximation’s error e as follows: DP ¼ 20 lgð1jejÞdB;
ð10Þ
DS ¼ 10 lgjejdB:
ð11Þ
3. Design examples Fig. 3. Approximation with linear combination of Chebyshev’s polynomials.
The properties of the filters are illustrated through two examples.
Fig. 4. Modulating functionarccosð:Þ and Chebyshev’s polynomial.
Fig. 5. Modulating functiontanhð:Þ and third basis function, m ¼10
2592
P. Apostolov / Signal Processing 91 (2011) 2589–2594
3.1. Low pass FIR filter Filter’s specification: cut-off frequency oc ¼0.2p rad/s; transition band width Doc ¼0.005p rad/s, sampling frequency osampl ¼2p rad/s; stop band attenuation DS Z 60 dB, power of the polynomial m¼4; weighting function W¼1. Parameter’s value b can be determined approximately
b¼
osampl DS 730: 32:9Doc
ð12Þ
The polynomial’s coefficients are determined with the Remez’ algorithm: a1 ¼0.5; a2 ¼0.5628; a3 ¼0; a4 ¼ 0.0628; a5 ¼0, and approximation’s error e ¼9.4108e 7. The exact values of the pass band ripple and stop band attenuation are calculated by Eqs. (10) and (11); DP¼ 8.1741e 6 dB; DS ¼ 60.2637 dB. Eq. (9) is the filter’s magnitude response. It should be noted that the coefficients of the impulse response are not obtained directly from polynomial’s coefficients as by the other methods, as the argument of the function cosð:Þ contains other, non-linear function. The impulse response can be determined with IFFT of the 2N samples of the magnitude response. The filter’s design is done using a method, known as ‘‘frequency sampling filter’’ [22,23]. It is based on FFT in 2N samples, where N is an integer positive number. For the purposes of the design a window function with 2N 1 values of the magnitude response (9) are calculated. The signal filtration is performed with convolution between FFT of the input signal and the window function. The filter’s structure is shown in Fig. 7. 1+ 2 1+ 1 1- 11 1- 2
2 1 - 1 - 2
0 0
pass / 0.5 stop /
1
Fig. 6. Approximation’s error e depending of the parameter b by fixed transition band.
IN
ADC
Buffer 2N
It should be noted that a larger number of samples should be determined to realize a narrow transition band. In this particular case the appropriate value is N ¼14. Therefore the window function consists of 8192 values, which determine the filter’s length. Fig. 8 shows the filter’s magnitude response. It is obvious that the characteristic is very close to the rectangular shape of the ideal low pass response. The magnitude response is maximally flat. It is seen that large parts of pass band and stop band are constant: 1 e and e. Since e is very small number (approximately 1e 6), it can be assumed that most of the values of the magnitude response in the pass band are equal to one, and in the stop band they are zeros. This circumstance leads to significant reduction of the calculation process of the filter. For the performance of the signal’s filtration, it is necessary to execute convolution only in the frequency band slightly wider than the transition band. Of course, it must include the values of the actual transition band. The rest of the pass band frequencies are transferred directly to output buffer (since it is not necessary to multiply by one), and those of the stop band are equal to zero. Fig. 9 shows the magnitude response of the same filter with described calculations’ minimization. Thus the filter of the considered example is realized with only 64 multiplications. 3.2. Band pass FIR filter Filter’s specification: Low cut-off frequency oc1 ¼ 0.15p rad/s; high cut-off frequency oc2 ¼0.25p rad/s; transition band width for both cut-off frequencies Doc ¼0.005p rad/s, sampling frequency osampl ¼2p rad/s; stop band attenuation for both stop bands DS Z60 dB, weighting function W¼1. A band pass filter with direct approximation of the ideal function of a band pass filter with arbitrary specification cannot be designed with the considered method. This is due to the nature of the Third basis function, which compresses its oscillations with non-linear dependence. Filters with 10–15 dB attenuation in the stop band can be designed with direct approximation. The Remez’ algorithm loses convergence for high values of the stop band attenuation. It is appropriate that the band pass filter be designed with window function, which is a combination of low pass and high pass filters of fourth order of the polynomial, with cut-off frequencies equal to those of the band pass filter’s specification. Fig. 10 shows the magnitude response of the band pas filter, designed in the described manner.
FFT 2N
Buffer 2N
Window Function 2N–1 Fig. 7. Filter’s structure.
IFFT 2N
OUT
P. Apostolov / Signal Processing 91 (2011) 2589–2594
2593
Fig. 8. Low pass FIR filter–magnitude response.
Fig. 9. Magnitude response with minimization of the calculations.
4. Conclusion FIR filters with responses close to ideal filter’s response can be designed with the proposed method. A basis function that compresses the oscillations of the function in the transition band is proposed. The method gets its name from this compression. The compression ratio is determined by the parameter b. The approximation is
carried out using LN norm that most closely approximates the transition band. The approximating polynomial is obtained by Remez’ algorithm, which has fast convergence and low computational complexity. In this case, the low degree of the polynomial (four) involves iterative solving of a system of six linear equations. In other methods they are much more. The method has very low computational cost and design complexity, which is an
2594
P. Apostolov / Signal Processing 91 (2011) 2589–2594
Fig. 10. Band pass FIR filter.
important advantage. With four degree polynomial FIR filter with arbitrary specifications in a fixed transition band can be designed. The approximation error e depends on the parameter b, and defines the pass band ripple and the stop band attenuation. No other polynomial of 4th or lower degree with better approximating properties is known till date. When the width of the transition band is equal to zero and parameter b ¼N, the polynomial’s graph coincides with the ideal response’s rectangular shape. Of course, this is impossible in reality. The theoretical design possibilities are limited by the computer calculation accuracy. The implementation of the proposed method makes sense in the design of filters with extreme, close to ideal response. The repeated reduction of the computational operations (Fig.7) is effective in approximations with values of e close to zero and a very narrow transition band. The proposed method may be a good alternative in several applications in the FIR filters’ design. References [1] B. Porat, A Course in Digital Signal Processing, New York, Wiley, 1997. [2] D.W. Tufts, J.T. Francis, Designing digital low-pass filters— comparison of some methods and criteria, IEEE Trans. Audio Electroacoust. AE-18 (1970) 487–494. [3] V.R. Algazi, M. Suk, On the frequency weighted least-squares design of finite duration filters, IEEE Trans. Circuits Syst. CAS-22 (12) (Dec. 1975) 943–953. [4] C.S. Burrus, A.W. Soewito, R.A. Gopinath, Least squared error FIR filter design with transition bands, IEEE Trans. Signal Process. 40 (6) (Jun. 1992) 1327–1340. [5] W.C. Kellogg, Time domain design of nonrecursive least meansquare digital filters, IEEE Trans. Audio Electroacoust. AE-20 (2) (Jun. 1972) 155–158. [6] Y.C. Lim, S.R. Parker, Discrete coefficient FIR digital filter design based upon an LMS criteria, IEEE Trans. Circuits Syst. CAS-30 (10) (Oct. 1983) 723–739. [7] M. Okuda, M. Ikehara, S. Takahashi, Fast and stable least-squaresapproach for the design of linear phase FIR filters, IEEE Trans. Signal Process. 46 (6) (Jun. 1998) 1485–1493.
[8] E.Z. Psarakis, G.V. Moustakides, An L2 based method for the design of one dimensional zero phase FIR filters, IEEE Trans. Circuits Syst. II: Analog Digital Signal Process. 44 (July 1997) 591–601. [9] I.W. Selesnick, M. Lang, C.S. Burrus, Constrained least square design of FIR filters without specified transition bands, IEEE Trans. Signal Process. SP-44 (Aug. 1996) 1879–1892. [10] I.W. Selesnick, M. Lang, C.S. Burrus, A modified algorithm for constrained least square design of multiband FIR filters without specified transition bands, IEEE Trans. Signal Process. SP-46 (Feb. 1998) 497–501. [11] Y.M. Law, C.W. Kok, Constrained eigenfilter design without specified transition bands. IEEE Trans. Circuits Syst. II 14-21. [12] J.W. Adams, J.L. Sullivan, Peak constrained least-squares optimization, IEEE Trans. Signal Process. SP-46 (Feb. 1998) 306–321. [13] P.P. Vaidyanathan, T.Q. Nguyen, Eigenfilters: a new approach to least squares FIR filter design and applications including Nyquist filters, IEEE Trans. Circuits Syst. CAS-34 (1) (Jan. 1987) 11–23. [14] J.W. Adams, FIR digital filters with least-squares stopbands subject to peak-gain constraints, IEEE Trans. Circuits Syst. Theory 38 (4) (Apr. 1991) 376–388. [15] E.Z. Psarakis, A weighted L2-based method for the design of arbitrary one-dimensional FIR digital filters, Signal Process. 86 (2006) 937–950. [16] Liron D. Grossmann, Yonina C. Eldar, An L1-method for the design of linear-phase FIR digital filters, Signal Process. 55 (11) (2007) 5253–5265. [17] T.W. Parks, J.H. McClellan, A program for the design of linear phase FIR digital filters, IEEE Trans. Audio Electroacoust. AU-20 (August 1972) 196–199. [18] J.H. McClellan, T.W. Parks, A unified approach to the design of optimum FIR linear-phase digital filter, IEEE Trans. Circuit Theory CT-20 (Nov. 1973) 697–701. [19] L.R. Rabiner, J.H. McClellan, T.W. Parks, FIR digital filter design techniques using weighted Chebyshev approximations, Proc. IEEE 63 (Apr. 1975) 595–610. [20] P.S. Apostolov, Mathematical approximations in the technical communications, DSc thesis, Sofia, 2010. [21] E.Ya. Remez, Fundamentals of numerical methods for Chebyshev approximations, Naukova Dumka, Kiev, 1969. [22] S.P. Harris, E.C. Ifeachor, Automatic design of frequency sampling filters by hybrid Genetic Algorithm Techniques, IEEE Transaction Signal Process. 46 (12) (1988) 3304–3314. [23] P.A. Lynn, Frequency sampling filters with integer multipliers, in: R.E. Bogner, A.G. Constantenides (Eds.), Introduction to Digital Filtering, Wiley, New York, 1975.