數(shù)字信號處理實(shí)驗7

1、 Laboratory Exercise 7DIGITAL FILTER DESIGN7.1DESIGN OF IIR FILTERSProject 7.1Estimation of IIR Filter OrderAnswers:Q7.1 The normalized passband edge angular frequency Wp is -0.2The normalized stopband edge angular frequency Ws is -0.4The desired passband ripple Rp is -0.5dBThe desired stopband ripp

2、le Rs is -40dB(1) Using these values and buttord we get the lowest order for a Butterworth lowpass filter to be - 8The corresponding normalized passband edge frequency Wn is - 0.2469 or 0.2469pi(2) Using these values and cheb1ord we get the lowest order for a Type 1 Chebyshev lowpass filter to be -5

3、The corresponding normalized passband edge frequency Wn is - 0.2000(3) Using these values and cheb2ord we get the lowest order for a Type 2 Chebyshev lowpass filter to be -5N, Wn = cheb2ord(0.2,0.4,0.5,40).The corresponding normalized passband edge frequency Wn is - 0.4000(4) Using these values and

4、ellipord we get the lowest order for an elliptic lowpass filter to be 4N, Wn = ellipord(0.2,0.4,0.5,40).From the above results we observe that the Elliptic filter has the lowest order meeting the specifications.Q7.2 The normalized passband edge angular frequency Wp is - 0.6000The normalized stopband

5、 edge angular frequency Ws is - 0.3429The desired passband ripple Rp is -1dBThe desired stopband ripple Rs is -50dB(1) Using these values and buttord we get the lowest order for a Butterworth highpass filter to be -8N, Wn = buttord(Wp,Ws,Rp,Rs).The corresponding normalized passband edge frequency Wn

6、 is - 0.5647(2) Using these values and cheb1ord we get the lowest order for a Type 1 Chebyshev highpass filter to be 5N,Wn = cheb1ord(Wp,Ws,Rp,Rs).The corresponding normalized passband edge frequency Wn is - 0.6000(3) Using these values and cheb2ord we get the lowest order for a Type 2 Chebyshev hig

7、hpass filter to be -5N,Wn = cheb2ord(Wp,Ws,Rp,Rs).The corresponding normalized passband edge frequency Wn is 0.3429(4) Using these values and ellipord we get the lowest order for an elliptic highpass filter to be -4N,Wn = ellipord(Wp,Ws,Rp,Rs).The corresponding normalized passband edge frequency Wn

8、is Wn = 0.6000,From the above results we observe that the Elliptic filter has the lowest order meeting the specifications.Project 7.2 IIR Filter DesignA copy of Program P7_1 is given below:%程序p7.1%巴特沃斯帶阻濾波器的設(shè)計ws=0.4 0.6;wp=0.3 0.7;rp=0.4;rs=50;%估計濾波器的階數(shù)N1,wn1=buttord(wp,ws,rp,rs);%設(shè)計濾波器num,den=butte

9、r(N1,wn1,'stop');%顯示傳輸函數(shù)disp('分子系數(shù)是');disp(num);disp('分母系數(shù)是');disp(den);%計算增益響應(yīng)g,w=gain(num,den);%繪制增益響應(yīng)plot(w/pi,g);gridaxis(0 1 -60 5);xlabel('omega/pi');ylabel('增益，dB');title('巴特沃斯帶阻濾波器的設(shè)計');Answers:Q7.5 The coefficients of the Butterworth bandstop

10、transfer function generated by running Program P7_1 are as follows:分子系數(shù)是 Columns 1 through 6 0.0330 0.0000 0.2972 0.0000 1.1889 0.0000 Columns 7 through 12 2.7741 0.0000 4.1611 0.0001 4.1611 0.0000 Columns 13 through 18 2.7741 0.0000 1.1889 0.0000 0.2972 0.0000 Column 19 0.0330分母系數(shù)是 Columns 1 throug

11、h 6 1.0000 0.0000 2.6621 0.0000 4.1451 0.0001 Columns 7 through 12 4.1273 0.0001 2.8977 0.0000 1.4381 0.0000 Columns 13 through 18 0.5027 0.0000 0.1178 0.0000 0.0167 0.0000 Column 19 0.0011The exact expression for the transfer function is The gain response of the filter as designed is given below:Fr

12、om the plot we conclude that the design 符合the specifications. The plot of the unwrapped phase response and the group delay response of this filter is given below:Here is the program to find and plot the unwrapped phase response and group delay:% 程序 Q7_5B% 巴特沃斯帯阻濾波器的設(shè)計% Plot the unwrapped phase and t

13、he group delay.Ws = 0.4 0.6; Wp = 0.2 0.8; Rp = 0.4; Rs = 50;% 估計濾波器節(jié)數(shù)N1, Wn1 = buttord(Wp, Ws, Rp, Rs);% 設(shè)計濾波器num,den = butter(N1,Wn1,'stop');%顯示傳輸函數(shù)，繪制不卷繞相位wp = 0:pi/1023:pi;wg = 0:pi/511:pi;Hz = freqz(num,den,wp);Phase = unwrap(angle(Hz);figure(1);plot(wp/pi,Phase);grid;% axis(0 1 a b);xl

14、abel('omega /pi'); ylabel('Unwrapped Phase (rad)');title('Unwrapped Phase Response of a Butterworth Bandstop Filter');% 全延時GR = grpdelay(num,den,wg);figure(2);plot(wg/pi,GR);grid;%axis(0 1 a b);xlabel('omega /pi'); ylabel('Group Delay (sec)');title('Group

15、Delay of a Butterworth Bandstop Filter');Project 7.4Estimation of FIR Filter OrderAnswers:Q7.13The estimated order of a linear-phase lowpass FIR filter with the following specifications: wp = 2 kHz, ws = 2.5 kHz, dp = 0.005, ds = 0.005, and FT = 10 kHz obtained using kaiord is -The purpose of th

16、e command ceil is - To round the estimated order up to the next largest integer; the order has to be integer, so if the formula returns a fraction it needs to be rounded up to the next whole number. The purpose of the command nargin is -To detect if kaiord has been called with fourarguments or with

17、five. If five, its assumed that all the frequencies are analog and thatthe last argument is the sampling frequency. If four, then the sampling frequencydefaults to 2, implying that the other frequency arguments are in units of cycles per sample.Project 7.5FIR Filter DesignAnswers:Q7.20The MATLAB pro

18、gram to design and plot the gain and phase responses of a linear-phase FIR filter using fir1 is shown below. The filter order is estimated using kaiserord. The output data are the filter coefficients. % 程序Q7_20% 線性相位FIR低通濾波器的設(shè)計% 滿足Q7.13的參數(shù)% 寫出分子系數(shù)% 傳輸函數(shù)% - 增益響應(yīng)% - 相位響應(yīng)% - 繪制未卷繞相位響應(yīng)%clear;Fp = 2*103;

19、Fs = 2.5*103;FT = 10*103;Rp = 0.005;Rs = 0.005;%估計階數(shù)N = kaiord(Fp,Fs,Rp,Rs,FT)Wp = 2*Fp/FT; % 規(guī)劃角頻率Ws = 2*Fs/FT; Wn = Wp + (Ws - Wp)/2;h = fir1(N,Wn);disp('Numerator Coefficients are ');disp(h);g, w = gain(h,1); figure(1);plot(w/pi,g);grid;%axis(0 1 -60 5);xlabel('omega /pi'); ylabel

20、('Gain in dB');title('Gain Response');% 頻率響應(yīng)w2 = 0:pi/511:pi;Hz = freqz(h,1,w2);MagH = abs(Hz);T1 = 1.005*ones(1,length(w2);T2 = 0.995*ones(1,length(w2);T3 = 0.005*ones(1,length(w2);figure(4);plot(w2/pi,MagH,w2/pi,T1,w2/pi,T2,w2/pi,T3);grid;% 相位figure(2);Phase = angle(Hz);plot(w2/pi,

21、Phase);grid;xlabel('omega /pi'); ylabel('Phase (rad)');title('Phase Response');figure(3);UPhase = unwrap(Phase);plot(w2/pi,UPhase);grid;xlabel('omega /pi'); ylabel('Unwrapped Phase (rad)');title('Unwrapped Phase Response');The coefficients of the lowpa

22、ss filter corresponding to the specifications given in Question 7.20 are as shown below 0.0010 -0.0004 -0.0015 0.0000 0.0024 0.0010 -0.0038 -0.0032 0.00490.0071 -0.0050 -0.0128 0.0026 0.0202 0.0038 -0.0284 -0.0166 0.03660.0404 -0.0436 -0.0909 0.0483 0.3129 0.4498 0.3129 0.0483 -0.0909-0.0436 0.0404

23、0.0366 -0.0166 -0.0284 0.0038 0.0202 0.0026 -0.0128-0.0050 0.0071 0.0049 -0.0032 -0.0038 0.0010 0.0024 0.0000 -0.0015-0.0004 0.0010The generated gain and phase responses are given below: From the gain plot we observe that the filter as designed _DOES NOT_ meet the specifications.AS SHOWN in the two

24、detail plots above, with N=46, neither the passband spec at wp= 0.4 (normalized frequency) nor the stopband spec at ws = 0.5 (normalized frequency) are met. So this design DOES NOT meet the spec.The filter order that meets the specifications is N=66Q7.23The MATLAB program to design and plot the gain

25、 and phase responses of a linear-phase FIR filter using fir1 and kaiser is shown below. The filter order N is estimated using Eq. (7.37) and the parameter b is computed using Eq. (7.36). The output data are the filter coefficients.% 程序 Q7_23% 使用凱澤共識設(shè)計FIR低通濾波器% - 增益響應(yīng)%clear;% 參數(shù)滿足Q7.23.Wp = 0.31;Ws =

26、 0.41;Wn = Wp + (Ws-Wp)/2;As = 50;Ds = 10(-As/20);Dp = Ds; if As > 21N = ceil(As-7.95)*2/(14.36*(abs(Wp-Ws)+1)elseN = ceil(0.9222*2/abs(Wp-Ws)+1)end% (7.36) 計算Bif As > 50BTA = 0.1102*(As-8.7);elseif As >= 21BTA = 0.5842*(As-21)0.4+0.07886*(As-21);elseBTA = 0;endWin = kaiser(N+1,BTA);h = fir

27、1(N,Wn,Win);% 分子系數(shù)disp('Numerator Coefficients are ');disp(h);%增益響應(yīng)g, w = gain(h,1);figure(1);plot(w/pi,g);grid;axis(0 1 -80 5);xlabel('omega /pi'); ylabel('Gain in dB');title('Gain Response');%頻率響應(yīng)w2 = 0:pi/511:pi;Hz = freqz(h,1,w2);figure(2);Phase = angle(Hz);plot(w2/pi,Phase);grid;xlabel('omega /pi'); ylabel('Phase (rad)');title('Phase Response');figure(3);UPhase = u