This preview shows page 1. Sign up to view the full content.
Unformatted text preview: Filtering References: 1. Biosignal and Biomedical Image Processing MATLA B‐Based Applica>ons, JOHN L. SEMMLOW 2. Bores Signal Processing ‐ Introduc>on to DSP ‐ Filtering hKp://
www.bores.com/courses/intro/ﬁlters/4_window.htm The goal of ﬁltering • To reshape the spectrum to one’s advantage – Ex: to improve the SNR • It can be considered as a process of selec>ng, or suppressing, certain frequency components of a signal • The general concept of ﬁltering is linear processing The characteris>cs of signal and noise • Signals: most signal are narrow band • Noise : broadband • The broadest band noise is white noise For example The Time waveform of a recorded signal FFT Filtering Frequency based analysis • Applied on periodic or aperiodic signals • It can’t be applied on transient Reponses of inﬁnite length, ex: step func>on or system with nonzero ini>al condi>ons Analog domain: Laplace transform (s, s=σ+jϖ), a representa>on of complex frequency in place of jw in the Fourier transform Digital domain: Z‐transform (z, |z|ejw) Three types of ﬁlters based on their ﬁltering characteris>cs passband the band of frequency components
that are allowed to pass stopband the band of frequency components
that are suppressed passband ripple the maximum amount by which
attenuation in the passband may
deviate from nominal gain stopband attenuation the minimum amount by which
frequency components in the
stopband are attenuated the stopband aKenua>on is formally speciﬁed as the aKenua>on to the top of the ﬁrst sidelobe of the ﬁlter's frequency response. The groups of ﬁlters based on their approach: • Finite impulse response (FIR) ﬁlters (digital ﬁltering) • Inﬁnite impulse response (IIR) ﬁlters (similar to analog ﬁltering) Digital transfer func>on • Is the most useful applica>on of z‐transform line in its ability to deﬁne the digital equivalent of a transfer func>on • It is deﬁned as: X
a Z‐transform The z transform is deﬁned as a sum of signal values x[n] mul>plied by powers of z Z: arbitrary complex number • A interval shia of n samples, or associated >me shia of nTs seconds • Every data sample in sequence x(n) is associated with a unique power of z which deﬁnes a sample’s posi>on in the sequence • The >me shiaing property of z−n can be formally stated as: • The >me shiaing characteris>c of Z‐transform can be used to deﬁne a unit delay process, Z‐1 • Z‐transform: a unit delay process • Ex: The output is the same as input but shiaed (or delayed) by one data sample An earlier signal value can be generated from the present one, because the z transform of x[n‐1] is just the z transform of x[n] mul>plied by(Z‐1): • Applying the operator z‐1 to an input value (say xn) gives the previous input (xn‐1): z‐1 xn = xn‐1 Suppose we have an input sequence x0 = 5 x1 = ‐2 x2 = 0 x3 = 7 x4 = 10 Then z‐1 x1 = x0 = 5 Then z‐1 x2 = x1 = ‐2 z‐1 x3 = x2 = 0 z‐1 x0 = x‐1 (usually =0) • Similarly, applying the operator z‐1 to another input value (say Yn) gives the previous input (Yn‐1): z‐1 yn = yn‐1 • Applying the delay operator z‐1 twice produces a delay of two sampling intervals: z‐1 (z‐1 xn) = z‐1 xn‐1 = xn‐2 z‐2 z‐2 xn = xn‐2 This nota>on can be extended to delays of three or more sampling intervals, the appropriate power of z‐1 being used. • now applying this nota>on in the descrip>on of a recursive digital ﬁlter A second order digital ﬁlter, given it symmetry form by : a0yn + a1yn‐1 + a2yn‐2 = b0xn + b1xn‐1 + b2xn‐2 yn-1 = z-1 yn yn-2 = z-2 yn xn-1 = z-1 xn xn-2 = z-2 xn (a0 + a1z‐1 + a2z‐2) yn = (b0 + b1z‐1 + b2z‐2) xn yn / xn = (b0 + b1z‐1 + b2z‐2) / (a0 + a1z‐1 + a2z‐2) the general form of the transfer func;on for a second‐order recursive (IIR) ﬁlter For a non‐recursive ﬁlter (FIR) b0 is regarded as 1, and all the other b coeﬃcients are zero Then yn / xn = a0 + a1z‐1 + a2z‐2 Digital ﬁlters are also characterized by their response to an impulse: a signal consis>ng of a single value followed by zeroes: • The impulse response is an indica>on of how long the ﬁlter takes to seKle into a steady state( the ﬁlter's stability ) • The impulse response that con>nues oscilla>ng in the long term indicates the ﬁlter may be prone to instability. Matlab Implementa>on •
• Digital transfer func>on Func>on ﬁlter y= ﬁlter (b,a, x) y: output, b,a: coeﬃcient of H(z), x: input Example Fs=1000; data length=512 The impulse responses of H(Z) Exercise Fs=1000; data length=512 Func>on freqz • Digital ﬁlter frequency response • [H,W] = FREQZ(b,a,N) returns the N‐point complex frequency response vector H and the N‐point frequency vector W in radians/sample of the ﬁlter: jw ‐jw ‐jmw jw B(e) b(1) + b(2)e + .... + b(m+1)e H(e) = ‐‐‐‐ = ‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐‐ jw ‐jw ‐jnw A(e) a(1) + a(2)e + .... + a(n+1)e given numerator and denominator coeﬃcients in vectors B and A. The frequency response is evaluated at N points equally spaced around the upper half of the unit circle. If N isn't speciﬁed, it defaults to 512. • [H,F] = freqz(b,a,N,Fs) and [H,F] = freqz(B,A,N, Fs) return frequency vector F (in Hz), where Fs is the sampling frequency (in Hz). • Filtering can be done directly in the frequency domain, by opera>ng on the signal's frequency spectrum A digital ﬁlter is made up from previous inputs and previous outputs involved with two convolu>ons • Since ﬁltering is a frequency selec>ve process, the important thing about a digital ﬁlter is its frequency response The ﬁlter's frequency response • The ﬁlter's frequency response can be calculated from its ﬁlter equa>on: ? Coeﬃcient J: image part, square root of minus one to calculate the ﬁlter coeﬃcients having ﬁrst deﬁned the desired frequency response. So we are faced with an inverse problem Prac>cally, it requires the least possible amount of computa>on ‐> using the smallest number of coeﬃcients. To simplify the calcula>on by excluding the possibility of feedback of ﬁlter The impulse response of the above ﬁlter must necessarily be ﬁnite in dura>on Finite impulse response (FIR) ﬁlter • The impulse response of the above ﬁlter must necessarily be ﬁnite in dura>on • The filter's frequency response is also simplified,
because all the bottom half is omitted this frequency response is just the Fourier transform of the ﬁlter coeﬃcients Therefore, the coeﬃcients for an FIR ﬁlter can be calculated simply by taking the inverse Fourier transform of the desired frequency response BUT... • to deﬁne a sharp ﬁlter needs closely spaced frequency samples ‐ so a lot of them • so the inverse Fourier transform will give us a lot of ﬁlter coeﬃcients • but we don't want a lot of ﬁlter coeﬃcients So, a beKer way for calcula>ng FIR ﬁlter coeﬃcients based on throwing away the small ones: • To hypothesize that it would be s>ll ok by removing lots of ﬁlter coeﬃcients (small values) • To specify the desired frequency response using lots of samples • To calculate the inverse Fourier transform • But, this gives us a lot of ﬁlter coeﬃcients • so to truncate the ﬁlter coeﬃcients to give us less • then to calculate the Fourier transform of the truncated set of coeﬃcients to see if it s>ll matches our requirement But, • trunca>ng the ﬁlter coeﬃcients means the ﬁlter's frequency response can only be deﬁned coarsely • Luckily, we already know a way to sharpen up the frequency spectrum of a truncated signal, by applying a window func>on. • So aaer trunca>on, we can apply a window func>on to sharpen up the ﬁlter's frequency response A beKer recipe for calcula>ng FIR ﬁlter coeﬃcients • To pretend we don't mind lots of ﬁlter coeﬃcients • To specify the desired frequency response using lots of samples • To calculate the inverse Fourier transform • But, this gives us a lot of ﬁlter coeﬃcients • So to truncate the ﬁlter coeﬃcients to give us less • To apply a window func>on to sharpen up the ﬁlter's frequency response • Then to calculate the Fourier transform of the truncated set of coeﬃcients to see if it s>ll matches our requirement This is called the window method of FIR ﬁlter design. BUT... • The stopband attenuation is fixed in most
window func>ons An appropriate choice of window func>on is an important step for designing the FIR ﬁlter Matlab implementa>on • FIR ﬁlter design: – Two stage – Three stage Two stage FIR ﬁlter design • Determine the coeﬃcient, b by know order (number of coeﬃcient) • Func>on ﬁr1 – FIR ﬁlter design using the window method – b=ﬁr1(n,wn,’aype’, window) – Wn: 0‐1; normalized by fs/2 – fytpe: high‐> high pass stop‐> stopband ﬁlter non‐sepciﬁed‐> may low or band pass according to the wn Window: window length= n+1 Example • Design a window‐based bandpass ﬁlter with cutoﬀ frequencies of 5 and 15 Hz. Assume a sampling frequency of 100 Hz; Filter order = 128, data length=512; Matlab code for example 2 Three stage of FIR ﬁlter (for matlab 6.0 only) • To determine the order • To determine the coeﬃcient • To apply the designed ﬁlter on the signals • To determine the order func>on remezorder • To determine the coeﬃcient func>on remez • To apply the designed ﬁlter on the signals func>on ﬁlter Inﬁnite impulse response (IIR) ﬁlter. • A ﬁlter is subjected to an impulse and its output need not necessarily become zero aaer the impulse has run through the summa>on • The impulse response of such a ﬁlter can be inﬁnite in dura>on. IIR • The lter's frequency response can be simpliﬁed by subs>tu>ng a new variable, z Complex numbers can be drawn using an Argand diagram The complex variable z is shown as a vector on the Argand diagram Unit circle • Z, traces a circle of radius 1 on the Argand diagram. • The map of values in the z plane is called the transfer func>on H(z). • The transfer func>on H(z) evaluated around the unit circle on the Argand diagram of z • The shape of the H(z) is determined by the posi>ons of its poles and zeroes: if z is inside the unit circle, the transient terms will die away if z is on the unit circle, oscilla>ons will be in a steady state if z is outside the unit circle, the transient terms will increase Zero and Pole =0, zeros =0, poles, inﬁnitely large, • The frequency response can be determined by tracing around the unit circle of the z plane: – poles cause bumps – zeroes cause dips – the closer to the unit circle, the sharper the feature Impulse invariance • No one knows how to design a IIR ﬁlter • The beKer way to do is to fall back on the analogue ﬁlter • Then to calculate and sampled its impulse responses • To use the results as the ﬁlter coeﬃcient • The above process is called as the impulse invariance • Two problems for using this method – Aliasing and Frequency resolu>on MATLAB implementa>on • func>on ﬁlqilt • It improves the phase performance of IIR ﬁlters by running the algorithm in both forward and reverse direc>ons • y = ﬁlKilt(b,a,x) Two stage FIR ﬁlter design • Yule–Walker recursive ﬁlter is the only IIR ﬁlter that is not supported by an order‐selec>on rou>ne (i.e., a three‐stage design process). • [b,a] = yulewalk(n,f,m) • n: order • f are relaSve to fs/2,the ﬁrst point in f must be zero and the last point 1. • m is a vector of the desired ﬁlter gains at the frequencies speciﬁed in f Example • Design an 12th‐order Yule–Walker bandpass ﬁlter with cutoﬀ frequencies of 0.25 and 0.5. Plot the frequency response of this ﬁlter and compare with that produced by the FIR ﬁlter ﬁr2 of the same order (fs=256, data length=256) Exercise • Design an 14 th‐order Yule–Walker bandpass ﬁlter with cutoﬀ frequencies of 0.3 and 0.8. Plot the frequency response of this ﬁlter and compare with that produced by the FIR ﬁlter ﬁr2 of the same order (fs=256, data length=256) Three stage IIR ﬁlter design • 1. to determine the order – BuKerworth, Chebyshev, and ellip>c—are supported by order selec>on rou>nes – [n,wn] = buUord(wp, ws, rp, rs); BuUerworth ﬁlter wp : the passband frequency rela>ve to fs/2, ws : the stopband frequency rela>ve to fs/2 rp : the passband ripple in db, rs : the stopband ripple in db • 2. to determine the coeﬃcient [b,a] – [n,wn] = buUord(wp, ws, rp, rs); BuUerworth ﬁlter • 3. Filtering the signal with the designed IIR ﬁlter Example • Example 4.8 Plot the frequency response curves (in db) obtained from an 8th‐order lowpass ﬁlter using the BuKerworth, Chebyshev Type I (func>on cheby1) and II (cheby2), and ellip>c ﬁlters (ellip). Use a cutoﬀ frequency of 200 Hz and assume a sampling frequency of 2 kHz. For all ﬁlters, the passband ripple should be less than 3 db and the minimum stopband aKenua>on should be 60 db. Homework • 1. Find and plot the frequency spectrum and the impulse response of a digital linear process having the digital transfer func>on: Fs=1000; data length=512 2. a) To Load the ecg data from previous homework b) To construct a noise consisSng of two closely spaced sinusoids of 200 and 230 Hz with SNR of −8 db and −12 db respec>vely. c) To Plot the magnitude spectrum using the FFT. d) To apply an 24 coeﬃcient FIR bandpass window type ﬁlter to the data using the ﬁr1 MATLAB rouSne. e) To plot the bandpass ﬁltered data. • Use the ﬁrst stage IIR design rou>nes, buUord, cheby1ord, cheby2ord, and elliptord to ﬁnd the ﬁlter order required for a lowpass ﬁlter that aUenuates 40 db/octave. (An octave is a doubling in frequency). Assume a cutoﬀ frequency of 200 Hz and a sampling frequency of 2 kHz. a) To Load the ecg data from previous homework b) To construct a noise consisSng of two closely spaced sinusoids of 60 Hz and a zero min unit variance random noise with SNR of −8 db and −12 db respec>vely. c) Use the three step IIR ﬁlter design method, with passband ripple<1.5dB, stopband aKenua>on> 60dB • Plot the spectra of the signals without and with the IIR ﬁltering ...
View Full Document
- Spring '09