數字信號處理實驗報告

1、實驗報告課程名稱 數字信號處理 姓名： 學號：實驗二、基于MATLAB的數字濾波器的設計1實驗目的1）熟悉IIR數字濾波器設計的原理與方法。2）熟悉用窗函數法設計FIR數字濾波器設計的原理與方法。3）掌握IIR和FIR數字濾波器的計算機仿真方法。4）通過觀察對實際心電信號的濾波作用，獲得數字濾波的感性知識。2實驗原理IIR濾波器常用的一種設計方法是借助模擬濾波器的設計方法進行的，即，先設計模擬濾波器得到傳輸函數 ，然后將 按某種方法轉換成數字濾波器的系統(tǒng)函數 。轉換方法主要有脈沖響應不變法和雙線性變換法。脈沖響應不變法就是使數字濾波器的單位脈沖響應模仿模擬濾波器的沖激響應，由于時域采樣后造成頻

2、域的周期延拓，可能會造成頻域的混疊，因此不適合高通、帶阻濾波器的設計。另一種常用的方法是雙線性變換法。雙線性變換法為了克服s平面與z平面的非一一對應關系，采用了非線性頻率壓縮方法：首先將原s平面壓縮到s1平面的一條橫帶內，再將s1平面映射到z平面上。雖然雙線性變換法的模擬頻率與數字頻率間是非線性的關系，但消除了頻率的混疊。FIR濾波器常用的一種設計方法是窗口法。如果希望濾波器的頻率響應為 ，則其對應的單位脈沖響應為 。窗函數設計法的基本原理是用有限長單位脈沖響應序列h(n)逼近hd(n)。由于hd(n)往往是無限長序列，且是非因果的，所以用窗函數w(n)將hd(n)截斷，并進行加權處理，得到h

3、(n)= hd(n)w(n)，h(n)就作為實際設計的FIR數字濾波器的單位脈沖響應序列，其頻率響應函數為 ，式中N為所選窗函數w(n)的長度。用窗函數法設計的濾波器性能取決于窗函數的類型及窗口長度。設計中，要根據對阻帶最小衰減和過渡帶寬度的要求選擇合適的窗函數類型和窗口長度。選定窗函數類型和窗口長度后，求出單位脈沖響應和頻率響應。頻率響應是否滿足要求，要進行驗算，如果不滿足要求，則要重新選擇窗函數類型和長度，再次驗算，直至滿足要求為止。如果要求線性相位特性，則單位脈沖響應還必須滿足偶對稱或奇對稱特性。根據單位脈沖響應的奇偶對稱性和長度的奇偶性又將線性相位FIR濾波器分成四類，要根據所設計的濾

4、波特性正確選擇其中一類。3、 實驗內容1） 設計給定指標的Butterworth數字低通濾波器。設計指標參數為：在通帶內0，0.2頻率區(qū)間上，最大衰減小于1dB；在阻帶0.3，內，最小衰減大于15dB。并畫出濾波器的頻響特性曲線。2） 設計給定指標的FIR數字低通濾波器。3） 用所設計的數字濾波器對實際心電信號進行仿真濾波處理，分別打印出濾波前后的心電信號時域和頻域波形。4、 實驗步驟1） 復習有關Butterworth模擬濾波器設計和用雙線性變換法設計IIR數字濾波器的內容，用雙線性變換法設計數字濾波器的系統(tǒng)函數。2） 復習用窗函數法設計FIR數字濾波器的內容，并設計數字濾波器的系統(tǒng)函數。3

5、） 用Matlab中的函數畫出心電信號的時域和頻域波形。4） 用Matlab中的函數設計滿足要求的數字濾波器，畫出濾波器的頻響曲線。5） 用設計出的濾波器對心電信號進行濾波處理，畫出濾波后的信號時域與頻域波形。6） 用C語言編寫濾波器仿真程序，對心電信號進行濾波處理，將結果與步驟5的結果進行比較。5、實驗結果1）設計給定指標的Butterworth數字低通濾波器:代碼：clcN,wc= buttord(0.2,0.3,1,15);b,a= butter(N,wc)設計結果：b = 0.0007 0.0044 0.0111 0.0148 0.0111 0.0044 0.0007a = 1.000

6、0 -3.1836 4.6222 -3.7795 1.8136 -0.4800 0.05442）設計給定指標的FIR數字低通濾波器:代碼：clc;wp=0.2*pi;ws=0.3*pi;deltaw=ws-wp;N0=ceil(6.6*pi/deltaw);%計算海明窗長度N=N0+mod(N0+1,2);wd=(hamming(N); wc=(wp+ws)/2;t=(N-1)/2;n=0:1:N-1;hd=sin(wc*(n-t+eps)./(pi*(n-t+eps);h=hd.*wdA,w=freqz(h);mag=abs(A);db=20*log10(mag+eps)/max(mag);

7、plot(w/pi,db);title(幅度響應(db);xlabel(w/pi);ylabel(20log|H(ejw|(db);axis(0,1,-80,5);頻率響應：產生心電圖信號：代碼：x=-4,-2,0,-4,-6,-4,-2,-4,-6,-6,-4,-4,-6,-6,-2,6,12,8,0,-16,-38,-60,-84,-90,-66,-32,-4,2,-4,8,12,12,10,6,6,6,4,0,0,0,0,0,-2,-4,0,0,0,-2,-2,0,0,-2,-2,-2,-2,0,2,-4,-2,0,-2,-4,-4,2,0,0,-2,-4,-2,0,0,-2,-4,-2

8、,0,0,-4,-4,-2,-2,-4,-6,-6,-4,-4,8,-10,-8,-6,-6,-8,-12,-10,-8,-8,-10,-12,-10,-8,-8,-10,-10,-8,-6,-6,-8,-8,-4,-2,-4,-4,-4,0,0,-2,-4,-2,-2,0,-4;tt=0:119;figure(2);plot(tt,x);y=fft(x);figure(3);plot(abs(y);時域波形：頻域波形：經過IIR濾波和FIR濾波：代碼：N1,wc= buttord(0.2,0.3,1,15);b,a= butter(N1,wc);Y1=filter(b,a,x);figure

9、(4);plot(Y1);Y2=conv(x,h);figure(5);plot(Y2);F1=fft(Y1);figure(6);plot(abs(F1);F2=fft(Y2);figure(7);plot(abs(F2);經過IIR濾波后時域波形：經IIR濾波后的頻域波形：經FIR濾波后的時域波形：經FIR濾波后的頻域波形： FIR濾波器的C語言代碼：RunIntegerFIRFilter(double *FiltCoeff, int NumTaps, double *Signal, double *FilteredSignal, int NumSigPts) int j, k; long

10、 long y; long long One = 1; double TwoToN = One (RegisterSize-1); / Minus one because of the sign bit. long long Reg256; long long IntFirCoeff256; for(j=0; j256; j+)Regj = 0; / Init the delay registers. for(j=0; jNumTaps; j+)IntFirCoeffj = RndWholePart(TwoToN * FirCoeffj); / Round the coefficient to the nearest int for(j=0; j1; k-)Regk-1 = Regk-2; Reg0 = R