Biosignal and Biomedical Image Processing MATLAB-Based Applications Muya phần 4

94 Chapter 4 replaced by the filter coefficients, b(n). Hence, FIR filters can be implemented using either convolution or MATLAB’s filterroutine. Eq. (8) indicates that the filter coefficients (or weights) of an FIR filter are the same as the impulse response of the filter. Since the frequency response of a process having an im-pulse response h(n) is simply the Fourier transform of h(n), the frequency re-sponse of an FIR filter having coefficients b(n) is just the Fourier transform of b(n): N−1 X(m) = b(n) e(−j2π mn/N) (9) n=0 Eq. (9) is a special case of Eq. (5) when the denominator equals one. If b(n) generally consists of a small number of elements, this equation can some-times be determined manually as well as by computer. The inverse operation, going from a desired frequency response to the coefficient function, b(n), is known as filter design. Since the frequency re-sponse is the Fourier transform of the filter coefficients, the coefficients can be found from the inverse Fourier transform of the desired frequency response. This design strategy is illustrated below in the design of a FIR lowpass filter based on the spectrum of an ideal filter. This filter is referred to as a rectangular window filter* since its spectrum is ideally a rectangular window. FIR Filter Design The ideal lowpass filter was first introduced in Chapter 1 as a rectangular win-dow in the frequency domain (Figure 1.7). The inverse Fourier transform of a rectangular window function is given in Eq. (25) in Chapter 2 and repeated here with a minor variable change: sin[2πfcTs(n − L/2)] π(n − L/2) (10) where fc is the cutoff frequency; Ts is the sample interval in seconds; and L is the length of the filter. The argument, n − L/2, is used to make the coefficient function symmetrical giving the filter linear phase characteristics. Linear phase characteristics are a desirable feature not easily attainable with IIR filters. The coefficient function, b(n), produced by Eq. (10), is shown for two values of fc in Figure 4.4. Again, this function is the same as the impulse response. Unfortu- *This filter is sometimes called a window filter, but the term rectangular window filter will be used in this text so as not to confuse the filter with a window function as described in the last chapter. This can be particularly confusing since, as we show later, rectangular window filters also use window functions! TLFeBOOK Digital Filters 95 FIGURE 4.4 Symmetrical weighting function of a rectangular filter (Eq. (10) trun-cated at 64 coefficients. The cutoff frequencies are given relative to the sampling frequency, fs, as is often done in discussing digital filter frequencies. Left: Low-pass filter with a cutoff frequency of 0.1fs/2 Hz. Right: Lowpass cutoff frequency of 0.4fs/2 Hz. nately this coefficient function must be infinitely long to produce the filter char-acteristics of an ideal filter; truncating it will result in a lowpass filter that is less than ideal. Figure 4.5 shows the frequency response, obtained by taking the Fourier transform of the coefficients for two different lengths. This filter also shows a couple of artifacts associated with finite length: an oscillation in the frequency curve which increases in frequency when the coefficient function is longer, and a peak in the passband which becomes narrower and higher when the coefficient function is lengthened. Since the artifacts seen in Figure 4.5 are due to truncation of an (ideally) infinite function, we might expect that some of the window functions described in Chapter 3 would help. In discussing window frequency characteristics in Chapter 3, we noted that it is desirable to have a narrow mainlobe and rapidly diminishing sidelobes, and that the various window functions were designed to make different compromises between these two features. When applied to an FIR weight function, the width of the mainlobe will influence the sharpness of the transition band, and the sidelobe energy will influence the oscillations seen TLFeBOOK 96 Chapter 4 FIGURE 4.5 Freuquency characteristics of an FIR filter based in a weighting func-tion derived from Eq. (10). The weighting functions were abruptly truncated at 17 and 65 coefficients. The artifacts associated with this truncation are clearly seen. The lowpass cutoff frequency is 100 Hz. in Figure 4.5. Figure 4.6 shows the frequency characteristics that are produced by the same coefficient function used in Figure 4.4 except that a Hamming window has been applied to the filter weights. The artifacts are considerably diminished by the Hamming window: the overshoot in the passband has disap-peared and the oscillations are barely visible in the plot. As with the unwin-dowed filter, there is a significant improvement in the sharpness of the transition band for the filter when more coefficients are used. The FIR filter coefficients for highpass, bandpass, and bandstop filters can be derived in the same manner from equations generated by applying an inverse FT to rectangular structures having the appropriate associated shape. These equations have the same general form as Eq. (10) except they include additional terms: sin[π(n − L/2)] sin[2πfcTs(n − L/2)] π(n − L/2) π(n − L/2) Highpass (11) sin[2πfHT(n − L/2)] sin[2πfLTs(n − L/2)] π(n − L/2) π(n − L/2) Bandpass (12) TLFeBOOK Digital Filters 97 FIGURE 4.6 Frequency characteristics produced by an FIR filter identical to the one used in Figure 4.5 except a Hamming function has been applied to the filter coefficients. (See Example 1 for the MATLAB code.) sin[2πfLT(n − L/2)] sin[π(n − L/2)] π(n − L/2) π(n − L/2) sin[2πfHTs(n − L/2)] π(n − L/2) Bandstop (13) An FIR bandpass filter designed using Eq. (12) is shown in Figure 4.7 for two different truncation lengths. Implementation of other FIR filter types is a part of the problem set at the end of this chapter. A variety of FIR filters exist that use strategies other than the rectangular window to construct the filter coef-ficients, and some of these are explored in the section on MATLAB implemen-tation. One FIR filter of particular interest is the filter used to construct the derivative of a waveform since the derivative is often of interest in the analysis of biosignals. The next section explores a popular filter for this operation. Derivative Operation: The Two-Point Central Difference Algorithm The derivative is a common operation in signal processing and is particularly useful in analyzing certain physiological signals. Digital differentiation is de- TLFeBOOK 98 Chapter 4 FIGURE 4.7 Frequency characteristics of an FIR Bandpass filter with a coefficient function described by Eq. (12) in conjuction with the Blackman window function. The low and high cutoff frequencies were 50 and 150 Hz. The filter function was truncated at 33 and 129 coefficients. These figures were generated with code similar to that in Example 4.2 below, except modified according to Eq. (12) fined as x/ t and can be implemented by taking the difference between two adjacent points, scaling by 1/Ts, and repeating this operation along the entire waveform. In the context of the FIR filters described above, this is equiva-lent to a two coefficient filter, [−1, +1]/Ts, and this is the approach taken by MATLAB’s dervroutine. The frequency characteristic of the derivative opera-tion is a linear increase with frequency, Figure 4.8 (dashed line) so there is considerable gain at the higher frequencies. Since the higher frequencies fre-quently contain a greater percentage of noise, this operation tends to produce a noisy derivative curve. Figure 4.9A shows a noisy physiological motor response (vergence eye movements) and the derivative obtained using the dervfunction. Figure 4.9B shows the same response and derivative when the derivative was calculated using the two-point central difference algorithm. This algorithm acts TLFeBOOK ... - --nqh--