數(shù)字信號(hào)處理實(shí)驗(yàn)指導(dǎo)書

1、實(shí)驗(yàn)1 時(shí)域離散信號(hào)的產(chǎn)生一、實(shí)驗(yàn)?zāi)康膌 學(xué)會(huì)運(yùn)用MATLAB產(chǎn)生常用離散時(shí)間信號(hào)。二、實(shí)驗(yàn)涉及的matlab子函數(shù)1、square功能：產(chǎn)生矩形波調(diào)用格式：x=square(t)；類似于sin（t），產(chǎn)生周期為2*pi，幅值為+1的方波。x=square(t，duty)；產(chǎn)生制定周期的矩形波，其中duty用于指定脈沖寬度與整個(gè)周期的比例。2、rand功能：產(chǎn)生rand隨機(jī)信號(hào)。調(diào)用格式：x=rand（n，m）；用于產(chǎn)生一組具有n行m列的隨機(jī)信號(hào)。三、實(shí)驗(yàn)原理在時(shí)間軸的離散點(diǎn)上取值的信號(hào)，稱為離散時(shí)間信號(hào)。通常，離散時(shí)間信號(hào)用x（n）表示，其幅度可以在某一范圍內(nèi)連續(xù)取值。由于信號(hào)處理所用的設(shè)備

2、主要是計(jì)算機(jī)或?qū)Ｓ玫男盘?hào)處理芯片，均以有限的位數(shù)來表示信號(hào)的幅度，因此，信號(hào)的幅度也必須“量化”，即取離散值。我們把時(shí)間和幅度上均取離散值的信號(hào)稱為時(shí)域離散信號(hào)或數(shù)字信號(hào)。在MATLAB中，時(shí)域離散信號(hào)可以通過編寫程序直接生成，也可以通過對(duì)連續(xù)信號(hào)的等間隔抽樣獲得。 下面介紹常用的時(shí)域離散信號(hào)及其程序。1、單位抽樣序列MATLAB源程序?yàn)?) function x,n = impuls (n0,n1,n2)% Generates x(n) = delta(n-n0); n=n0 處建立一個(gè)單位抽樣序列% x,n = impuls (n0,n1,n2)if (n0 < n1) | (n0

3、> n2) | (n1 > n2)error('arguments must satisfy n1 <= n0 <= n2')endn = n1:n2;x = zeros(1,(n0-n1), 1, zeros(1,(n2-n0);將上述文件存為：impuls.m，在命令窗口輸入n0=0,n1=-10,n2=11;x,n=impuls (n0,n1,n2); stem(n,x,filled)2) n1=-5;n2=5;n0=0;n=n1:n2;x=n=n0; stem(n,x,'filled','k');axis(n1,n

4、2,1.1*min(x),1.1*max(x);title('單位脈沖序列');xlabel('時(shí)間(n)');ylabel('幅度x(n)');3) n1=-5;n2=5;k=0; n=n1:n2; nt=length(n); %求n點(diǎn)的個(gè)數(shù)nk=abs(k-n1)+1; %確定k在n序列中的位置x=zeros(1,nt); %對(duì)所有樣點(diǎn)置0x(nk)=1; %對(duì)抽樣點(diǎn)置1stem(n,x,'filled','k');axis(n1,n2,0,1.1*max(x);title('單位脈沖序列');

5、xlabel('時(shí)間(n)');Ylabel('幅度x(n)');2、單位階躍序列MATLAB源程序?yàn)椋?) n1=-2;n2=8;n0=0; n=n1:n2; %生成離散信號(hào)的時(shí)間序列x=n>=n0; %生成離散信號(hào)x(n)stem(n, x,'filled','k'); %繪制脈沖桿圖，且圓點(diǎn)處用實(shí)芯圓表示axis(n1,n2,0,1.1*max(x);title('單位階躍序列');xlabel('時(shí)間(n)');Ylabel('幅度x(n)');2) n1=-2;n2=

6、8;k=0; n=n1:n2; nt=length(n); %求n點(diǎn)的個(gè)數(shù)nk=abs(k-n1)+1; %確定k在n序列中的位置x=zeros(1,nk-1),ones(1,nt-nk+1); %對(duì)所有樣點(diǎn)置0stem(n,x,'filled','k');axis(n1,n2,0,1.1*max(x);title('單位階躍序列');xlabel('時(shí)間(n)');ylabel('幅度x(n)');3、正弦序列x（n）=Um例、已知一時(shí)域周期性正弦信號(hào)的頻率為1HZ，振幅幅度為1V，在窗口中顯示兩個(gè)周期的信號(hào)波形

7、，并對(duì)該信號(hào)的一個(gè)周期進(jìn)行32點(diǎn)采樣獲得離散信號(hào)。顯示連續(xù)信號(hào)和采樣獲得離散信號(hào)波形。MATLAB源程序?yàn)椋篺=1;Um=1;nt=2; %輸入信號(hào)頻率、振幅和顯示周期數(shù)N=32; T=1/f; %N為采樣點(diǎn)數(shù)，T為窗口顯示時(shí)間 dt=T/N; %采樣時(shí)間間隔n=0:nt*N-1; %離散信號(hào)的時(shí)間序列tn=n*dt; %時(shí)間序列樣點(diǎn)在時(shí)間軸上的位置x=Um*sin(2*f*pi*tn);subplot(2,1,1);plot(tn,x); %顯示原信號(hào)axis(0,nt*T,1.1*min(x) 1.1*max(x);ylabel('x(t)');subplot(2,1,2)

8、;stem(tn,x); %顯示經(jīng)采樣的信號(hào)axis(0,nt*T,1.1*min(x),1.1*max(x);ylabel('x(n)');4、矩形序列將square表示式中的t換成n，且n取整數(shù)，則可以獲得矩形序列。例、一個(gè)連續(xù)的周期性矩形信號(hào)頻率為5kHZ，信號(hào)幅度為0-2V之間，脈沖寬度與周期的比例為1:4,且要求在窗口上顯示其2個(gè)周期的信號(hào)波形，并對(duì)信號(hào)的一個(gè)周期進(jìn)行16點(diǎn)采樣來獲得離散信號(hào)，顯示原連續(xù)信號(hào)與采樣獲得的離散信號(hào)。MATLAB源程序?yàn)椋篺=5000;nt=2; N=16;T=1/f; dt=T/N; n=0:nt*N-1;tn=n*dt;x=square

9、(2*f*pi*tn,25)+1; %產(chǎn)生時(shí)域信號(hào),且幅度在02V之間 subplot(2,1,1);stairs(tn,x,'k'); axis(0 nt*T 1.1*min(x) 1.1*max(x);ylabel('x(t)');subplot(2,1,2);stem(tn,x,'filled','k'); axis(0 nt*T 1.1*min(x) 1.1*max(x);ylabel('x(n)');注意：直接用square子函數(shù)產(chǎn)生的信號(hào)波形，其幅度為-11之間。5、rand函數(shù)在實(shí)際系統(tǒng)的研究和處理中

10、，常常需要產(chǎn)生隨機(jī)信號(hào)，MATLAB提供的rand函數(shù)可以為我們生成隨機(jī)信號(hào)。例、生成一組41點(diǎn)構(gòu)成的連續(xù)隨機(jī)信號(hào)和與之相應(yīng)的隨機(jī)序列。MATLAB源程序?yàn)椋簍n=0:40;N=length(tn);x=rand(1,N);subplot(1,2,1),plot(tn,x,'k');subplot(1,2,2),stem(tn,x,'filled','k');四、實(shí)驗(yàn)任務(wù)1、產(chǎn)生離散序列：（1）f（n）= （顯示-3<n<4）（2）f（n）=u（n-1） （顯示-5<n<5）2、一個(gè)連續(xù)的周期性正弦信號(hào)頻率為50HZ，信號(hào)

11、幅度在02V之間，在窗口上顯示2個(gè)周期的信號(hào)波形，對(duì)信號(hào)的一個(gè)周期進(jìn)行16點(diǎn)采樣來獲得離散信號(hào)，顯示原連續(xù)信號(hào)和采樣獲得的離散信號(hào)波形。3、一個(gè)連續(xù)的周期性方波信號(hào)頻率為200HZ,信號(hào)幅度在02V之間，在窗口上顯示兩個(gè)周期的信號(hào)波形，用4Khz的頻率對(duì)連續(xù)信號(hào)進(jìn)行采樣，顯示原連續(xù)信號(hào)和采樣獲得的離散信號(hào)波形。實(shí)驗(yàn)2 離散序列的基本運(yùn)算一、實(shí)驗(yàn)?zāi)康膌 學(xué)會(huì)運(yùn)用MATLAB進(jìn)行離散序列的運(yùn)算，并掌握程序的編寫方法。二、實(shí)驗(yàn)涉及的matlab子函數(shù)1、find功能：尋找非零元素的索引號(hào)。調(diào)用格式：find（n>=min(n1)）&(n<=max(n1)）；在符合關(guān)系運(yùn)算條件的范

12、圍內(nèi)尋找非零元素的索引號(hào)。2、fliplr功能：對(duì)矩陣行元素進(jìn)行左右翻轉(zhuǎn)。調(diào)用格式：x1= fliplr（x）；將x的行元素左右翻轉(zhuǎn)，賦給變量x1。三、實(shí)驗(yàn)原理離散序列的時(shí)域運(yùn)算包括信號(hào)的相加、相乘，信號(hào)的時(shí)域變換包括信號(hào)的移位、反折、尺度變換等。在MATLAB中，離散序列的相加、相乘等運(yùn)算是兩個(gè)向量之間的運(yùn)算，因此參加運(yùn)算的兩個(gè)序列向量必須具有相同的維數(shù)，否則應(yīng)進(jìn)行相應(yīng)的處理。下面介紹各種離散序列的時(shí)域運(yùn)算和時(shí)域變換的性質(zhì)。1、 序列的移位x1（n）=x（n-m）例、x0=u（n），（顯示-10<n<10）x1=u（n+6）（顯示-10<n<10）,x2= u（n-4

13、）（顯示-10<n<10）MATLAB源程序?yàn)椋簄1=-10;n2=10;k0=0;k1=-6;k2=4; n=n1:n2; %生成離散信號(hào)的時(shí)間序列x0=n>=k0; %生成離散信號(hào)x0(n)x1=(n-k1)>=0; %生成離散信號(hào)x1(n)x2=(n-k2)>=0; %生成離散信號(hào)x2(n)subplot(3,1,1),stem(n, x0,'filled','k'); axis(n1,n2,1.1*min(x0),1.1*max(x0);ylabel('u(n)');subplot(3,1,2),stem(n

14、, x1,'filled','k'); axis(n1,n2,1.1*min(x1),1.1*max(x1);ylabel('u(n+6)');subplot(3,1,3),stem(n, x2,'filled','k'); axis(n1,n2,1.1*min(x2),1.1*max(x2);ylabel('u(n-4)');例、已知x（n）=，求x（n-2）和x（n+2）在-210區(qū)間的波形。MATLAB源程序?yàn)椋簄=-2:10;n0=2;n1=-2;x=2*sin(2*pi*n/10); %建

15、立原信號(hào)x(n)x1=2*sin(2*pi*(n-n0)/10); %建立x(n-2)信號(hào)x2=2*sin(2*pi*(n-n1)/10); %建立x(n+2)信號(hào)subplot(3,1,1),stem(n,x,'filled','k'); ylabel('x(n)');subplot(3,1,2),stem(n,x1,'filled','k'); ylabel('x(n-2)');subplot(3,1,3),stem(n,x2,'filled','k'); yla

16、bel('x(n+2)');2、 序列的相加情況1 兩序列具有相同的維數(shù)例、求x(n)= (0<n<10)MATLAB源程序?yàn)椋簄1=0;n2=10;n01=2;n02=4;n=n1:n2;x1=(n-n01)=0;x2=(n-n02)=0;x3=x1+x2;subplot(3,1,1);stem(n,x1,'filled','k');axis(n1,n2,1.1*min(x1),1.1*max(x1);ylabel('(n-2)');subplot(3,1,2);stem(n,x2,'filled',

17、'k');axis(n1,n2,1.1*min(x2),1.1*max(x2);ylabel('(n-4)');subplot(3,1,3);stem(n,x3,'filled','k');axis(n1,n2,1.1*min(x3),1.1*max(x3);ylabel('(n-2)+(n-4)');情況2 兩序列具有不同的維數(shù)例、x1（n）=u（n+2） （-4<n<6）X2(n)=u(n-4) (-5<n<8)求x（n）= x1（n）+ X2(n)=u(n-4)MATLAB源程序?yàn)椋簄

18、1=-4:6;n01=-2; x1=(n1-n01)>=0; %建立x1信號(hào)n2=-5:8;n02=4;x2=(n2-n02)>=0; %建立x2信號(hào)n=min(n1,n2):max(n1,n2); %為x信號(hào)建立時(shí)間序列nN=length(n); %求時(shí)間序列n的點(diǎn)數(shù)Ny1=zeros(1,N); %新建一維N列的y1全0數(shù)組y2=zeros(1,N); %新建一維N列的y2全0數(shù)組y1(find(n>=min(n1)&(n<=max(n1)=x1; %為y1賦值y2(find(n>=min(n2)&(n<=max(n2)=x2; %為y2

19、賦值x=y1+y2;subplot(3,1,1),stem(n1,x1,'filled','k');ylabel('x1(n)');axis(min(n),max(n),1.1*min(x),1.1*max(x);subplot(3,1,2),stem(n2,x2,'filled','k');ylabel('x2(n)');axis(min(n),max(n),1.1*min(x),1.1*max(x);subplot(3,1,3),stem(n,x,'filled','k&

20、#39;);ylabel('x(n)');axis(min(n),max(n),1.1*min(x),1.1*max(x);3、 序列相乘例、x1（n）= （-4<n<10）X2(n)=u(n+1) (-2<n<6) 求x（n）= x1（n）xX2(n)MATLAB源程序?yàn)椋簄1=-4:10; x1=3*exp(-0.25*n); %建立x1信號(hào)n2=-2:6;n02=-1;x2=(n2-n02)>=0; %建立x2信號(hào)n=min(n1,n2):max(n1,n2); %為x信號(hào)建立時(shí)間序列nN=length(n); %求時(shí)間序列n的點(diǎn)數(shù)Ny1=z

21、eros(1,N); %新建一維N列的y1全0數(shù)組y2=zeros(1,N); %新建一維N列的y2全0數(shù)組y1(find(n>=min(n1)&(n<=max(n1)=x1; %為y1賦值y2(find(n>=min(n2)&(n<=max(n2)=x2; %為y2賦值x=y1.*y2;subplot(3,1,1),stem(n1,x1,'filled','k');ylabel('x1(n)');axis(min(n),max(n),1.1*min(x1),1.1*max(x1);subplot(3,1,

22、2),stem(n2,x2,'filled','k');ylabel('x2(n)');axis(min(n),max(n),1.1*min(x2),1.1*max(x2);subplot(3,1,3),stem(n,x,'filled','k');ylabel('x(n)');axis(min(n),max(n),1.1*min(x),1.1*max(x);4、 序列反折例、x（n）= （-4<n<4）求x（-n）。MATLAB源程序?yàn)椋簄=-4:4;x=exp(-0.3*n);x1=

23、fliplr(x);n1=-fliplr(n);subplot(1,2,1),stem(n,x,'filled','k');title('x(n)');subplot(1,2,2),stem(n1,x1,'filled','k');title('x(-n)');5、 序列的尺度變換例、x（n）=sin（2*pi*n），求x（2n）和x（n/2）MATLAB源程序?yàn)椋簄=0:40;tn=n./20; %每周期取20個(gè)點(diǎn)x=sin(2*pi*tn); %建立原信號(hào)x(n)x1=sin(2*pi*tn*2

24、); %建立x(2n)信號(hào)x2=sin(2*pi*tn/2); %建立x(n/2)信號(hào)subplot(3,1,1),stem(tn,x,'filled','k'); ylabel('x(n)');axis(0,2,1.1*min(x),1.1*max(x);subplot(3,1,2),stem(tn,x1,'filled','k'); ylabel('x(2n)');axis(0,2,1.1*min(x),1.1*max(x);subplot(3,1,3),stem(tn,x2,'fill

25、ed','k'); ylabel('x(n/2)');axis(0,2,1.1*min(x),1.1*max(x); 四、實(shí)驗(yàn)任務(wù)1、x（n）= （-5<n<5）2、x（n）=u（n-2）+u（n+2） （-5<n<5）3、已知x（n）=3cos（2*pi*n/10）,顯示x（n-3）和x（n+3）在020區(qū)間的波形。4、已知x1=exp（- n /16），x2=5sin（2*pi*n/10）,顯示x1乘以x2在024區(qū)間的波形。5、已知x（n）=nsin（n），顯示020區(qū)間的波形 y1(n)=x(n-3),y2（n）=x（-n

26、），y3=-x(n),y4=x(-n+3),y5（n）=x（n/2）.實(shí)驗(yàn)3 模擬原型濾波器的設(shè)計(jì)一、實(shí)驗(yàn)?zāi)康膌 學(xué)會(huì)運(yùn)用MATLAB設(shè)計(jì)模擬低通濾波器原型的設(shè)計(jì)方法。二、實(shí)驗(yàn)涉及的matlab子函數(shù)1、buttord功能：確定巴特沃斯濾波器的階數(shù)和3dB截止頻率。調(diào)用格式：N,Wc= buttord(wp,ws,Rp,As,s)，計(jì)算巴特沃斯濾波器的階數(shù)和3dB截止頻率。Rp為通帶最大衰減，As為阻帶最小衰減。2、che1ord功能：確定切比雪夫1型濾波器的階數(shù)和通帶截止頻率。調(diào)用格式：N,Wp= buttord(wp,ws,Rp,As,s)3、che2ord功能：確定切比雪夫2型濾波器的階

27、數(shù)和阻帶截止頻率。調(diào)用格式： N,Wn= buttord(wp,ws,Rp,As,s) 4、ellipord功能：確定橢圓濾波器的階數(shù)和通帶截止頻率。調(diào)用格式： N,Wn= buttord(wp,ws,Rp,As,s)。5、buttap功能：巴特沃斯模擬低通濾波器原型（即歸一化的濾波器）。調(diào)用格式：z，p，k= buttap(n)，設(shè)計(jì)巴特沃斯模擬低通濾波器原型，其系統(tǒng)函數(shù)為6、cheb1ap功能：切比雪夫1型模擬低通濾波器原型。調(diào)用格式：z，p，k= cheb1ap (n，Rp)，設(shè)計(jì)切比雪夫1型模擬低通濾波器原型，其通帶最大衰減為Rp，系統(tǒng)函數(shù)為7、cheb2ap功能：切比雪夫2型模擬低通

28、濾波器原型。調(diào)用格式：z，p，k= cheb2ap (n，As)，設(shè)計(jì)切比雪夫2型模擬低通濾波器原型，其阻帶最小衰減為As，系統(tǒng)函數(shù)為8、ellipap功能：橢圓模擬低通濾波器原型。調(diào)用格式：z，p，k= ellipap (n，Rp，As)，設(shè)計(jì)橢圓模擬低通濾波器原型，其系統(tǒng)函數(shù)為9、poly功能：求某向量指定根所對(duì)應(yīng)的特征多項(xiàng)式。調(diào)用格式：P=poly（），多項(xiàng)式P是一個(gè)特征多項(xiàng)式，的元素是多項(xiàng)式P的根。例如：>>a=1 2 ；P=poly(a)；則P = 1 -3 210、poly2strPa= poly2str(a,s)例如：>> P= 1,-3,2; PA=po

29、ly2str(P,'s')則PA = s2 - 3 s + 211、pzmap功能：顯示連續(xù)系統(tǒng)的零極點(diǎn)分布圖。調(diào)用格式：pzmap（b，a）；繪制由行向量b和a構(gòu)成的系統(tǒng)函數(shù)確定的零極點(diǎn)分布圖pzmap(p，z)；繪制由列向量z確定的零點(diǎn)和p確定的極點(diǎn)構(gòu)成的零極點(diǎn)分布圖。三、實(shí)驗(yàn)原理由于IIR數(shù)字濾波器是在已知的歸一化的低通模擬濾波器的基礎(chǔ)上設(shè)計(jì)的，主要包括巴特沃斯、切比雪夫、橢圓低通濾波器，因此把這些低通濾波器稱為濾波器原型。下面介紹各種濾波器的設(shè)計(jì)。例1、 進(jìn)行巴特沃斯濾波器原型的設(shè)計(jì)，獲得任意階數(shù)N的歸一化的系統(tǒng)函數(shù)公式和幅頻響應(yīng)。MATLAB源程序?yàn)椋?巴特沃斯模擬濾

30、波器原型n=input('N= '); %輸入濾波器階數(shù)N%計(jì)算n階模擬低通原型，得到左半平面零極點(diǎn)z0,p0,k0=buttap(n); b0=k0*poly(z0) %求濾波器系數(shù)b0a0=poly(p0) %求濾波器系數(shù)a0h,w=freqs(b0,a0); %顯示系統(tǒng)的頻率特性plot(w,abs(h),'r'), %畫幅頻特性圖axis(0,5,0,1.1),ylable('幅度');xlabel('f（HZ）'); Pb=poly2str(b0,'p') %給出b0決定的關(guān)于p多項(xiàng)式Pa=poly2st

31、r(a0,'p') %給出a0決定的關(guān)于p多項(xiàng)式輸入階數(shù)N=2則運(yùn)行結(jié)果為b0 = 1a0 = 1.0000 1.4142 1.0000Pb = 1Pa = p2 + 1.4142 p + 1即歸一化的2階巴特沃斯濾波器原型的系統(tǒng)函數(shù)為例2、 通過模擬濾波器原型設(shè)計(jì)一個(gè)巴特沃斯模擬低通濾波器，要求通帶截止頻率2khz，通帶最大衰減1db，阻帶截止頻率5khz，阻帶最小衰減20db。MATLAB源程序?yàn)椋?巴特沃斯模擬濾波器fp=2000;wp=2*pi*fp %輸入濾波器的通帶截止頻率fs=5000;ws=2*pi*fs %輸入濾波器的阻帶截止頻率Rp=1;As=20; %輸入

32、濾波器的通阻帶衰減指標(biāo)%計(jì)算濾波器的階數(shù)和3dB截止頻率N,wc=buttord(wp,ws,Rp,As,'s') %計(jì)算n階模擬低通原型z0,p0,k0=buttap(n) %字母后加0表明這是原型濾波器的各指標(biāo)，而不是所求的濾波器的b0=k0*poly(z0) %求歸一化濾波器分子系數(shù)b0a0=poly(p0) %求歸一化濾波器分母系數(shù)a0H,w=freqs(b0,a0); %求歸一化系統(tǒng)的頻率特性dbH=20*log10(abs(H)/max(abs(H); %將歸一化系統(tǒng)的幅頻特性化為分貝值，注意此時(shí)的分貝值為負(fù)的，為了使分貝圖和幅頻響應(yīng)圖一致subplot(2,2,1

33、),plot(w*wc/(2*pi),abs(H),grid %畫所求濾波器的幅頻響應(yīng)圖axis(0,6000,0,1.1);ylabel('幅度');xlabel('f(Hz)');subplot(2,2,2),plot(w*wc/(2*pi),angle(H),grid %畫所求濾波器的相頻響應(yīng)圖axis(0,6000,-4,4);ylabel('相位');xlabel('f(Hz)');subplot(2,2,3),plot(w*wc/(2*pi),dbH),grid %畫所求濾波器的幅頻響應(yīng)分貝圖axis(0,6000,-

34、30,2);ylabel('幅度(dB)');xlabel('f(Hz)');subplot(2,2,4),plot(p0*wc,'x'); %畫所求濾波器的極點(diǎn)圖，用x號(hào)標(biāo)出axis square,axis equal,grid on %使x、y軸等比例顯示%以下程序可不寫wx0=wp,ws/wc; %設(shè)置頻率向量Hx=freqs(b0,a0,wx0); %計(jì)算通帶阻帶該兩點(diǎn)的頻率特性，即H(jwp)、H(jws)dbHx=-20*log10(abs(Hx)/max(abs(H) %將H(jwp)、H(jws)化為分貝值，即通帶最大衰減、阻帶最

35、小衰減運(yùn)行結(jié)果為wp = 1.2566e+004ws = 3.1416e+004N = 4wc = 1.7689e+004z0 = p0 = -0.3827 + 0.9239i （此為歸一化濾波器的極點(diǎn)，第四幅圖所畫的為所求濾波器的極點(diǎn)p0*wc） -0.3827 - 0.9239i -0.9239 + 0.3827i -0.9239 - 0.3827ik0 = 1b0 = 1 a0 = 1.0000 2.6131 3.4142 2.6131 1.0000dbHx = 0.2730 20.0000 （實(shí)際求出來的濾波器的通帶最大衰減、阻帶最小衰減，可見阻帶指標(biāo)剛好，通帶指標(biāo)有富余）所求濾波器的

36、系統(tǒng)函數(shù)為例3、 進(jìn)行切比雪夫1型濾波器原型的設(shè)計(jì)，獲得任意階數(shù)N的歸一化的系統(tǒng)函數(shù)公式和幅頻響應(yīng)。MATLAB源程序?yàn)椋簄=input('N= '); %由使用者輸入濾波器階數(shù)NRp=input('Rp='); %輸入濾波器通帶衰減指標(biāo)%計(jì)算n階模擬低通原型，得到左半平面零極點(diǎn)z0,p0,k0=cheb1ap(n,Rp); b0=k0*poly(z0) %求濾波器系數(shù)b0a0=poly(p0) %求濾波器系數(shù)a0h,w=freqs(b0,a0); %顯示系統(tǒng)的頻率特性plot(w,abs(h),'r'),axis(0,5,0,1.1),ylab

37、el('幅度');xlabel('f（HZ）'); title('切比雪夫1型原型模擬濾波器');Pb=poly2str(b0,'p') %給出b0決定的關(guān)于s多項(xiàng)式Pa=poly2str(a0,'p') %給出a0決定的關(guān)于s多項(xiàng)式輸入N=5,Rp=1運(yùn)行結(jié)果為b0 = 0.1228a0 = 1.0000 0.9368 1.6888 0.9744 0.5805 0.1228Pb = 0.12283Pa = p5 + 0.93682 p4 + 1.6888 p3 + 0.9744 p2 + 0.58053 p +

38、0.12283即5階通帶最大衰減1db的切比雪夫1型原形濾波器系統(tǒng)函數(shù)為例4、設(shè)計(jì)一個(gè)切比雪夫1型模擬低通濾波器，要求通帶截止頻率2khz，通帶最大衰減1db，阻帶截止頻率5khz，阻帶最小衰減20db。MATLAB源程序?yàn)椋?切比雪夫型模擬濾波器fp=2000;wp=2*pi*fp %輸入濾波器的通帶截止頻率fs=5000;ws=2*pi*fs %輸入濾波器的阻帶截止頻率Rp=1;As=20; %輸入濾波器的通阻帶衰減指標(biāo)%計(jì)算濾波器的階數(shù)和通帶截止頻率n,wp=cheb1ord(wp,ws,Rp,As,'s')%計(jì)算n階模擬低通原型，得到左半平面零極點(diǎn)z0,p0,k0=ch

39、eb1ap(n,Rp)b0=k0*poly(z0) %求濾波器系數(shù)b0a0=poly(p0) %求濾波器系數(shù)a0H,w=freqs(b0,a0); %求系統(tǒng)的頻率特性dbH=20*log10(abs(H)/max(abs(H); %化為分貝值subplot(2,2,1),plot(w*wp/(2*pi),abs(H),'k'),grid axis(0,6000,0,1.1);ylabel('幅度');xlabel('f(Hz)');subplot(2,2,2),plot(w*wp/(2*pi),angle(H),'k'),grid

40、axis(0,6000,-4,4);ylabel('相位');xlabel('f(Hz)');subplot(2,2,3),plot(w*wp/(2*pi),dbH,'k'),grid axis(0,6000,-50,2);ylabel('幅度(dB)');xlabel('f(Hz)');subplot(2,2,4),plot(p0*wp,'xk') %顯示系統(tǒng)的極點(diǎn)分布axis square,axis equal,grid on %使X,Y軸等比例顯示wx0=wp,ws/wp; %設(shè)置頻率向量Hx

41、=freqs(b0,a0,wx0); %計(jì)算該兩點(diǎn)的頻率特性dbHx=-20*log10(abs(Hx)/max(abs(H) %化為分貝值運(yùn)行結(jié)果為wp = 1.2566e+004ws = 3.1416e+004n = 3wp = 1.2566e+004z0 = p0 = -0.2471 + 0.9660i -0.4942 -0.2471 - 0.9660ik0 = 0.4913b0 = 0.4913a0 = 1.0000 0.9883 1.2384 0.4913dbHx = 0.9996 28.9442 （可見實(shí)現(xiàn)同樣的指標(biāo)，階數(shù)比巴特沃斯低，但性能卻更好）例5、進(jìn)行橢圓型濾波器原型的設(shè)計(jì)

42、，獲得任意階數(shù)N的歸一化的系統(tǒng)函數(shù)公式和幅頻響應(yīng)。MATLAB源程序?yàn)椋簄=input('N= '); %由使用者輸入濾波器階數(shù)NRp=input('Rp='); %輸入濾波器通帶衰減指標(biāo)As=input('As=');%計(jì)算n階模擬低通原型，得到左半平面零極點(diǎn)z0,p0,k0=ellipap(n,Rp,As); b0=k0*poly(z0) %求濾波器系數(shù)b0a0=poly(p0) %求濾波器系數(shù)a0h,w=freqs(b0,a0); %顯示系統(tǒng)的頻率特性plot(w,abs(h),'r'),axis(0,5,0,1.1),yl

43、abel('幅度');xlabel('f（HZ）'); title('橢圓型模擬濾波器原型');Pb=poly2str(b0,'p') %給出b0決定的關(guān)于s多項(xiàng)式Pa=poly2str(a0,'p') %給出a0決定的關(guān)于s多項(xiàng)式輸入N= 4 Rp=1 As=20運(yùn)行結(jié)果為b0 = 0.1000 0 0.5422 0 0.5256a0 = 1.0000 0.9038 1.6765 0.8688 0.5897Pb = 0.1 p4 + 0.5422 p2 + 0.52559Pa = p4 + 0.90377 p3

44、+ 1.6765 p2 + 0.86881 p + 0.58972例6、設(shè)計(jì)一個(gè)切比雪夫1型模擬低通濾波器，要求通帶截止頻率2khz，通帶最大衰減2db，阻帶截止頻率5khz，阻帶最小衰減50db。MATLAB源程序?yàn)椋?橢圓模擬濾波器fp=2000;wp=2*pi*fp %輸入濾波器的通帶截止頻率fs=5000;ws=2*pi*fs %輸入濾波器的阻帶截止頻率Rp=2;As=50; %輸入濾波器的通阻帶衰減指標(biāo)%計(jì)算濾波器的階數(shù)和通帶截止頻率n,wn=ellipord(wp,ws,Rp,As,'s') %當(dāng)不知道對(duì)哪個(gè)頻率歸一化時(shí)，用wn表示，最后在相應(yīng)的指標(biāo)處去除歸一化即可

45、%計(jì)算n階模擬低通原型，得到左半平面零極點(diǎn)z0,p0,k0=ellipap(n,Rp,As) b0=k0*poly(z0) %求濾波器系數(shù)b0Pb=poly2str(b0,'s') %給出b0決定的關(guān)于s多項(xiàng)式a0=(poly(p0) %求濾波器系數(shù)a0Pa=poly2str(a0,'s') %給出a0決定的關(guān)于s多項(xiàng)式H,w=freqs(b0,a0); %求系統(tǒng)的頻率特性dbH=20*log10(abs(H)/max(abs(H); %化為分貝值subplot(2,2,1),plot(w*wn/(2*pi),abs(H),grid %去除歸一化axis(0,1

46、0000,-0.1,1.1);ylabel('幅度');xlabel('f(Hz)');subplot(2,2,2),plot(w*wn/(2*pi),angle(H),grid %去除歸一化axis(0,10000,-4,4);ylabel('相位');xlabel('f(Hz)');subplot(2,2,3),plot(w*wn/(2*pi),dbH),grid axis(0,10000,-100,2);ylabel('幅度(dB)');xlabel('f(Hz)');subplot(2,2,

47、4),pzmap(b0,a0); %此為歸一化原型濾波器的零極點(diǎn)，不是所求濾波器的零極點(diǎn)wx0=wp,ws/wn; %設(shè)置頻率向量Hx=freqs(b0,a0,wx0); %計(jì)算該兩點(diǎn)的頻率特性dbHx=-20*log10(abs(Hx)/max(abs(H) %化為分貝值運(yùn)行結(jié)果為wp = 1.2566e+004ws = 3.1416e+004n = 4wn = 1.2566e+004 （wn=wp可見橢圓濾波器是對(duì)通帶截止頻率進(jìn)行歸一化的）z0 = 0 - 4.2226i 0 + 4.2226i 0 - 1.8719i 0 + 1.8719ip0 = -0.2673 - 0.4381i -

48、0.2673 + 0.4381i -0.0879 - 0.9678i -0.0879 + 0.9678ik0 = 0.0032b0 = 0.0032 0 0.0675 0 0.1976Pb = 0.0031623 s4 + 0.067465 s2 + 0.19757a0 = 1.0000 0.7104 1.3017 0.5512 0.2487Pa = s4 + 0.71038 s3 + 1.3017 s2 + 0.5512 s + 0.24872dbHx = 1.9993 50.2300四、實(shí)驗(yàn)任務(wù)設(shè)計(jì)一個(gè)模擬低通濾波器，通帶截止頻率為6khz，通帶最大衰減為1db，阻帶截止頻率為15khz，

49、阻帶最小衰減為30db。要求：分別實(shí)現(xiàn)符合以上指標(biāo)的巴特沃斯、切比雪夫1型、橢圓濾波器，繪制幅頻特性、幅頻特性分布圖、相頻特性、零極點(diǎn)分布圖，并列寫系統(tǒng)函數(shù)表示式，注意：是所求系統(tǒng)的系統(tǒng)函數(shù)，而不是歸一化的系統(tǒng)函數(shù)。實(shí)驗(yàn)4 模擬域的頻率變換法一、實(shí)驗(yàn)?zāi)康膌 學(xué)會(huì)運(yùn)用MATLAB設(shè)計(jì)模擬低通、帶通、高通、帶阻濾波器的設(shè)計(jì)方法。二、實(shí)驗(yàn)涉及的matlab子函數(shù)1、lp2lp功能：低通到低通模擬濾波器變換。調(diào)用格式：bt,at= lp2lp(b,a,w0)，將系統(tǒng)函數(shù)表示的截止頻率為1rad/s的模擬低通濾波器原型變換為截止頻率為w0的低通濾波器。2、lp2hp功能：低通到高通模擬濾波器變換。調(diào)用格

50、式：bt,at= lp2hp(b,a,w0)，將系統(tǒng)函數(shù)表示的截止頻率為1rad/s的模擬低通濾波器原型變換為截止頻率為w0的高通濾波器。3、lp2bp功能：低通到帶通模擬濾波器變換。調(diào)用格式：bt,at= lp2bp(b,a,w0,bw)，將系統(tǒng)函數(shù)表示的截止頻率為1rad/s的模擬低通濾波器原型變換為中心頻率為w0、帶寬為bw的帶通濾波器。如果被設(shè)計(jì)的濾波器低端截止頻率為w1，高端截止頻率為w2，則W0=sqrt（w1*w2），bw=w2-w14、lp2bs功能：低通到帶阻模擬濾波器變換。調(diào)用格式：bt,at= lp2bs(b,a,w0,bw)，將系統(tǒng)函數(shù)表示的截止頻率為1rad/s的模擬

51、低通濾波器原型變換為中心頻率為w0、帶寬為bw的帶阻濾波器。如果被設(shè)計(jì)的濾波器低端截止頻率為w1，高端截止頻率為w2，則W0=sqrt（w1*w2），bw=w2-w1。三、實(shí)驗(yàn)原理下面舉例說明用頻率變換法設(shè)計(jì)各種模擬濾波器的過程。例1、通過頻率變換法設(shè)計(jì)一個(gè)巴特沃斯模擬低通濾波器，要求通帶截止頻率2khz，通帶最大衰減1db，阻帶截止頻率5khz，阻帶最小衰減20db。MATLAB源程序?yàn)椋?將巴特沃斯模擬濾波器原型變換成實(shí)際模擬低通濾波器fp=2000;wp=2*pi*fp; %輸入實(shí)際濾波器的通帶截止頻率fs=5000;ws=2*pi*fs; %輸入實(shí)際濾波器的阻帶截止頻率Rp=1;As=

52、20; %輸入濾波器的通阻帶衰減指標(biāo)%計(jì)算濾波器的階數(shù)和3dB截止頻率n,wc=buttord(wp,ws,Rp,As,'s')%計(jì)算n階模擬低通原型，得到左半平面零極點(diǎn)z0,p0,k0=buttap(n); b0=k0*poly(z0) %求歸一化的濾波器系數(shù)b0a0=poly(p0) %求歸一化的濾波器系數(shù)a0 H,w0=freqs(b0,a0); %求歸一化的濾波器頻率特性dbH=20*log10(abs(H)/max(abs(H); %幅度化為分貝值%變換為實(shí)際模擬低通濾波器bt,at=lp2lp(b0,a0,wc); %從歸一化低通變換到實(shí)際低通Ht,wt=freqs

53、(bt,at); %求實(shí)際系統(tǒng)的頻率特性dbHt=20*log10(abs(Ht)/max(abs(Ht); %幅度化為分貝值%實(shí)際模擬低通頻率特性作圖subplot(2,2,1),plot(wt/2/pi,abs(Ht);axis(0,6000,0,1.1);ylabel('幅度');xlabel('f(Hz)');gridsubplot(2,2,2),plot(wt/2/pi,dbHt);axis(0,6000,-30,1);ylabel('幅度dB');xlabel('f(Hz)');gridsubplot(2,2,3),p

54、lot(wt/2/pi,angle(Ht);axis(0,6000,-4,4);ylabel('相位');xlabel('f(Hz)');grid運(yùn)行結(jié)果為n = 4wc = 1.7689e+004b0 = 1a0 = 1.0000 2.6131 3.4142 2.6131 1.0000例2、通過頻率變換法設(shè)計(jì)一個(gè)巴特沃斯模擬高通濾波器，要求通帶截止頻率5khz，通帶最大衰減1db，阻帶截止頻率2khz，阻帶最小衰減20db。MATLAB源程序?yàn)椋?將巴特沃斯模擬濾波器原型變換成實(shí)際模擬高通濾波器fp=5000;wp=2*pi*fp; %輸入實(shí)際濾波器的通帶截止頻率fs=2000;ws=2*pi*fs; %輸入實(shí)際濾波器的阻帶截止頻率Rp=1;As=20; %輸入濾波器的通阻帶衰減指標(biāo)%計(jì)算濾波器的階數(shù)和3dB截止頻率n,wc=buttord(wp,ws,Rp,As,'s')%計(jì)算n階模擬低通原型，得到左半平面零極點(diǎn)z0,p0,k0=buttap(n); b0=k0*poly(z0) %求歸一化的濾波器系數(shù)b0a0=pol