Filter Design: FIR Filters

Massachusetts Institute of Technology Department of Electrical Engineering and Computer Science 6.341: Discrete-Time Signal Processing OpenCourseWare 2006 Lecture 9 Filter Design: FIR Filters Reading: Sections 7.2 - 7.6 in Oppenheim, Schafer & Buck (OSB). Unlike discrete-time IIR filters which are generally obtained by transforming continuoustime IIR systems, FIR filters are almost always implemented in discrete time: a desired magni­ tude response is approximated directly in the DT domain. In addition, a linear phase constraint is often imposed during FIR filter design. In this lecture, we will look at two FIR filter design techniques: by windowing and through the Parks-McClellan Algorithm. FIR Filter Design by Windowing In designing FIR filter, given the frequency response Hd (ejω ) and impulse response hd [n] of an ideal system, we would like to approximate the infinitely long hd [n] with a finite sequence h[n], where h[n] = 0 except for 0 ≤ n ≤ M . Consider an ideal lowpass filter whose frequency response is finite and rectangular. A possible approximation error criterion can be defined as: � π 1 |Hd (ejω ) − H(ejω )|2 dω . �= 2π −π To minimize �, use Parseval’s theorem: � π 1 |Hd (ejω ) − H(ejω )|2 dω �= 2π −π = ∞ � |hd [n] − h[n]|2 = n=−∞ h[n] = |hd [n] − h[n]|2 + n=0 � ⇒ M � hd [n] 0 ≤ n ≤ M 0 otherwise � |hd [n] − h[n]|2 n=Z\[0,M ] . The optimal FIR approximation using the mean square error criterion gives the truncation of the ideal impulse response. We can represent h[n] as the product of the idal impulse response with a finite-duration rectangular window w[m]: � 1 0≤n≤M h[n] = hd [n]w[n], w[n] = 0 otherwise � π 1 ⇒ H(ejω ) = Hd (ejθ )W (ej(ω−θ) )dθ . 2π −π 1 OSB Figure 7.19(a) depicts this convolution process implied by truncation of the ideal impulse response. The resulting magnitude approximation as the sinc pulse W (ej(ω−θ) ) slides pass the ideal frequency response Hd (ejω ) is shown in Figure 7.19(b). When W (ej(ω−θ) ) moves across the discontinuity of Hd (ejω ), a transition band results and ripples occur on both sides. As shown in OSB Figure 7.23, the main lobe of the window frequency response controls transition bandwidth Δω ≈ 2π(2/(M + 1)). The mainlobe is defined as the region between the first zero crossings on either side of the origin. It is desirable to have W (ejω ) as concentrated in frequency as possible. On the other hand, sidelobes control passband and stopband ripples. The larger the area under the sidelobes, the larger the ripples. Passband and stopband ripples are approximately equal over a wide range of frequencies. For any given window length, we would like W (ejω ) to be “most like an impulse,” with a narrow main lobe and low sidelobes. For example, OSB Figure 7.20 shows a magnitude response of a rectangular window where M = 7. Although increasing the rectangular window length M narrows the lobe widths, areas under each lobe remains constant, thus the ripples occur more rapidly, but with the same amplitude. To reduce the height of the sidelobes, we will need to taper the rectangular window smoothly to zero at the ends. Some commonly used windows are displayed in OSB Figure 7.21. Their properties are quantitatively compared in OSB Table 7.1. Corresponding analytical definitions can be found in Section 7.2.1. Comparison among these commonly used windows shows that at a given length, the rectan­ gular window has the narrowest mainlobe; it gives the sharpest transition at the discontinuity of Hd (ejω ). The Bartlett window is triangular, with a slope discontinuity at its center, while the Hanning, Hamming and Blackman windows are smoother. By tapering the window smoothly to zero, the sidelobes can be reduced in amplitudes, but the trade-off is larger main lobes. All of these windows are symmetric about M/2, hence their frequency responses have generalized linear phase. The Kaiser windows, as defined below, are a family of near optimal windows that allow controlled trade-offs between the sidelobe amplitudes and mainlobe widths. By varying the β parameter, it is possible to approximately duplicate the other types of windows using a Kaiser window. The Kaiser window family: w[n] = � � � ⎧ 2 ⎪ ) ⎨ I0 β 1 − ( n−α α ⎪ ⎩ 0 I0 (β) 0≤n≤M otherwise where M , β : shape parameter 2 I0 (·): zeroth-order modified Bessel Function of the first kind α= The definition implies that the Kaiser class of windows is parameterized by length (M + 1) and shape (β). The following figure illustrates its general form. As β increases, I0 (β) also increases to give a smoother window. 2 Kaiser Window Given a set of filter specifications, the values of M and β needed can be determined numerically. As in OSB Figure 7.23, define: Transition bandwidth: Δω = ωs − ωp Peak approximation error: δ (passband: 1 ± δ stopband: ± δ) Stopband attenuation: A = −20 log10 δ Then the values of β, M needed to achieve A are ⎧ A > 50 ⎨ 0.1102(A − 8.7) 0.5842(A − 21)0.4 + 0.07886(A − 21) 21 ≤ A ≤ 50 β= ⎩ 0 A < 21 , M≥ A−8 . 2.285Δω See OSB Section 7.3 for examples of the Kaiser window method. Optimum FIR Filter Approximation: The Parks-McClellan Algorithm Although the rectangular windowing method provides the best mean-square approximation to a desired frequency response, just because it is optimal does not mean it is good. FIR filters designed with windows exhibit oscillatory behavior around the discontinuity of the ideal fre­ quency response and does not allow separate control of the passband and stopband ripples. An alternative FIR design technique is the Parks-McClellan algorithm which is based on polynomial approximations. Filter Design as Polynomial Approximation Consider the DTFT of a causal FIR system of length M + 1: jω H(e ) = M � h[n]e −jωn n=0 = M � n=0 3 h[n]xn |x=e−jω . Let’s further constrain our system to a type I generalized linear-phase filter defined by a sym­ metric impulse response with M even: h[n] = h[M − n] 0 ≤ n ≤ M, M even. Its corresponding frequency response is H(ejω ) = M � h[n]e−jωn = Ae (ejω )ejωM/2 , n=0 where M/2 � jω Ae (e ) = he [n]ejωn he [n] = h[n + M/2] = he [−n] n=−M/2 M/2 = he [0] + � he [n](e−jωn + ejωn ) n=1 M/2 = he [0] + 2 � he [n] cos ωn n=1 Here h[n] is just he [n] with an integer delay of M/2. Rewrite the cos ωn terms using Chebyshev polynomials: Tn (x) = cos(n cos−1 x) ⇒ L=M/2 jω Ae (e ) = � k=0 cos ωn = Tn (cos ω) M/2 k ak (cos ω) = � ak xk |x=cos ω = P (x)|x=cos ω . k=0 Similar procedures can also be carried out for type II, III, IV systems (See Section 7.4 of OSB). In particular, when M is odd: �ω � jω −jωM/2 H(e ) = e cos P (cos ω) 2 The Alternation Theorem With Ae (ejω ) written as a trignometric polynomial, Parks and McClellan applied the Alterna­ tion Theorem from approximation theory to the filter design problem. The theorem is stated in Section 7.4 of OSB (p.489). Although the alternation theorem does not lead to any optimal filter design directly, what it allows us to do is to determine whether a filter response written as P (x) is uniquely optimal under the minimax error criterion, i.e., if it minimizes the maximum weighted approximation error. The Parks-McClellan Algorithm uses the alternation theorem to iteratively determine an optimal equiripple approximation. 4 OSB Figures 7.35 and 7.36 show a typical example of an optimal lowpass filter frequency response for L = 7 and its equivalent polynomial approximation function P (x) . The alternating points are ω1 , ω2 , . . . , ω8 , π, where Ep (x) = Ep (cos ω) alternates between its maximum and minimum for L + 2 times. Such an approximation is said to be equiripple. The weighting function in this case satisfies: � 1 K cos ωp ≤ cos ω ≤ 1 , Wp (cos ω) = 1 −1 ≤ cos ω ≤ ωs where K = δ1 /δ2 . P (x) turns into a hyperbolic cosine function outside of Fp = [−1, 1]. Since we are only interested in ω ∈ [0, π], we do not care about the behavior of P (x) outside of Fp . Comparison between the Figures 7.35 and 7.36 shows that for lowpass filters (or any piecewise constant filters), the alternations in P (x) can be counted directly from the frequency response Ae (ejω ). In OSB Figure 7.37 are other examples of optimum lowpass filter approximations for L = 7. All four satisfy the alternation theorem. Nonetheless, the uniqueness of optimal filters implied by the alternation theorem is not violated because ωp and/or ωs are different for each. By comparison, OSB Figure 7.38 and 7.39 show approximations that are not optimal according to the alternation theorem. There are 8 alternations in each case. More generally, it can be shown that for type I lowpass filters, the following properties are implied by the alternation theorem: • The maximum possible number of alternations of the error is L + 3. • Alternations will always occur at ωp and ωs . • The passband and stopband have equal ripples, and the transition band is monotonic. All points with zero slope inside the passband and all points with zero slope inside the stopband (for 0 < ω < ωp and ωs < ω < π) will correspond to alternations, i.e., the filter is equiripple, except possibly at ω = 0, π. See Section 7.4 of OSB for proofs. So far we have focused on only optimal lowpass FIR filters, which have two approximation bands. For multi-band filters such as the bandpass filter, the alternation theorem still holds, but it does imply a different set of properties on the approximating polynomial. For example, the optimal approximation can have more than L + 3 alternations and have ripples (local maxima/minima) in the transition band. OSB Figure 7.47 is such an example. This system is optimal since the alternation theorem is satisfied, but the ripple in its transition band makes it undesirable in practical applications. Sections 7.5 and 7.6 of OSB give more examples of FIR equiripple approximations and further discuss practical implementation issues that should be considered during the filter design process. Parks-McClellan Algorithm The alternation theorem provides us with a tool for efficiently finding the optimum filter in the minimax sense. For a type I lowpass filter, the following set of equations is satisfied by the 5 optimum filter according to the alternation theorem: W (ωi )[Hd (ejωi ) − Ae (ejωi )] = (−1)i+1 δ i = 1, 2, . . . , (L + 2) , where δ is the optimum error. Parks and McClellan found that given fixed ωp , ωs , an iterative algorithm finds the optimal approximating Ae (ejωi ): Algorithm 1. Initial guess of (L + 2) extremal frequencies ω0 , ω1 , . . . , ωL+1 . 2. Solve for polynomial coefficients and (L + 2) equations in (L + 2) unknowns. 3. Evaluate Ae (ejωi ) or E(ω) on a dense set of frequencies and choose a new set of extremal frequen­ cies. 4. Check whether extremal frequencies changed. If yes go to step 2. If no, algorithm has converged. This algorithm is described in more detail and illustrated by a flowchart in OSB Section 7.4.3. Keep in mind that the optimality of filters designed by the PM algorithm is only in the minimax error sense, and the equiripple characteristic may be less desirable when compared with filters designed by other techniques. 6