第一篇:數字信號處理實驗報告完整版
實驗 1
利用 T DFT 分析信號頻譜
一、實驗目的
1.加深對 DFT 原理的理解。
2.應用 DFT 分析信號的頻譜。
3.深刻理解利用 DFT 分析信號頻譜的原理,分析實現過程中出現的現象及解決方法。
二、實驗設備與環境
計算機、MATLAB 軟件環境 三、實驗基礎理論
T 1.DFT 與 與 T DTFT 的關系
有限長序列 的離散時間傅里葉變換 在頻率區間的 N 個等間隔分布的點 上的 N 個取樣值可以由下式表示:
212 /0()|()()0 1Nj knjNk NkX e x n e X k k N??? ?????? ? ? ? ?? 由上式可知,序列 的 N 點 DFT ,實際上就是 序列的 DTFT 在 N 個等間隔頻率點 上樣本。
2.利用 T DFT 求 求 DTFT
方法 1 1:由恢復出的方法如下:
由圖 2.1 所示流程可知:
101()()()Nj j n kn j nNn n kX e x n e X k W eN? ? ?? ? ?? ? ???? ??? ?? ?? ?? ?? ?? ? ? 由上式可以得到:
IDFT DTFT
第二篇:數字信號處理實驗報告
南京郵電大學
實 驗 報 告
實驗名稱_____熟悉MATLAB環境 ___ 快速傅里葉變換及其應用 ____IIR數字濾波器的設計_ FIR數字濾波器的設計
課程名稱 數字信號處理A
班級學號_______09002111___________ 姓 名 王都超
開課時間 2011/2012學年,第 二 學期
實驗一
熟悉MATLAB環境
一、實驗目的
(1)熟悉MATLAB的主要操作命令。(2)學會簡單的矩陣輸入和數據讀寫。(3)掌握簡單的繪圖命令。
(4)用MATLAB編程并學會創建函數。(5)觀察離散系統的頻率響應。
二、實驗內容
(1)數組的加、減、乘、除和乘方運算。輸入A=[1 2 3 4],B=[3,4,5,6],求
C=A+B,D=A-B,E=A.*B,F=A./B,G=A.^B。并用stem語句畫出A、B、C、D、E、F、G。
D =
-2 E =
F =
0.3333
0.5000
0.6000
0.6667 G =
243
4096(2)用MATLAB實現下列序列: a)x(n)?0.8n 0?n?1
5n=0:1:15;x1=0.8.^n;a=(0.2+3*i)*n;stem(x1)b)x(n)?e(0.2?3j)n 0?n?15
n=0:1:15;x2=exp(a);a=(0.2+3*i)*n;stem(x2)
c)x(n)?3cos(0.125?n?0.2?)?2sin(0.25?n?0.1?)
0?n?15
(4)繪出下列時間函數的圖形,對x軸、y軸以及圖形上方均須加上適當的標注: a)x(t)?sin(2?t)0?t?10s
b)x(t)?cos(100?t)sin(?t)0?t?4s t=0:0.01:4;x=cos(100*pi*t).*sin(pi*t);plot(t,x, 'r-');xlabel('t'),ylabel('x(t)'),title('cos')
(6)給定一因果系統H(z)?(1?頻響應和相頻響應。
2z?1?z?2)/(1?0.67z?1?0.9z?2),求出并繪制H(z)的幅
(7)計算序列{8-2-1 2 3}和序列{2 3-1-3}的離散卷積,并作圖表示卷積結果。
(8)求以下差分方程所描述系統的單位脈沖響應h(n), 0?n?50
y(n)?0.1y(n?1)?0.06y(n?2)?x(n)?2x(n?1)
實驗過程與結果(含實驗程序、運行的數據結果和圖形); clear all;N=50;a=[1-2];b=[1 0.1-0.06];x1=[1 zeros(1,N-1)];n=0:1:N-1;h=filter(a,b,x1);stem(n,h)axis([-1 53-2.5 1.2])
實驗二
快速傅里葉變換及其應用
一、實驗目的
(1)在理論學習的基礎上,通過本實驗,加深對FFT的理解,熟悉MATLAB中的有關函數。(2)應用FFT對典型信號進行頻譜分析。
(3)了解應用FFT進行信號頻譜分析過程中可能出現的問題,以便在實際中正確應用FFT。(4)應用FFT實現序列的線性卷積和相關。
二、實驗內容
實驗中用到的信號序列 a)高斯序列
??(n?p)q?xa(n)??e?0?20?n?15 其他
b)衰減正弦序列
?e?ansin(2?fn)xb(n)??0?0?n?15其他
c)三角波序列 ?n?xc(n)??8?n?0?0?n?34?n?7 其他
d)反三角波序列
?4?n?xd(n)??n?4?0?0?n?34?n?7 其他
(1)觀察高斯序列的時域和幅頻特性,固定信號xa(n)中參數p=8,改變q的值,使q分別等于2,4,8,觀察它們的時域和幅頻特性,了解當q取不同值時,對信號序列的時域幅頻特性的影響;固定q=8,改變p,使p分別等于8,13,14,觀察參數p變化對信號序列的時域及幅頻特性的影響,觀察p等于多少時,會發生明顯的泄漏現象,混疊是否也隨之出現?記錄實驗中觀察到的現象,繪出相應的時域序列和幅頻特性曲線。
(3)觀察三角波和反三角波序列的時域和幅頻特性,用N=8點FFT分析信號序列xc(n)和觀察兩者的序列形狀和頻譜曲線有什么異同?繪出兩序列及其幅頻特性xd(n)的幅頻特性,曲線。
在xc(n)和xd(n)末尾補零,用N=32點FFT分析這兩個信號的幅頻特性,觀察幅頻特性發生了什么變化?兩種情況的FFT頻譜還有相同之處嗎?這些變化說明了什么?
(5)用FFT分別實現xa(n)(p=8,q=2)和xb(n)(a=0.1,f=0.0625)的16點循環卷積和線性卷積。
n=0:15;p=8;q=2;
xa=exp(-(n-p).^2/q);subplot(2,3,1);stem(n,xa,'.');title('xa波形');
Xa=fft(xa,16);subplot(2,3,4);stem(abs(Xa),'.');
title('Xa(k)=FFT[xa(n)]的波形 ');A=1;f=0.0625;a=0.1;
xb=exp(-a*n).*sin(2*pi*f*n);subplot(2,3,2);stem(n,xb,'.');title('xb波形');Xb=fft(xb,16);subplot(2,3,5);stem(abs(Xb),'.');
title('Xb(k)=FFT[xb(n)]的波形 ');
實驗過程與結果(含實驗程序、運行的數據結果和圖形);
實驗三 IIR數字濾波器的設計
一、實驗目的
(1)掌握雙線性變換法及脈沖響應不變法設計IIR數字濾波器的具體設計方法及其原理,熟悉用雙線性變換法及脈沖響應不變法設計低通、高通和帶通IIR數字濾波器的計算機編程。
(2)觀察雙線性變換及脈沖響應不變法設計的濾波器的頻域特性,了解雙線性變換法及脈沖響應不變法的特點。
(3)熟悉巴特沃思濾波器、切比雪夫濾波器和橢圓濾波器的頻率特性。
二、實驗內容(1)P162 例4.4 設采樣周期T=250?s(采樣頻率fs=4kHz),分別用脈沖響應不變法和雙線性變換法設計一個三階巴特沃思低通濾波器,其3dB邊界頻率為fc=1kHz。
脈沖響應不變法: fc=1000;fs=4000;OmegaC=2*pi*fc;[B,A]=butter(3, OmegaC,'s');[num1,den1]=impinvar(B,A,fs);[h1,w]=freqz(num1,den1);f = w/pi*fs/2;plot(f,abs(h1));
雙線性變換法: fc=1000;fs=4000;
OmegaC=2*fs*tan(pi*fc/fs);[B,A]=butter(3, OmegaC,'s');[num2,den2]=bilinear(B,A,fs);[h2,w]=freqz(num2,den2);f = w/pi*fs/2;plot(f,abs(h2));
同一圖中畫兩條曲線: fc=1000;fs=4000;OmegaC=2*pi*fc;[B,A]=butter(3, OmegaC,'s');[num1,den1]=impinvar(B,A,fs);[h1,w]=freqz(num1,den1);f = w/pi*fs/2;
OmegaC=2*fs*tan(pi*fc/fs);[B,A]=butter(3, OmegaC,'s');[num2,den2]=bilinear(B,A,fs);[h2,w]=freqz(num2,den2);f = w/pi*fs/2;plot(f,abs(h1),'r-.');hold on;plot(f,abs(h2),'g-');
(選做)(2)fc=0.2kHz,?=1dB,fr=0.3kHz,At=25dB,T=1ms;分別用脈沖響應不變法及雙線性變換法設計一巴特沃思數字低通濾波器,觀察所設計數字濾波器的幅頻特性曲線,記錄帶寬和衰減量,檢查是否滿足要求。比較這兩種方法的優缺點。
實驗過程與結果(含實驗程序、運行的數據結果和圖形);
實驗四
FIR數字濾波器的設計
一、實驗目的
(1)掌握用窗函數法,頻率采樣法及優化設計法設計FIR濾波器的原理及方法,熟悉相應的計算機編程;
(2)熟悉線性相位FIR濾波器的幅頻特性和相頻特性;
(3)了解各種不同窗函數對濾波器性能的影響。
二、實驗內容
(1)生成一個長度為20的矩形窗,畫出其時域和幅頻特性曲線。n=0:1:19;N=20;win(1:20)=1;[H,w]=freqz(win,1);subplot(2,1,1);stem(n,win)subplot(2,1,2);plot(w,abs(H));
(2)用矩形窗設計一個21階的線性相位低通FIR數字濾波器,截止頻率Wc=0.25π,求出濾波器系數,并繪出濾波器的幅頻特性。修改程序,分別得到階次為N=41,61的濾波器,并顯示其各自的幅頻曲線。
a)在上面所得的幾幅圖中,在截止頻率兩邊可以觀察到幅頻響應的擺動行為。請問波紋的數量與濾波器脈沖響應的長度之間有什么關系?
b)最大波紋的高度與濾波器脈沖響應的長度之間有什么關系?
實驗過程與結果(含實驗程序、運行的數據結果和圖形); 21階的線性相位低通FIR數字濾波器: Wc=0.25*pi;N=21;M=(N-1)/2;
%位移量
for n=0:(N-1)
if(n== fix(M))
%中間的點單獨算
hd(n+1)=Wc/pi;
else
hd(n+1)=sin(Wc*(n-M))/(pi*(n-M));end;end;win=boxcar(N);%%%不同窗函數
h=hd.*win';[H,w]=freqz(h,1);n=0:1:N-1;subplot(3,1,1);stem(n,h)subplot(3,1,2);plot(w,abs(H));subplot(3,1,3);plot(w,angle(H));
41階的線性相位低通FIR數字濾波器: Wc=0.25*pi;N=41;M=(N-1)/2;
%位移量
for n=0:(N-1)
if(n== fix(M))
%中間的點單獨算
hd(n+1)=Wc/pi;
else
hd(n+1)=sin(Wc*(n-M))/(pi*(n-M));end;end;win=boxcar(N);%%%不同窗函數
h=hd.*win';[H,w]=freqz(h,1);n=0:1:N-1;subplot(3,1,1);stem(n,h)subplot(3,1,2);plot(w,abs(H));subplot(3,1,3);plot(w,angle(H));
61階的線性相位低通FIR數字濾波器: Wc=0.25*pi;N=61;M=(N-1)/2;
%位移量
for n=0:(N-1)
if(n== fix(M))
%中間的點單獨算
hd(n+1)=Wc/pi;
else
hd(n+1)=sin(Wc*(n-M))/(pi*(n-M));end;end;win=boxcar(N);%%%不同窗函數
h=hd.*win';[H,w]=freqz(h,1);n=0:1:N-1;subplot(3,1,1);stem(n,h)subplot(3,1,2);plot(w,abs(H));subplot(3,1,3);plot(w,angle(H));
數字信號處理實驗小結及心得體會:
通過這次實驗,我對MATLAB語言有了一定的認識,雖然還不能完全用MATLAB獨立編寫程序,但對這種語言環境有了新的了解。我知道了一般的加減乘除在MATLAB中不同的意義。知道輸入、輸出語句怎么形成。通過快速傅里葉變換及其應用的實驗,加深了我對FFT的理解,還有對各典型信號的頻譜分析,改變參數后時域和幅頻特性的變化。IIR數字濾波器的設計讓我知道了巴特沃思濾波器和切比雪夫濾波器的頻率特性,還有雙線性變換及脈沖響應不變法設計的濾波器的頻率特性。做這個實驗的時候程序有點困難,很多細節問題不能考慮清楚,導致圖形出不來。FIR數字濾波器的設計出來的是三種窗的圖形,通過三種窗的比較,我了解了他們各自的特點,幅頻和相頻特性。我在這次實驗中的收獲很大,接觸了很多新的知識,但在實驗寫程序時,我發現自己還有很多不足。很多程序寫不完全。這是自己今后要加強的地方。
第三篇:數字信號處理實驗報告
JIANGSU
UNIVERSITY OF TECHNOLOGY
數字信號處理實驗報告
學院名稱: 電氣信息工程學院
專 業:
班 級: 姓 名: 學 號: 指導老師: 張維璽(教授)
2013年12月20日
實驗一 離散時間信號的產生
一、實驗目的
數字信號處理系統中的信號都是以離散時間形態存在的,所以對離散時間信號的研究是數字信號的基本所在。而要研究離散時間信號,首先需要產生出各種離散時間信號。使用MATLAB軟件可以很方便地產生各種常見的離散時間信號,而且它還具有強大繪圖功能,便于用戶直觀地處理輸出結果。
通過本實驗,學生將學習如何用MATLAB產生一些常見的離散時間信號,實現信號的卷積運算,并通過MATLAB中的繪圖工具對產生的信號進行觀察,加深對常用離散信號和信號卷積和運算的理解。
二、實驗原理
離散時間信號是指在離散時刻才有定義的信號,簡稱離散信號,或者序列。離散序列通常用x(n)來表示,自變量必須是整數。常見的離散信號如下:(1)單位沖激序列δ(n)
如果δ(n)在時間軸上延遲了k個單位,得到δ(n-k),即長度為N的單位沖激序列δ(n)可以通過下面的MATLAB命令獲得。
n=-(N-1):N-1 x=[zeros(1,N-1)1 zeros(1,N-1)]; stem(n,x)延遲K個采樣點的長度為N的單位沖激序列δ(n-k)(k n=0:N-1 y=[zeros(1,M)1 zeros(1,N-M-1)]; stem(n,y) (2)單位階躍序列u(n) 如果u(n)在時間軸上延遲了k個單位,得到u(n-k),即長度為N的單位階躍序列u(n)可以通過下面的MATLAB命令獲得。 n=-(N-1):N-1 x=[zeros(1,N-1)ones(1,N)]; stem(n,x)延遲的單位階躍序列可以使用類似于單位沖激序列的方法獲得。(3)矩形序列 矩形序列有一個重要的參數,就是序列的寬度N。矩形序列與u(n)之間的關系為矩形序列等= u(n)— u(n-N)。 因此,用MATLAB表示矩形序列可利用上面的單位階躍序列組合而成。(4)正弦序列x(n) 這里,正弦序列的參數都是實數。與連續的正弦信號不同,正弦序列的自變量n必須為整數??梢宰C明,只有當2π/w為有理數時,正弦序列具有周期性。 長度為N的正弦序列x(n)可以通過下面的MATLAB命令獲得。n=0:N-1 x=A*cos(2*pi*f*n/Fs+phase)(5)單邊實指數序列x(n) 長度為N的實指數序列x(n)可以通過下面的MATLAB命令實現。n=0:N-1 x=a.^n stem(n,x)單邊指數序列n的取值范圍為n>=0。當|a|>1時,單邊指數序列發散;當|a|<1時,單邊指數序列收斂。當a>0時,該序列均取正值;當a<0時,序列在正負擺動。 (6)負指數序列x(n) 當a=0時,得到虛指數序列x(n)。 與連續負指數信號一樣,我們將負指數序列實部和虛部的波形分開討論,得到如下結論: 1)當a>0時,負指數序列x(n)的實部和虛部分別是按指數規律增長的正弦振蕩序列; 2)當a<0時,負指數序列x(n)的實部和虛部分別是按指數規律衰減的正弦振蕩序列; 3)當a=0時,負指數序列x(n)即為虛指數序列,其實部和虛部分別是等幅的正弦振蕩序列; 長度為N的實指數序列x(n)可以通過下面的MATLAB命令實現。n=0:N-1 x=exp((a.+j*w)*n)stem(n,real(x))或 stem(n,imag(x)) 三、實驗內容及分析 ?1n?01、編制程序產生單位沖激序列??n???“?并繪出其圖及??n?”學號后兩位0n?0?形。程序:(1)N=4; n=-(N-1):N-1; x=[zeros(1,N-1)1 zeros(1,N-1)];stem(n,x); title('單位沖激序列'); grid on; (2)N=6; M=1;%學號01 n=-(N-1):N-1; y=[zeros(1,N-M+1)1 zeros(1,N-M-1)];stem(n,y); title('單位沖激序列');grid on; 分析:在上圖的基礎上向右平移了1個單位。 ?1n?02、編制程序產生單位階躍序列u?n???、u?n?“學號后兩位”?及 0n?0?u?n??u?n?“學號后兩位”?,并繪出其圖形。程序: 4 (1)N=5; n=-(N-1):N-1; x=[zeros(1,N-1)ones(1,N)];stem(n,x); title('單位階躍序列');grid on; (2)N=6; M=1;%學號01 n=-(N-1):N-1; x=[zeros(1,N-M+1)ones(1,N-M)];stem(n,x); title('單位階躍序列');grid on; 分析:在上圖的基礎上平移了1個單位.(3)N=6; M=1;%學號01 n=-(N-1):N-1; x=[zeros(1,N-1)ones(1,N)];y=[zeros(1,N-M+1)ones(1,N-M)];z=x-y;stem(n,z); title('單位階躍序列');grid on; 2?? 3、編制程序產生正弦序列x?n??cos?2?n?、x?n??cos??n?及 ?學號后兩位?x?n??sin?2n?并繪出其圖形。 程序:(1)N=5; A=1; w=2*pi;phi=0;n=0:0.05:N-1;x=A*cos(w*n+phi);stem(n,x);title('余弦信號');grid on; 分析:該序列具有周期性,且輸出為余弦信號.(2)N=5; A=1; w=2*pi/1;%學號01 phi=0;n=0:0.05:N-1;x=A*cos(w*n+phi);stem(n,x);title('余弦信號');grid on; ; 分析:該序列具有周期性,且輸出為余弦信號.(3)N=5; A=1; w=2*pi;phi=0; n=0:0.05:N-1;x=A*sin(w*n+phi);stem(n,x);title('正弦信號');grid on; 分析:該序列具有周期性,且輸出為正弦信號.4、編制程序產生復正弦序列x?n??e(2?j學號后兩位)n,并繪出其圖形。N=3; n=0:0.2:N-1; w=1;%學號01 x=exp((2+j*w)*n);subplot(2,1,1) stem(n,real(x)),title('實部');grid on;subplot(2,1,2) stem(n,imag(x)),title('虛部');grid on; 5、編制程序產生指數序列x?n??an,并繪出其圖形。其中a=學號后兩位、a=1/“學號后兩位”。 (1)N=10; n=0:N-1; a=1;%學號01 x=a.^n;stem(n,x);title('指數序列');grid on; (2)N=10; n=0:N-1; a=1;%學號01 x=a.^(-n);stem(n,x);title('指數序列');grid on; 實驗三 離散時間信號的頻域分析 一、實驗目的 信號的頻域分析是信號處理中一種有效的工具。在離散信號的頻域分析中,通常將信號表示成單位采樣序列的線性組合,而在頻域中,將信號表示成復變量或的線性組合。通過這樣的表示,可以將時域的離散序列映射到頻域以便于進一步的處理。 在本實驗中,將學習利用MATLAB計算離散時間信號的DTFT和DFT,并加深對其相互關系的理解。 二、實驗原理 (1)DTFT和DFT的定義及其相互關系。 (2)使用到的MATLAB命令有基于DTFT離散時間信號分析函數以及求解序列的DFT函數。 三、實驗內容及分析 (1)編程計算并畫出下面DTFT的實部、虛部、幅度和相位譜。 X(e)?jw0.0518?0.1553e1?1.2828ex(n)?cos?jw?jw?0.1553e?j2w?1.0388e?j2w?0.0518e?j3w?0.3418e?j3w (2)計算32點序列 5?n16,0≦n≦31的32點和64點DFT,分別繪出幅度譜圖形,并繪出該序列的DTFT圖形。 3-1 clear; x=[0.0518,-0.1553,0.1553,0.0518];y=[1,1.2828,1.0388,0.3418];w=[0:500]*pi/500 H=freqz(x,y,w); magX=abs(H);angX=angle(H);realX=real(H);imagX=imag(H);subplot(221);plot(w/pi,magX);grid; xlabel('frequency in pi unit');ylabel('magnitude');title('幅度 part');axis([0 0.9 0 1.1]); subplot(223);plot(w/pi,angX);grid; xlabel('frequency in pi unit');ylabel('radians');title('相位 part');axis([0 1-3.2 3.2]); subplot(222);plot(w/pi,realX);grid; xlabel('frequency in pi unit');ylabel('real part');title('實部 part');axis([0 1-1 1]); subplot(224);plot(w/pi,imagX);grid; xlabel('frequency in pi unit');ylabel('imaginary');title('虛部 part');axis([0 1-1 1.1]); 3-2 N=32;n=0:N-1; xn=cos(5*pi*n/16);k=0:1:N-1;Xk=fft(xn,N);subplot(2,1,1);stem(n,xn);subplot(2,1,2);stem(k,abs(Xk));title('32點');figure N=64;n=0:N-1; xn=cos(5*pi*n/16);k=0:1:N-1;Xk=fft(xn,N);subplot(2,1,1);stem(n,xn);subplot(2,1,2);stem(k,abs(Xk));title('64點'); (1) (2) 實驗四 離散時間LTI系統的Z域分析 一、實驗目的 本實驗通過使用MATLAB函數對離散時間系統的一些特性進行仿真分析,以加深對離散時間系統的零極點、穩定性,頻率響應等概念的理解。學會運用MATLAB分析離散時間系統的系統函數的零極點;學會運用MATLAB分析系統函數的零極點分布與其時域特性的關系;學會運用MATLAB進行離散時間系統的頻率特性分析。 二、實驗原理 離散時間系統的系統函數定義為系統零狀態響應的Z變化與激勵的Z變化之比。 在MATLAB中系統函數的零極點可通過函數roots得到,也可借助函數tf2zp得到,tf2zp的語句格式為 [Z,P,K]=tf2zp(B,A)其中,B與A分別表示H(z)的分子與分母多項式的系數向量。它的作用是將H(z)的有理分式表示式轉換為零極點增益形式。 若要獲得系統函數H(z)的零極點分布圖,可直接應用zplane函數,其語句格式為 Zplane(B,A) 其中,B與A分別表示H(z)的分子和分母多項式的系數向量。它的作用是在z平面上畫出單位圓、零點與極點。 離散系統中z變化建立了時域函數h(n)與z域函數H(z)之間的對應關系。因此,z變化的函數H(z)從形式可以反映h(n)的部分內在性質。可根據系統的傳遞函數H(z)求單位沖激響應h(n)的函數impz、filter等。 利用系統的頻率響應,可以分析系統對各種頻率成分的響應特性,并推出系統的特性(高通、低通、帶通、帶阻等)。 MATLAB提供了求離散時間系統頻響特性的函數freqz,調用freqz的格式主要有兩種。一種形式為 [H,w]= reqz(B,A,N)其中,B與A分別表示H(z)分子和分母多項式的系數向量;N為正整數,默認值為512;返回值w包含[0,π]范圍內的N個頻率等分點;返回值H則是離散時間系統頻率響應在0~π范圍內N個頻率處的值。另一種形式為 [H,w]= freqz(B,A,N,‘whole’) 與第一種方式不同之處在于角頻率的范圍由[0,π]擴展到[0,2π]。 三、實驗內容與結果分析 已知LTI離散時間系統,要求由鍵盤實現系統參數輸入,并繪出幅頻和相頻響應曲線和零極點分布圖,進而分析系統的濾波特性和穩定性。 (一)程序 b=[0.0528,0.797,0.1295,0.1295,0.797,0.0528]; a=[1,-1.8107,2.4947,-1.8801,0.9537,-0.2336];w=[0:20:500]*pi/500; x1=0.0528+0.797*exp(-1*j*w)+0.1295*exp(-2*j*w)+0.1295*exp(-3*j*w)+0.797*exp(-4*j*w)+0.0528*exp(-5*j*w); x2=1-1.8107*exp(-1*j*w)+2.4947*exp(-2*j*w)+1.8801*exp(-3*j*w)+0.9537*exp(-4*j*w)+0.2336*exp(-5*j*w);x22=x2+(x2==0)*eps;x=x1./x22;magx=abs(x); angx=angle(x).*180/pi; subplot(2,2,3);zplane(b,a);title('零極點圖');subplot(2,2,2);stem(w/pi,magx);title('幅度部分');ylabel('振幅');subplot(2,2,4);stem(w/pi,angx); xlabel('以pi為單位的頻率');title('相位部分');ylabel('相位'); (二)波形圖 圖4-1 幅頻、相頻響應曲線、零極點分布圖 實驗六 IIR數字濾波器的設計 一、實驗目的 從理論上講,任何的線性是不變(LTI)離散時間系統都可以看做一個數字濾波器,因此設計數字濾波器實際就是設計離散時間系統。數字濾波器你包括IIR(無限沖激響應)和FIR(有限沖激響應)型,在設計時通常采用不同的方法。 本實驗通過使用MATLAB函數對數字濾波器進行設計和和實現,要求掌握IIR數字巴特沃斯濾波器、數字切比雪夫濾波器的設計原理、設計方法和設計步驟;能根據給定的濾波器指標進行濾波器設計;同時也加深學生對數字濾波器的常用指標和設計過程的理解。 二、實驗原理 在IIR濾波器的設計中,常用的方法是:先根據設計要求尋找一個合適的模擬原型濾波器,然后根據一定的準則將此模擬原型濾波器轉換為數字濾波器。 IIR濾波器的階數就等于所選的模擬原型濾波器的階數,所以其階數確定主要是在模擬原型濾波器中進行的。 IIR數字濾波器的設計方法如下:(1)沖激響應不變法。(2)雙線性變化法。 一般來說,在要求時域沖激響應能模仿模擬濾波器的場合,一般使用沖激響應不變法。沖激響應不變法一個重要特點是頻率坐標的變化是線性的,因此如果模擬濾波器的頻率響應帶限于折疊頻率的話,則通過變換后濾波器的頻率響應可不失真地反映原響應與頻率的關系。 與沖激響應不變法比較,雙線性變化的主要優點是靠頻率的非線性關系得到s平面與z平面的單值一一對應關系,整個值對應于單位圓一周。所以從模擬傳遞函數可直接通過代數置換得到數字濾波器的傳遞函數。 MATLAB提供了一組標準的數字濾波器設計函數,大大簡化了濾波器的設計工程。 (1)butter。 (2)cheby1、cheby2。 三、實驗內容及分析 利用MATLAB編程方法或利用MATLAB中fdatool工具設計不同功能的IIR數字濾波器。 1、基于chebyshev I型模擬濾波器原型使用沖激不變轉換方法設計數字濾波器,要求參數為通帶截止頻率?p?0.4?;通帶最大衰減Ap?1dB;阻帶截止頻率?s?0.4?;阻帶最小衰減As?35dB。 程序: wp=0.2*pi; %通帶邊界頻率 ws=0.4*pi; %阻帶截止頻率 rp=1; %通帶最大衰減 rs=35; %阻帶最小衰減 Fs=1000; %?ùéè3é?ù??3?1000hz [N,Wn]=cheb1ord(wp,ws,rp,rs,'s'); [Z,P,K]=cheby1(N,rp,Wn,'s');[H,W]=zp2tf(Z,P,K); figure(1);freqs(H,W);[P,Q]=freqs(H,W);figure(2);plot(Q*Fs/(2*pi),abs(P));grid on; xlabel('頻率/Hz');ylabel('幅度'); 2、基于Butterworth型模擬濾波器原型使用雙線性變換方法設計數字濾波器的,要求參數為截止頻率?p?0.4?;通帶最大衰減Ap?1dB;阻帶截止頻率?s?0.25?;阻帶最小衰減AS?40dB。程序: wp=0.4*pi;ws=0.25*pi;rp=1;rs=40;fs=500;ts=1/fs;wp1=wp*ts;ws1=ws*ts; wp2=2*fs*tan(wp1/2);ws2=2*fs*tan(ws1/2); [N,Wn]=buttord(wp2,ws2,rp,rs,'s');[Z,P,K]=buttap(N);[Bap,Aap]=zp2tf(Z,P,K);[b,a]=lp2lp(Bap,Aap,Wn);[bz,az]=bilinear(b,a,fs);[H,W]=freqz(bz,az);subplot(2,1,1);plot(W/pi,abs(H));grid on;xlabel('頻率')ylabel('幅度')subplot(2,1,2); plot(W/pi,20*log10(abs(H)));grid on;xlabel('頻率');ylabel('幅度(dB)'); 實驗七 FIR數字濾波器的設計 一、實驗目的 掌握用窗函數設計FIR數字濾波的原理及其設計步驟;熟悉線性相位數字濾波器的特性。學習編寫數字濾波器的設計程序的方法,并能進行正確編程;根據給定的濾波器指標,給出設計步驟。 二、實驗原理 如果系統的沖激響應h(n)為已知,則系統的輸入輸出關系為 y(n)=x(n)*h(n) 對于低通濾波器,只要設計出低通濾波器的沖激響應函數,就可以由式得到系統的輸出了。 但是將h(n)作為濾波器的脈沖響應有兩個問題:一是它是無限長的;二是它是非因果的。對此,采取兩項措施:一是將h(n)截短;二是將其右移。 設計時,要根據阻帶的最小衰減和過渡帶寬度來選擇恰當的窗函數類型和窗口長度N。常用的窗函數有矩形窗、海明窗和布萊克曼窗等。 窗函數設計FIR濾波器步驟如下: (1)給定理想頻率響應的幅頻特性和相頻特性; (2)求理想單位脈沖響應,在實際計算中,可對理想頻率響應采樣。(3)根據過渡帶寬度和阻帶最小衰減,確定窗函數類型和窗口長度N;(4)求FIR濾波器單位脈沖響應; (5)分析幅頻特性,若不滿足要求,可適當改變窗函數形式或長度N,重復上述設計過程,以得到滿意的結果。 三、實驗內容及分析 1、分別用海明窗和布萊克曼窗設計一個48階的FIR帶通濾波器,通帶為Wn??0.450.55?。程序1:海明窗設計 N=48; Window=hamming(N+1);w1=0.45;w2=0.55;ws=[w1,w2]; b=fir1(N,ws/pi,Window);freqz(b,1,512);title('海明窗');grid on; 程序2:萊克曼窗設計 N=48; Window=blackman(N+1);w1=0.45;w2=0.55;ws=[w1,w2]; b=fir1(N,ws/pi,Window);freqz(b,1,512);title('布萊克曼窗');grid on; 2、用矩形窗設計一個線性相位高通濾波器。其中He??jw?e?j????????00.3????? 0???0.3?程序: N=9; alpha=(N-1)/2;Wc=0.7*pi;n=(0:8);i=n-alpha;i=i+(i==0)*eps; h=(-1).^n.*sin((i).*Wc)./((i).*pi);%矩形窗函數設計的系統脈沖響應 w=(0:1:500)*2*pi/500; H=h*exp(-j*n'*w);%矩形窗函數設計的頻響 magH=abs(H);% 矩形窗函數設計的振幅 subplot(211);stem(n,h); axis([0,8,-0.4,0.4]);title('矩形窗設計h(n)');line([0,10],[0,0]);xlabel('n');ylabel('h');subplot(212);plot(w/pi,magH); xlabel('以pi為單位的頻率');ylabel('H振幅');axis([0,2,0,1.7]);title('矩形窗設計振幅譜'); 實驗心得體會: 這次實驗使我進一步加深了對MATLAB軟件的使用。從上次的信號系統實驗的初步使用到這一次的深入了解,有了更深刻的認識。對這種語言環境也有了新的了解。 在實驗的過程中,我對數字濾波器的整個過程有了很好的理解和掌握。IIR數字濾波器的設計讓我知道了巴特沃思濾波器和切比雪夫濾波器的頻率特性,還有雙線性變換及脈沖響應不變法設計的濾波器的頻率特性。做這兩個實驗的時候程序有點困難,但經過細心的改寫圖形最終出來了。FIR數字濾波器的設計出來的是兩種窗的圖形,通過兩種窗的比較,我了解了他們各自的特點,幅頻和相頻特性。 最后,感謝張老師對我的諄諄教導! 數字信號處理 實驗報告 實驗一 序列的傅立葉變換 一、實驗目的 1.進一步加深理解DFS,DFT算法的原理;2.研究補零問題;3.快速傅立葉變換(FFT)的應用。 二、實驗步驟 1.復習DFS和DFT的定義,性質和應用; 2熟悉MATLAB語言的命令窗口、編程窗口和圖形窗口的使用;3利用提供的程序例子編寫實驗用程序;4.按實驗內容上機實驗,并進行實驗結果分析;5.寫出完整的實驗報告,并將程序附在后面。 三、實驗內容 1.周期方波序列的頻譜 試畫出下面四種情況下的的幅度頻譜, 并分析補零后,對信號頻譜的影響。 x(n)?cos(0.48?n)?cos(0.52?n)2.有限長序列x(n)的DFT(1)取x(n)(n=0:10)時,畫出x(n)的頻譜X(k)的幅度;(2)將(1)中的x(n)以補零的方式,使x(n)加長到(n:0~100)時,畫出x(n)的頻譜X(k)的幅度; (3)取x(n)(n:0~100)時,畫出x(n)的頻譜X(k)的幅度。利用FFT進行譜分析x(t)?2sin(4?t)?5cos(8?t)3.已知:模擬信號 以t=0.01n(n=0:N-1)進行采樣,求N點DFT的幅值譜。請分別畫出N=45;N=50;N=55;N=60時的幅值曲線。 四、實驗數據分析 1.周期方波序列的頻譜分析 首先定義一個功能函數dfs function[Xk]=dfs(xn,N)n=[0:1:N-1];k=[0:1:N-1];WN=exp(-j*2*pi/N);nk=n'*k;WNnk=WN.^nk;Xk=xn*WNnk;(1)L=5,N=20;%題1.(1)L=5;N=20;%對于(2),(3),(4)問,只要修改L,N的數值就好。n=1:N;xn=[ones(1,L),zeros(1,N-L)];Xk=dfs(xn,N);magXk=abs([Xk(N/2+1:N)Xk(1:N/2+1)]);k=[-N/2:N/2];figure(1)subplot(2,1,1);stem(n,xn);xlabel('n');ylabel('xtide(n)');title('DFS of SQ.wave:L=5,N=20');subplot(2,1,2);stem(k,magXk);axis([-N/2,N/2,0,16]);xlabel('k');ylabel('Xtide(k)'); (2)L=5,N=40; (3).L=5,N=60 (4)L=7,N=60; 結果分析:雖然周期序列不存在FT,但是一個周期序列可以利用其DFS系數X(k)表示它的頻譜分布規律,從以上各頻譜圖可以看出,隨著補零點數的增加,周期序列的諧波次數越來越多,其頻譜的包絡線越來越平滑連續,更能反映幅度值隨時間的變化。 2.有限長序列的DFT(1) %題2-(1)n=0:10;xn=cos(0.48*pi*n)+cos(0.52*pi*n);N=11;Xk=fft(xn,N); %序列x(n)的N點DFT k=0:N-1;wk=2*k/N; subplot(1,1,1);stem(wk,abs(Xk),'.');title('頻譜X(K)的幅度');xlabel('ω/π');ylabel('幅度'); (2)%題2-2 M=10;N=100;n=1:M;xn=cos(0.48*pi*n)+cos(0.52*pi*n);n1=[0:1:N-1];y1=[xn(1:1:M),zeros(1,N-M)];figure(1)subplot(2,1,1);stem(n1,y1);xlabel('n');ylabel('x(n)');title('序列x(n),0<=n<=100');axis([0,N,-2.5,2.5]);Y1=fft(y1);magY1=abs(Y1(1:1:N/2+1));k1=0:1:N/2;w1=2*pi/N*k1;subplot(2,1,2);title('x(n)的幅頻特性曲線');stem(w1/pi,magY1);axis([0,1,0,60]);xlabel('omega/pi');ylabel('|X(K)|'); (3) %Example2-3 M=10;N=100;n=0:M;xn=cos(0.48*pi*n)+cos(0.52*pi*n);n1=[0:1:N-1];y1=[xn(1:1:M),zeros(1,N-M)];figure(1)subplot(2,2,1);stem(n1,y1);xlabel('n');ylabel('x(n)');title('序列x(n),0<=n<=100');axis([0,N,-2.5,2.5]);YK=fft(y1);Y=abs(Y1(1:1:N/2+1));k1=0:1:N/2;w1=2*pi/N*k1;subplot(2,2,3);stem(w1/pi,Y);title('x(n)的幅頻特性曲線');axis([0,1,0,60]);xlabel('omega/pi');ylabel('|X(K)|');subplot(2,2,4);plot(angle(Y1));title('x(n)的相頻特性曲線');xlabel('omega/pi');ylabel('phi(omega)'); 結果分析:由上述仿真圖可得,隨著n取值范圍的增大,其頻譜在[0,2π]上的采樣間隔越來越小,采樣點越來越多。采樣點越多,其DFS頻譜越接近FT的頻譜。其相頻特性曲線呈現周期性變化。 3.問題三 %題3 %N=45 figure(1)subplot(2,2,1)N=45;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=fft(x,N);plot(q,abs(y))stem(q,abs(y))title('DFT N=45')%N=50 subplot(2,2,2)N=50;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=fft(x,N);plot(q,abs(y))stem(q,abs(y))title('DFT N=50')%N=55 subplot(2,2,3)N=55;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=fft(x,N);plot(q,abs(y))stem(q,abs(y))title('DFT N=55')%N=60 subplot(2,2,4)N=60;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=fft(x,N);plot(q,abs(y))stem(q,abs(y))title('DFT N=60') 結果分析:由上述仿真圖可得,隨著N取值的增大,其頻譜在[0,2π]上的采樣間隔越來越小,采樣點越來越多。采樣點越多,其DFS頻譜越接近FT的頻譜,幅值曲線越來越清晰,更能準確反應幅值隨時間的變化規律。 五、心得體會 對于周期序列的離散傅里葉變換,通過matlab的模擬,可以更好的了解掌握序列采樣間隔對其采樣頻譜的影響。在實驗過程中,學習如何使用matlab程序語言解決問題,是很有價值的。讓我對matlab的使用更為熟練。 實驗二 用雙線性變換法設計IIR數字濾波器 一、實驗目的 1.熟悉用雙線性變換法設計IIR數字濾波器的原理與方法; 2. 掌握數字濾波器的計算機仿真方法; 3、通過觀察對實際心電圖的濾波作用,獲得數字濾波器的感性知識。 二、實驗內容 1.用雙線性變換法設計一個巴特沃斯低通IIR濾波器,設計指標參數為:在通帶內頻率低于0.2π時,最大衰減小于1dB;在阻帶內[0.3π,π]頻率區間上,最小衰減大于15dB 2.以0.2π為采樣間隔,打印出數字濾波器在頻率區 間[0, 0.2π]上的幅值響應曲線。 3.用所設計的濾波器對實際的心電圖信號采樣序列 x(n)=[-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];)進行仿真濾波處理,并分別打印出濾波前后的心電圖信號波形圖,觀察總結濾波作用與效果。三.實驗步驟 1.復習有關巴特沃斯模擬濾波器設計和雙線性變換法設計IIR濾波器的內容 2.參考例子程序用MATLAB語言編寫仿真實驗用程序; 3.在通用計算機上運行仿真程序 4.寫出完整的實驗報告并回答思考題。四.實驗數據分析 1.巴特沃斯低通IIR濾波器的設計 由題可得,數字低通技術指標為 wp=0.2πrad, αp=1dB,ws=0.3πrad, αs=15dB 若T=1s,預畸變校正計算相應的模擬低通的技術指標為 Ωp=1dB,Ωs=15dB 通過計算可得階數N=5.3056,向上取整的N=6.Ωc=0.7663rad/s,這樣保證阻帶技術指標滿足要求,通帶指標有富余。 以下是通過matlab實現的巴特沃斯低通IIR濾波器的設計 %1 T=1;Fs=1/T;wpz=0.2;wsz=0.3;wp=2*tan(wpz*pi/2);ws=2*tan(wsz*pi/2);rp=1;rs=15;%預畸變校正轉換指標 [N,wc]=buttord(wp,ws,rp,rs,'s');%設計過渡模擬濾波器 [B,A]=butter(N,wc,'s');[Bz,Az]=bilinear(B,A,Fs);%用雙線性變換法轉換成數字濾波器 fk=0:1/512:1;wk=2*pi*fk;Hk=freqs(B,A,wk);figure(1);subplot(2,1,1);plot(fk,20*log10(abs(Hk)));grid on;title('模擬濾波器幅值響應曲線');xlabel('omega/pi');ylabel('幅度(dB)');axis([0,1,-100,5]);[Nd,wdc]=buttord(wpz,wsz,rp,rs);%調用buttord和butter直接設計數字濾波器 [Bz,Az]=butter(N,wdc);wk=0:pi/512:pi;Hz=freqz(Bz,Az,wk);subplot(2,1,2);plot(wk/pi,20*log10(abs(Hz)));grid on;title('數字濾波器幅值響應曲線');xlabel('omega/pi');ylabel('幅度(dB)');axis([0,1,-100,5]); 2.以0.2π為采樣間隔,打印出數字濾波器在頻率區 間[0, 0.2π]上的幅值響應曲線。 在第一問的基礎上加上下面的程序 %2 figure(2);freqz(Bz,Az,[0:0.02*pi:0.2*pi]) 3.在第一問程序的基礎上加上如下程序即可 figure(3);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];subplot(2,2,1);n=0:55;stem(n,x,'.');title('x(n)的脈沖響應');xlabel('n');ylabel('x(n)');A=0.09036;b1=[A,2*A,A];a1=[1,-1.2686,0.7051];h1=filter(b1,a1,x);[H1,w]=freqz(b1,a1,100);b2=[A,2*A,A];a2=[1,-1.0106,0.3583];h2=filter(b2,a2,h1);[H2,w]=freqz(b2,a2,100);b3=[A,2*A,A];a3=[1,-0.9044,0.2155];h3=filter(b3,a3,h2);[H3,w]=freqz(b3,a3,100);subplot(2,2,2);stem(n,h3,'.');xlabel('n');ylabel('y(n)');title('通過濾波器H1(z),H2(z),H3(z)后的y3(n)函數');subplot(2,2,3);H4=H1.*(H2);H=H4.*(H3);mag=abs(H);db=20*log10((mag+eps)/max(mag));plot(w/pi,db);xlabel('ω/π');ylabel('20log[Ha3(ejw)]');title('通過濾波器H1(z),H2(z),H3(z)后的對數頻率響應20log[Ha3(ejw)]函數');grid;figure(4);N=1024;n=0:N/2-1;Xk=fft(x,N);AXk=abs(Xk(1:N/2));f=(0:N/2-1)*Fs/N;f=f/Fs;subplot(211);plot(f,AXk);title('x(n)的頻譜');xlabel('f');ylabel('| X(k)|');axis([0,0.5,0,400]); 五、思考題 用雙線性變換設計數字濾波器的過程中,下面變換公式的T值的取值,對設計結果是否有影響?為什么? 21?z?1s? T1?z?1 雖然采用雙線性變換法設計數字濾波器不會產生混疊現象,T得取值可以任選。雙線性變換法具有非線性,T小一些,非線性的影響也就少一些。 六、心得體會 通過這次實驗,我學會了如何使用matlab語言來實現IIR低通濾波器的設計,同時對于其數字低通技術指標的就算有了進一步的認識體會。 在做實驗的過程中,面對種種困難,但是卻在解決困難的過程中收獲了很多。是一個很有意義的經歷。 實驗三 用窗函數法設計FIR數字濾波器 一.實驗目的 1.掌握用窗函數法設計FIR數字濾波器的原理和方法。2.熟悉線性相位FIR數字濾波器特性。3.了解各種窗函數對濾波特性的影響。二.實驗原理 如果所希望的濾波器的理想頻率響應函數為 Hd(e jω),則其對應的單位脈沖響應為 用窗函數w(n)將hd(n)截斷,并進行加權處理,得到: h(n)就作為實際設計的FIR數字濾波器的單位脈沖響應序列,其頻率響應函數H()為 如果要求線性相位特性,則h(n)還必須滿足: 根據上式中的正、負號和長度N的奇偶性又將線性相位FIR濾波器分成四類。要根據所設計的濾波特性正確選擇其中一類。例如,要設計線性相位低通特性,可選擇h(n)=h(N-1-n)一類,而不能選h(n)=-h(N-1-n)一類。三.實驗內容 1.用MATLAB產生各種窗函數 %窗函數 subplot(4,2,1)m=200;a=boxcar(m);%矩形窗 m=1:200;plot(m,a)title('矩形窗');subplot(4,2,2)m=200;b=bartlett(m)%三角窗 m=1:200;plot(m,b)title(' 三角窗');subplot(4,2,3)m=200;c=hanning(m);%漢寧窗 m=1:200;plot(m,c)title('漢寧窗');subplot(4,2,4)m=200;d=hamming(m);%哈明窗 m=1:200;title(' 哈明窗');plot(m,d)subplot(4,2,5)m=200;e=blackman(m);%布萊克曼窗 m=1:200;plot(m,e)title('布萊克曼窗');subplot(4,2,6)m=200;f=kaiser(m,7.865);% 凱塞窗 m=1:200;plot(m,f)title(' 凱塞窗');subplot(4,2,7)plot(m,a,'r*',m,b,'g+',m,c,'y*',m,d,'b.',m,e,'y.',m,f,'k.')title(' 各種窗函數'); 2.利用窗函數設計FIR濾波器 設計具有下列指標?p=0.25?,Rp=0.25dB,?s=0.3?,Rp=50dB的低通數字濾波器。 由已知條件Rp=50dB,通過計算可知應該選擇哈明窗,哈明窗的Rps=53dB,選擇窗函數時應該選擇Rps>=50的,因而選擇離50dB的哈明窗。 通過Bt=6.6π/N,Bt=ws-wp得到階數N。 以下是利用matlab來實現FIR濾波器的設計。 %窗函數法設計FIR數字濾波器 wp=pi*0.25;ws=pi*0.3;%Rs=50dB,所以選擇哈明窗 DB=ws-wp;%計算過渡帶寬度 N=ceil(6.6*pi/DB);%計算哈明窗所需長度N wc=(wp+ws)/2/pi;%計算理想低通濾波器通帶截止頻率(關于π歸一化)hn=fir1(N,wc);%調用fir1計算低通數字濾波器 %以下是繪圖部分 figure(2)M=1024;hk=fft(hn,M);n=0:N;subplot(1,2,1);stem(n,hn,'.');%繪制序列h(n)xlabel('n');ylabel('h(n)');title('數字濾波器h(n)');k=1:M/2;w=2*(0:M/2-1)/M;subplot(1,2,2);plot(w,20*log10(abs(hk(k))));%繪制損耗函數曲線 axis([0,1,-80,5]);xlabel('ω/π');ylabel('20lg|Hg(ω)|');title('h(n)的損耗函數曲線');grid on; 結果分析:根據損耗函數可得,利用哈明窗實現了題目所要求的FIR低通濾波器。滿足?p=0.25?,Rp=0.25dB,?s=0.3?,Rp=50dB。因而選擇哈明窗可以很好的滿足所需要求。 四、心得體會 本次實驗是使用窗函數來設計FIR濾波器,在MATLAB中只需要將給定的參數輸入到函數中即可馬上得到處結果。通過本次實驗,我對使用MATLAB快速設計濾波器的流程更為熟練,同時,也由衷地感嘆MATLAB的功能強大之處,它讓我們在設計時能節約大量的時間。 根據輸入的參數以及結果,使我對課本上的設計濾波器的知識更為了解。更加鞏固了理論知識。 懷化學院數學系實驗報告 實驗項目名稱:IIR數字濾波器的設計(1) 指 導老 師: 歐衛華 學 姓 實驗項目制定人:實驗項目審批人: 年月日 一、實驗目的掌脈沖相應不變法設計IIR-Butterworth數字濾波器的具體設計方法及原理。 二、實驗原理與方法 1.確定數字濾波器的性能指標:通帶臨界頻率fp、阻帶臨界頻率fs;通帶內的最大衰減Ap;阻帶內的最小衰減As;采樣周期T; 2.確定相應的數字角頻率,ωp=2πfp;ωr=2πfr; 3.根據Ωp和Ωs計算模擬低通原型濾波器的階數N,并求得低通原型的傳遞函 數Ha(s); 4.用上面的脈沖響應不變法公式代入Ha(s),求出所設計的傳遞函數H(z); 5.分析濾波器特性,檢查其是否滿足指標要求。 三、實驗內容及步驟 沖激響應不變法設計數字Butterworth低通濾波器 (1)、模擬濾波器的最小階數[N,wn]=buttord(wp,ws,rp,rs,'s'); (2)、設計模擬低通濾波器原型,[z,p,k]=buttap(N); (3)、將零極點形式轉換為傳遞函數形式,[Bap,Aap]=zp2tf(z,p,k); (4)、進行頻率變換,[b,a]=lp2lp(Bap,Aap,wn); (5)用脈沖相應不變法得到數字濾波器的系統函數[bz,az]=impinvar(b,a,fs); 四、實驗范例 用脈沖相應不變法設計一個Butterworth低通數字濾波器,使其特征逼近一個低通Butterworth模擬濾波器的下列性能指標,通帶截止頻率Wp=2*pi*2000rad/s,通帶波紋Rp小于3dB,阻帶邊界頻率為Ws=2*pi*3000rad/s阻帶衰減大于15dB,采樣頻率Fs=10000;z,假設一個信號x(t)=sin(2*pi*f1*t) +0.5*cos(2*pi*f2*t),其中f1=1000Hz,f2=4000Hz,試將原信號與通過該濾波器的輸出信號進行比較。 wp=2000*2*pi;%濾波器截止頻率 ws=3000*2*pi; rp=3;rs=15;%通帶波紋和阻帶衰減 fs=10000;%采樣頻率 Nn=128; [N,wn]=buttord(wp,ws,rp,rs,'s');%模擬濾波器的最小階數 [z,p,k]=buttap(N);%設計模擬低通濾波器原型 [Bap,Aap]=zp2tf(z,p,k);%將零極點形式轉換為傳遞函數形式 [b,a]=lp2lp(Bap,Aap,wn);%進行頻率變換 [bz,az]=impinvar(b,a,fs);%應用脈沖相應不變法得到數字濾波器的系統函數 figure(1); [h,f]=freqz(bz,az,Nn,fs);%畫出數字濾波器的幅頻特性和相頻特性 subplot(2,1,1),plot(f,20*log10(abs(h))); xlabel('頻率/Hz');ylabel('振幅/dB');grid on; subplot(2,1,2),plot(f,180/pi*unwrap(angle(h))); xlabel('頻率/Hz');ylabel('振幅/^o');grid on; figure(2); f1=1000;f2=4000;%輸入信號的頻率 N=100;%數據長度 dt=1/fs;n=0:N-1;t=n*dt;%采樣間隔和時間序列 x=sin(2*pi*f1*t)+0.5*cos(2*pi*f2*t);%濾波器輸入信號 subplot(2,1,1),plot(t,x),title('輸入信號')%畫出輸入信號 %y=filtfilt(bz,az,x); y1=filter(bz,az,x);%用上面設計的濾波器對輸入信號濾波 subplot(2,1,2),plot(t,y1,'r-'),title('輸出信號'),xlabel('時間/s');legend('filter') 五、實驗習題 用脈沖相應不變法設計一個Butterworth低通數字濾波器,通帶頻率為0= 六,實驗結果第四篇:數字信號處理實驗報告
第五篇:六 數字信號處理實驗報告--IIR數字濾波器設計