實(shí)驗(yàn)部分：matlab在數(shù)字信號(hào)處理中的應(yīng)用

1、Matlab在數(shù)字信號(hào)處理中的應(yīng)用（基礎(chǔ)）一、數(shù)據(jù)類型：1、整數(shù)：Matlab支持8位，16位，32位和64位的有符號(hào)和無符號(hào)整數(shù)數(shù)據(jù)類型。如：x=int8(50); %指定x的數(shù)據(jù)類型為int8. x=502、浮點(diǎn)數(shù)：matlab的默認(rèn)數(shù)據(jù)類型是雙精度類型（double），為了節(jié)省存蓄空間，matlab也支持單精度數(shù)據(jù)類型的數(shù)組。Realmin(single)Ans= 1.1755e-038Realmax(double) Ans= 2.2251e-308 3、復(fù)數(shù)：matlab中虛數(shù)單位由i或者j表示。 Z=6+7j另一種創(chuàng)建復(fù)數(shù)的方法可以通過complex()函數(shù)，complex()函數(shù)的

2、調(diào)用格式：C=complex(a,b)，返回結(jié)果c為復(fù)數(shù)，實(shí)部是a，虛部是b。二、數(shù)組的創(chuàng)建 1、一維數(shù)組的創(chuàng)建：創(chuàng)建一維行向量，只需要把所有數(shù)組元素用空格或者逗號(hào)分隔，并用方括號(hào)吧所有數(shù)組元素括起來即可。如用分號(hào)，即為列向量。 創(chuàng)建等差的一維數(shù)組：格式 Var=start-val:step:stop-val。如果步長是1，可以省略。 2、二維數(shù)組的創(chuàng)建；在創(chuàng)建二維數(shù)組時(shí)，用逗號(hào)或者空格區(qū)分同一行的不同元素，用分號(hào)或者軟回車區(qū)分不同的行。三、函數(shù)流程控制 1、順序結(jié)構(gòu)。 2、判斷語句（if-else if-else-end）.3、循環(huán)語句（for-end）四、作圖 1、二維圖：plot(x,y

3、,linespec)，linespec參數(shù)，用于對(duì)圖像外觀屬性的控制，包括線條的形狀，顏色和點(diǎn)的形狀，顏色。stem(x,y);繪制脈沖桿圖圖形。Stairs(x,y);繪制階梯圖圖形。 2、圖像子窗口：subplot(m,n,p)，將圖像分為mn個(gè)子區(qū)域，在第p個(gè)區(qū)域中繪制圖像。 3、坐標(biāo)軸：axis(xmin,xmax,ymin,ymax).指定當(dāng)前圖像中x軸和y軸的范圍。 4、圖形注釋：1）標(biāo)題：title(圖形名字)。2）坐標(biāo)軸名：xlabel(x軸的名稱)，ylabel(y軸的名稱)。特殊符號(hào)的輸入: alpha的輸入，則自動(dòng)轉(zhuǎn)變成，實(shí)驗(yàn)一、幾種典型離散時(shí)間序列Matlab中處理的數(shù)

4、組，將下標(biāo)放在變量后面的小括號(hào)內(nèi)，且約定從1開始遞增。例如：x=5,4,3,2,1,0，表示x(1)=5, x(2)=4, x(3)=3, x(4)=2, x(5)=1, x(6)=0。要表示一個(gè)下標(biāo)不由1開始的數(shù)組x(n)，一般應(yīng)采用兩個(gè)矢量，如：n=-3:5;x=1,-1,3,2,0,-2,-1,2,1;這表示一個(gè)含有9個(gè)點(diǎn)的矢量，n為一組時(shí)間矢量，對(duì)應(yīng)x有：x(-3)=1, x(-2)=-1.。 連續(xù)信號(hào)作圖使用plot()函數(shù)，繪制線性圖。離散信號(hào)作圖使用stem()函數(shù)，繪制脈沖桿圖。一些常用的函數(shù)：abs():求絕對(duì)值（幅值）。調(diào)用格式：y=abs(x)。length():取某一變

5、量的長度（采樣點(diǎn)數(shù)）。調(diào)用格式：N=length(n)，取n的點(diǎn)數(shù)，賦值給N。real()：取一個(gè)復(fù)數(shù)的實(shí)部，調(diào)用格式：x=real(h)；取復(fù)數(shù)h的實(shí)部，賦值給變量x。imag():取復(fù)數(shù)的虛部，調(diào)用格式：x=imag(h)；取復(fù)數(shù)h的實(shí)部，賦值給變量yx=sawtooth(t);類似于sin(t)，產(chǎn)生周期為2pi，幅值從-1到+1的鋸齒波。x=sowtooth(t,width);產(chǎn)生三角波，其中width(0<width<1)為標(biāo)量用于確定最大值的位置。x=square(t);產(chǎn)生類似于sin(t),周期為2pi，幅值我1的方波，x=square(t,duty),產(chǎn)生指定周期

6、的矩形波，其中duty用于指定脈沖寬度與整個(gè)周期的比例。rand(n,m);產(chǎn)生一組具有n行m列的隨機(jī)信號(hào)。1、 單位沖激序列：1）利用零序列：x=zeros(1,N)，生成一個(gè)1N維的零向量。2）利用邏輯關(guān)系表達(dá)式產(chǎn)生單位沖激序列：x=(n-n0)=0；只在n=n0的地方產(chǎn)生1.例：MATLAB程序如下：%采樣邏輯關(guān)系求脈沖序列。 n1=-5;n2=5;n0=0; n=n1:n2;x=n=n0;%作圖部分stem(n,x,filled);axis(n1,n2,0,1.1*max(x);title(單位脈沖序列)；xlabel(時(shí)間(n);ylabel(幅度x(n);%采樣零序列求脈沖序列。n

7、1=-5;n2=5; k=0;n=n1:n2;nt=length(n); %求采樣點(diǎn)n的個(gè)數(shù)（長度）。nk=abs(k-n1)+1;x=zeros(1,nt);x(nk)=1;%作圖同上。2、單位階躍序列： 1）利用1序列：x=ones(1,N)，產(chǎn)生一個(gè)1N維的全1向量。2）利用邏輯關(guān)系表達(dá)式產(chǎn)生單位階躍序列：x=(n-n0)>=0。Matlab程序： n=0:49; x=ones(1,50); close all;stem(n,x);title(單位階躍信號(hào)序列)；3、單位矩形序列： 1）x=ones(1,N)， 2）利用邏輯關(guān)系表達(dá)式產(chǎn)生：x=(n-n0>=0)&(n

8、-nf<=0)。matlab程序： N=10;n=0:49;x=sign(sign(N-1-n)+1);close all; %關(guān)閉所有打開的圖形窗口stem(n,x);注：sign(x), 符號(hào)函數(shù)，當(dāng)x大于0時(shí)值為1，當(dāng)x等于0時(shí)值為0，當(dāng)x小于0時(shí)值為-1.4、正弦序列： x=a*sin(omega*n+thwlta)； x=a*sin(2*pi*f0/Fs*n+thelta)；例：頻率為1.振幅為1的正弦信號(hào)，在窗口中顯示2個(gè)周期的信號(hào)波形，并對(duì)該信號(hào)的一個(gè)周期進(jìn)行32點(diǎn)采樣。獲得離散信號(hào)。做出連續(xù)信號(hào)和離散信號(hào)的圖形。MATLAB程序如下： f=1; Um-1; nt=2; %

9、頻率，振幅，周期的個(gè)數(shù)。 N=32; T=1/f; %采樣點(diǎn)數(shù)， 周期dt=T/N; %采樣時(shí)間間隔 n=0:nt*N-1; tn=n*dt;x=Um*sin(2*f*pi*tn);%作圖部分 subplot(2,1,1), plot(tn, x); axis(0, nt*T, 1.1*min(x), 1.1*max(x);ylabel(連續(xù)正弦信號(hào)x(t);subplot(2,1,2), stem(tn, x);axis(0, nt*T, 1.1*min(x), 1.1*max(x);ylabel(離散正弦序列x(n); 5、實(shí)指數(shù)序列 x(n)=an; 例：編寫產(chǎn)生a=1/2和a=2的實(shí)指

10、數(shù)連續(xù)信號(hào)和離散信號(hào)序列的程序。 MATLAB程序如下： n1=-10; n2=10; a1=1/2; a2=2; na1=-10:0; na2=0:10; x1=a1.na1; x2=a2.na2; %作圖部分 subplot(2,2,1), plot(na1,x1); %作連續(xù)圖形 title(實(shí)指數(shù)原信號(hào)（a<1）); subplot(2,2,2), stem(na1,x1,filled); %作離散圖形。 title(實(shí)指數(shù)序列(a<1);subplot(2,2,3), plot(na2,x2); %作連續(xù)圖形 title(實(shí)指數(shù)原信號(hào)（a>1）); subplot(

11、2,2,4), stem(na2,x2,filled); %作離散圖形。 title(實(shí)指數(shù)序列(a>1);6、復(fù)指數(shù)序列： x=exp(sigma+jomega)*n);7、矩形波序列：y=rectpuls(t,width).該函數(shù)產(chǎn)生一個(gè)幅度為1寬度為width，且以t=0為對(duì)稱軸的矩形脈沖信號(hào)，width的默認(rèn)值為1.y=square(t,DUTY),產(chǎn)生一個(gè)周期為2*pi,幅值為+1（-1）的周期性方波信號(hào)。其中DUTY表示信號(hào)的占空比。默認(rèn)值為0.5.例：矩形脈沖信號(hào)的波形圖： 2 (0<=t<=1) f(t)= 0 (t<1,t>1)MATLAB程序如

12、下：t=-0.5:0.01:3;t0=0.5;width=1;ft=2*rectpuls(t-t0,width);plot(t,ft);grid on;axis(-0.5,3,0.2,2.2);title(矩形脈沖信號(hào))；例：產(chǎn)生一個(gè)頻率為10HZ，占空比為30%的周期方波信號(hào)。MATLAB程序如下：t=0:0.001:3;y=square(2*pi*10*t,30);plot(t,y);aixs(0,0.3,-1.2,1.2);title(周期方波信號(hào))；實(shí)驗(yàn)二、序列的基本運(yùn)算 1、序列的加法和乘法 x=x1+x2; x=x1.*x2例：已知x1(n)=u(n+2) (-4<n<

13、6) x2(n)=u(n-4) (-5<n<8) 求：x(n)=x1(n)+x2(n) MATLAB程序如下：n1=-4:6; n01=-2;x1=(n1-n01)>=0;n2=-5:8;n02=4;x2=(n2-n02)>=0;%用0值來擴(kuò)展它們的序列號(hào)，變成相同的起點(diǎn)和終點(diǎn)，原來的值不變。n=min(n1,n2):max(n1,n2);N=length(n);y1=zeros(1,N); y2=zeros(1,N);y1(find(n>=min(n1)&(n<=max(n1)=x1;y2(find(n>=min(n2)&(n<

14、=max(n2)=x2;x=y1+y2;subplot(3,1,1),stem(n,y1); axis(min(n),max(n),1.1*min(x),1.1*max(x);subplot(3,1,2),stem(n,y2); axis(min(n),max(n),1.1*min(x),1.1*max(x);subplot(3,1,3),stem(n, x);axis(min(n),max(n),1.1*min(x),1.1*max(x);注：序列的乘法與上程序相同。2、 序列的翻轉(zhuǎn)： 翻轉(zhuǎn)運(yùn)算用fliplr()函數(shù)實(shí)現(xiàn)，設(shè)序列x(n)，樣值向量x和位置向量nx表示，則翻轉(zhuǎn)之后的序列y(n)

15、的樣值向量y和位置向量ny表示，則y=fliplr(x); ny=-fliplr(nx);3、 序列的移位設(shè)序列x(n)，樣值向量x和位置向量nx表示，移位n0之后的序列y(n)的樣值向量y和位置向量ny表示，則y=x;ny=nx+n0;例：已知一正弦信號(hào)：x(n)=2sin(2pi*n/10)求其移位信號(hào)x(n-2)在-2<n<10區(qū)間的序列波形？MATLAB程序如下： n=-2:10; n0=2；x=2*sin(2*pi*n/10); %建立原信號(hào)x(n)x1=2*sin(2*pi*(n-n0)/10); %建立x(n-2)信號(hào)subplot(2,1,1), stem(n,x,

16、filled);ylabel(x(n);subplot(2,1,2), stem(n,x1,filled);4、 離散序列卷積：MATLAB提供一個(gè)conv函數(shù)功能：進(jìn)行兩個(gè)序列間的卷積運(yùn)算。調(diào)用格式：y=conv(x,h),用于求取兩個(gè)有限長度序列x和h的卷積，y的長度等于x和h的長度之和減1。注意：conv函數(shù)默認(rèn)兩個(gè)信號(hào)的時(shí)間序列從n=0開始的。例：已知兩個(gè)信號(hào)序列； (0<n<20) (0<n<10)2、如果兩個(gè)信號(hào)不是從n=0開始的，則采用y , ny=conv_new（x , nx , h , nh）函數(shù)。其中x 是輸入序列，nx是它的序列號(hào)，h是另一個(gè)序列

17、，nh是它的序列號(hào)，y是卷積和，ny是它的序列號(hào)。%conv_new.m實(shí)現(xiàn)任意位置序列卷積運(yùn)算，返回值是卷積和的值y和時(shí)間向量ny。functiony,ny=conv_new(x,nx,h,nh)n1=nx(1)+nh(1);n2=nx(length(x)+nh(length(h);ny=n1:n2;y=conv(x,h);例：x=3,11,7,0,-1,4,2; nx=-3:3;h=2,3,0,-5,2,1,; nh=-1:4;y,ny=conv_new(x,nx,h,nh);subplot(3,1,1),stem(nx,x); axis(min(nx),max(nx),1.1*min(x

18、),1.1*max(x);subplot(3,1,2),stem(nh,h); axis(min(nh),max(nh),1.1*min(h),1.1*max(h);subplot(3,1,3),stem(ny, y);axis(min(ny),max(ny),1.1*min(y),1.1*max(y);實(shí)驗(yàn)三：離散系統(tǒng)的沖激響應(yīng)和階躍響應(yīng)1、 impz()：功能：求解數(shù)字系統(tǒng)的沖激響應(yīng)。調(diào)用格式：h,t=impz(b,a);求解數(shù)字系統(tǒng)的沖激響應(yīng)h，取樣點(diǎn)數(shù)為缺省值。h,t=impz(b,a,n);求解數(shù)字系統(tǒng)的沖激響應(yīng)h，取樣點(diǎn)數(shù)為n值。impz(b,a);在當(dāng)前窗口用stem(t,h )

19、函數(shù)出圖。2、 dstep():功能：求解數(shù)字系統(tǒng)的階躍響應(yīng)。調(diào)用格式：h,t=dstep(b,a);求解數(shù)字系統(tǒng)的階躍響應(yīng)h，取樣點(diǎn)數(shù)為缺省值。h,t=dstep(b,a,n);求解數(shù)字系統(tǒng)的階躍響應(yīng)h，取樣點(diǎn)數(shù)為n值。dstep(b,a)；在當(dāng)前窗口用stairs(t,h)函數(shù)出圖。3、filter子函數(shù)功能：對(duì)數(shù)字系統(tǒng)的輸入信號(hào)進(jìn)行濾波處理。因?yàn)橐粋€(gè)離散系統(tǒng)可以看作是一個(gè)濾波器，系統(tǒng)的輸出就是輸入經(jīng)過濾波器濾波的結(jié)果。調(diào)用格式：y=filtet(b , a, x)，對(duì)于由矢量b, a決定的數(shù)字系統(tǒng)（b和a分別表示系統(tǒng)函數(shù)H（z）對(duì)應(yīng)的分子項(xiàng)和分母項(xiàng)系數(shù)構(gòu)成的數(shù)組，而且分母系數(shù)要?dú)w一化處理

20、。）當(dāng)輸入信號(hào)為x時(shí)，對(duì)x中的數(shù)據(jù)進(jìn)行濾波，結(jié)果存于y中，長度取max(na , nb). y , zf=filter(b , a, x);除得到結(jié)果矢量y外，還得到x 的最終狀態(tài)矢量zf。 y=filter(b , a, x, zi)；可在zi中指定x的初始狀態(tài)。4、filtic子函數(shù)功能：為filter子函數(shù)選擇初始條件。調(diào)用格式：zi=filtic(b ,a ,y ,x)；求給定輸入x和y時(shí)的初始狀態(tài)。 zi=filtic(b , a, y)；求x=0，給定輸入y時(shí)的初始狀態(tài)。其中，x和y分別是表示過去的輸入和輸出。例：已知一個(gè)因果系統(tǒng)的差分方程為 6y(n)+2y(n-2)=x(n)+

21、3x(n-1)+3x(n-2)+x(n-3)滿足初始條件y(-1)=0,x(-1)=0.求系統(tǒng)的單位沖激響應(yīng)和階躍響應(yīng)。 將上述方程對(duì)y(n)項(xiàng)系數(shù)進(jìn)行歸一化，得到其系統(tǒng)函數(shù)分子和分母系數(shù) a0=1, a1=0, a2=1/3, a3=0 b0=1/6, b1=1/2, b2=1/2, b3=1/6 用impz()函數(shù)的MATLAB程序（取N=32點(diǎn)作圖） a=1,0,1/3, 0; b=1/6, 1/2, 1/2, 1/6;N=32; n=0:N-1;hn=impz(b,a,n);gn=dstep(b,a,n);subplot(1,2,1),stem(n,hn,filled); title(

22、系統(tǒng)的單位階躍響應(yīng))；ylabel(h(n);xlabel(n);axis(0,N-1,-1.1*min(hn),1.1*max(hn);subplot(1,2,2), stem(n,gn,k);title(系統(tǒng)的單位階躍響應(yīng))；ylabel(g(n);xlabel(n);axis(0,N-1,-1.1*min(gn),1.1*max(gn);用filter()函數(shù)的MATLAB程序（取N=32點(diǎn)作圖） a=1,0,1/3, 0; b=1/6, 1/2, 1/2, 1/6;xi=filtic(b,a,0,0);N=32; n=0:N-1;x1=n=0; %單位沖激信號(hào)hn=filter(b,a

23、,x1,xi);subplot(1,2,1),stem(n,hn,filled); title(系統(tǒng)的單位階躍響應(yīng))；ylabel(h(n);xlabel(n);axis(0,N-1,-1.1*min(hn),1.1*max(hn);x2=n>=0; %單位階躍信號(hào)gn=filter(b,a,x2,xi);subplot(1,2,2), stem(n,gn,k);title(系統(tǒng)的單位階躍響應(yīng))；ylabel(g(n);xlabel(n);axis(0,N-1,-1.1*min(gn),1.1*max(gn);注：1、hold() 控制當(dāng)前圖形是否刷新的雙向切換開關(guān)。hold on使當(dāng)前

24、軸及圖形保持而不被刷新，準(zhǔn)備接受此后將繪制的新曲線。hold off使當(dāng)前軸及圖形不再具備刷新的性質(zhì)。3、 pause() 暫停執(zhí)行文件，等待用戶按任意鍵繼續(xù)，pause(n) 在繼續(xù)執(zhí)行之前，暫停n秒。實(shí)驗(yàn)四： 離散LSI系統(tǒng)的時(shí)域響應(yīng)對(duì)于離散LSI系統(tǒng)的響應(yīng)，MATLAB為我們提供了多種求解方法：（1） 用conv子函數(shù)進(jìn)行卷積積分，求任意輸入的系統(tǒng)零狀態(tài)響應(yīng)。（2） 用dlsim子函數(shù)求任意輸入的系統(tǒng)零狀態(tài)響應(yīng)。（3） 用filter和filtic子函數(shù)求任意輸入的系統(tǒng)完全響應(yīng)。1、 dlsim子函數(shù)功能：求解離散系統(tǒng)的響應(yīng)。調(diào)用格式：y=dlsim(b , a , x),求輸入信號(hào)為x

25、時(shí)系統(tǒng)的響應(yīng)。其中，b和a分別表示系統(tǒng)函數(shù)H（z）中，由對(duì)應(yīng)的分子項(xiàng)和分母項(xiàng)系數(shù)構(gòu)成的數(shù)組，而且分母系數(shù)要?dú)w一化處理。例：（書本P96頁15題）已知一個(gè)用以下差分方程表示的線性移不變因果系統(tǒng)為： 當(dāng)激勵(lì)時(shí)，求系統(tǒng)的響應(yīng)。 2、 用filter和filtic子函數(shù)求LSI系統(tǒng)對(duì)任意輸入的響應(yīng)Filter子函數(shù)功能：對(duì)數(shù)字系統(tǒng)的輸入信號(hào)進(jìn)行濾波處理。調(diào)用格式：y=filtet(b , a, x)，對(duì)于由b和a決定的數(shù)字系統(tǒng)（b和a分別表示系統(tǒng)函數(shù)H（z）中，由對(duì)應(yīng)的分子項(xiàng)和分母項(xiàng)系數(shù)構(gòu)成的數(shù)組，而且分母系數(shù)要?dú)w一化處理。）當(dāng)輸入信號(hào)為x時(shí)，對(duì)x中的數(shù)據(jù)進(jìn)行濾波，結(jié)果存于y中，長度取max(na ,

26、 nb). y , zf=filter(b , a, x);除得到結(jié)果矢量y外，還得到x 的最終狀態(tài)矢量zf。 y=filter(b , a, x, zi)；可在zi中指定x的初始狀態(tài)。Filtic子函數(shù)功能：為filter子函數(shù)選擇初始條件。調(diào)用格式：zi=filtic(b ,a ,y ,x)；求給定輸入x和y時(shí)的初始狀態(tài)。 zi=filtic(b , a, y)；求x=0，給定輸入y時(shí)的初始狀態(tài)。其中，x和y分別是表示過去的輸入和輸出。例：已知一個(gè)因果系統(tǒng)的差分方程為： 6y(n)-2y(n-4)=x(n)-3x(n-2)+3x(n-4)-x(n-6)滿足初始條件y(-1)=0,x(-1)

27、=0，求系統(tǒng)的單位沖激響應(yīng)和單位階躍響應(yīng)，時(shí)間軸上N取32點(diǎn)作圖。 例：已知一個(gè)IIR數(shù)字低通濾波器的系統(tǒng)函數(shù)為： =輸入一個(gè)矩形信號(hào)序列 x=square(n/5) (-2<n<10)作業(yè)：書本p96頁16題。作業(yè)：1、已知離散線性時(shí)不變的系統(tǒng)函數(shù)，請分別用impz和dstep子函數(shù)，filter和filtic子函數(shù)兩種方法求解系統(tǒng)的沖激響應(yīng)和階躍響應(yīng)。 （1） （2） 2、一個(gè)LSI系統(tǒng)的系統(tǒng)函數(shù)表示式為： 滿足初始條件y(-1)=5,y(-2)=5，試用filtet和filtic子函數(shù)求此系統(tǒng)的輸入序列x(n)為下列信號(hào)時(shí)的零輸入，零狀態(tài)以及完全響應(yīng)。（1） x(n)=（2）

28、 x(n)=實(shí)驗(yàn)五： z變換及其應(yīng)用1、ztrans子函數(shù)功能：返回?zé)o限長序列函數(shù)x(n)的z變換。調(diào)用格式：z=ztrans(x)；求無限長序列函數(shù)x(n)的z變換X(z)，返回z變換的表達(dá)式。2、iztrans子函數(shù)功能：求函數(shù)X（z）的z反變換x(n).調(diào)用格式：x=iztrans(X(z)；求函數(shù)X（z）的z反變換x(n)，返回z反變換的表達(dá)式。3、syms子函數(shù) 功能：定義多個(gè)符號(hào)對(duì)象。 調(diào)用格式：syms a b wo;把字符a b wo定義為基本的符號(hào)對(duì)象。4、residuez子函數(shù)功能：有理多項(xiàng)式的部分分式展開。調(diào)用格式：r p c=residuez(b, a)；把b(z)/a

29、(z)展開成部分分式的形式。 b , a=residuez(r p c)；根據(jù)部分分式的r p c 數(shù)組，返回有理多項(xiàng)式。其中：b ,a為按降冪排列的多項(xiàng)式的分子和分母的系數(shù)組；r為余數(shù)數(shù)組；p為極點(diǎn)數(shù)組；c為無窮限多項(xiàng)式系數(shù)數(shù)組。有理多項(xiàng)式如下： X（z）=注意：利用ztrans()子函數(shù)時(shí)，它只給出了z變換的表達(dá)式，而沒有給出收斂域。另外，由于這一功能還不盡完美，因而有的序列的z變換還求不出來，z的反變換也存在同樣的問題。例：用部分分式法求解函數(shù) 的z反變換，寫出h(n)的表達(dá)式，并用圖形與impz求得的結(jié)果相比較。MATLAB程序： %求z的反變換 b=0,1,0; a=1,-12,36

30、; r,p,c=residuez(b,a) %由此可知，這個(gè)多項(xiàng)式含有重極點(diǎn)，多項(xiàng)式分解后表示為： H(z)=-0.1667/(1-6z-1)+0.1667/(1-6z-1)2 =-0.1667/(1-6z-1)+0.1667z/6*6z-1/(1-6z-1)2 根據(jù)時(shí)域位移性質(zhì)，可寫出z反變換公式 h(n)=-0.1667(6)nu(n)+0.1667/6 *(n+1)6n+1u(n+1).%作圖 N=8;n=0:N-1; h=r(1)*p(1).n.*n>=0+r(2).*(n+1).*p(2).n.*n-1>=0; subplot(1,2,2),stem(n,h);title

31、(用部分分式法求反變換h(n);%用沖擊響應(yīng)法求h(n) h1=impz(b,a,N); subplot(1,2,2),stem(n,h1); title(用impz()求反變換h(n);例：已知一個(gè)離散系統(tǒng)的系統(tǒng)函數(shù)，輸入序列求系統(tǒng)在變換域的響應(yīng)Y(z)及時(shí)間域的響應(yīng)y(n).MATLAB程序： syms z; X=z./(z-1); H=z.2./(z.2-1.5*z+0.5); Y=X.*H; y=iztrans(Y); 實(shí)驗(yàn)六： 離散系統(tǒng)的描述模型及其轉(zhuǎn)換系統(tǒng)傳遞函數(shù)（tf）模型系統(tǒng)零極點(diǎn)增益（zpk）模型極點(diǎn)留數(shù)（rpk）模型二次分式（sos）模型 1、tf2zp功能：將系統(tǒng)的傳遞函

32、數(shù)（tf）模型轉(zhuǎn)換為系統(tǒng)函數(shù)的零極點(diǎn)增益（zpk）模型調(diào)用格式：z,p,k=tf2zp(num, den)；輸入系統(tǒng)傳遞函數(shù)模型中分子（num）、分母（den）；多項(xiàng)式的系數(shù)向量，求系統(tǒng)函數(shù)的零極點(diǎn)增益模型中的零點(diǎn)向量z,極點(diǎn)向量p和增益系數(shù)k.其中z, p, k.為列量。2、zp2tf功能：將系統(tǒng)函數(shù)的零極點(diǎn)增益（zpk）模型轉(zhuǎn)換為系統(tǒng)的傳遞函數(shù)（tf）模型調(diào)用格式：num, den=zp2tf(z, p, k)；輸入系統(tǒng)函數(shù)的零極點(diǎn)增益模型中的零點(diǎn)向量z,極點(diǎn)向量p和增益系數(shù)k.其中z, p, k.為列量。求系統(tǒng)傳遞函數(shù)模型中分子（num）、分母（den），多項(xiàng)式的系數(shù)向量，。3、tf2s

33、os功能：將系統(tǒng)的傳遞函數(shù)（tf）模型轉(zhuǎn)換為系統(tǒng)函數(shù)的二次分式（sos）模型調(diào)用格式：sos, g=tf2sos(num, den)；輸入系統(tǒng)傳遞函數(shù)模型中分子（num）、分母（den）；多項(xiàng)式的系數(shù)向量，求系統(tǒng)函數(shù)的二次分式模型的系數(shù)矩陣sos，增益系數(shù)g.4、 sos2tf功能：將系統(tǒng)函數(shù)的二次分式（sos）模型轉(zhuǎn)換系統(tǒng)的傳遞函數(shù)（tf）模型調(diào)用格式：num, den=sos2tf(sos, g)；輸入系統(tǒng)函數(shù)的二次分式模型的系數(shù)矩陣sos，增益系數(shù)g（默認(rèn)值為1），求系統(tǒng)傳遞函數(shù)模型中分子（num）、分母（den），多項(xiàng)式的系數(shù)向量。5、sos2zp功能：將系統(tǒng)函數(shù)的二次分式（sos）模

34、型轉(zhuǎn)換系統(tǒng)函數(shù)的零極點(diǎn)增益（zpk）模型調(diào)用格式：z, p ,k=sos2zp(sos, g)；輸入系統(tǒng)函數(shù)的二次分式模型的系數(shù)矩陣sos，增益系數(shù)g（默認(rèn)值為1），求系統(tǒng)函數(shù)的零極點(diǎn)增益模型中的零點(diǎn)向量z,極點(diǎn)向量p和增益系數(shù)k.其中z, p, k.為列量。6、zp2sos功能：將系統(tǒng)函數(shù)的零極點(diǎn)增益（zpk）模型轉(zhuǎn)換為系統(tǒng)函數(shù)的二次分式（sos）模型調(diào)用格式：sos, g=zp2sos(z, p, k); 輸入系統(tǒng)函數(shù)的零極點(diǎn)增益模型中的零點(diǎn)向量z,極點(diǎn)向量p和增益系數(shù)k.其中z, p, k.為列量。，求系統(tǒng)函數(shù)的二次分式模型的系數(shù)矩陣sos，增益系數(shù)g.7、residuez子函數(shù)功能：有

35、理多項(xiàng)式的部分分式展開。調(diào)用格式：r p c=residuez(b, a)；把b(z)/a(z)展開成部分分式的形式。 b , a=residuez(r p c)；根據(jù)部分分式的r p c 數(shù)組，返回有理多項(xiàng)式。其中：b ,a為按降冪排列的多項(xiàng)式的分子和分母的系數(shù)組；r為余數(shù)數(shù)組；p為極點(diǎn)數(shù)組；c為無窮限多項(xiàng)式系數(shù)數(shù)組。例：已知離散時(shí)間系統(tǒng)的傳遞函數(shù) H(z)=10z-1/(1-3z-1+2z-2)求系統(tǒng)的零點(diǎn)向量z，極點(diǎn)向量p和增益系數(shù)k，并列出系統(tǒng)函數(shù)的零-極點(diǎn)增益模型。 MATLAB程序： num=0,10,0; den=1,-3,2; z,p,k=tf2zp(num,den);根據(jù)結(jié)果

36、，零-極點(diǎn)增益模型的系統(tǒng)函數(shù)為： H(z)=10(z/(z-2)(z/(z-1).實(shí)驗(yàn)七： 離散系統(tǒng)的零極點(diǎn)分析1、zplane功能：顯示離散系統(tǒng)的零極點(diǎn)分布圖調(diào)用格式：zplane(z, p)；繪制由列向量z確定的零點(diǎn)，列向量p確定的極點(diǎn)構(gòu)成的零極點(diǎn)分布圖。 zplane(b, a)；繪制由行向量b和a構(gòu)成的系統(tǒng)函數(shù)確定的零極點(diǎn)分布圖。 hz, hp, ht=zplane(z, p)；執(zhí)行后可得到3個(gè)句柄向量，hz為零點(diǎn)線句柄，hp為極點(diǎn)線句柄，ht為坐標(biāo)軸，單位圓及文本對(duì)象的句柄。2、roots功能：求多項(xiàng)式的根調(diào)用格式：r=roots(a)；由多項(xiàng)式的分子或分母系數(shù)向量求根向量。其中，多

37、項(xiàng)式的分子或分母系數(shù)按降冪排列，得到的跟向量為列向量。例（p65）：已知系統(tǒng)的零極點(diǎn)增益模型分別為： ， 求這些系統(tǒng)的零極點(diǎn)分布圖以及系統(tǒng)的沖激響應(yīng)，判斷系統(tǒng)的穩(wěn)定性。結(jié)論：當(dāng)極點(diǎn)處于單位圓內(nèi)，系統(tǒng)的沖激響應(yīng)曲線隨著頻率的增大而收斂；當(dāng)極點(diǎn)處于單位圓上，系統(tǒng)的沖激響應(yīng)曲線為等幅震蕩；當(dāng)極點(diǎn)處于單位圓外，系統(tǒng)的沖激響應(yīng)曲線隨著頻率的增大而發(fā)散。實(shí)驗(yàn)八： 離散系統(tǒng)的頻率響應(yīng)一、離散系統(tǒng)頻率響應(yīng)的基本概念已知穩(wěn)定系統(tǒng)傳遞函數(shù)的零極點(diǎn)增益（zpk）模型為： 則系統(tǒng)的頻率響應(yīng)函數(shù)為： 其中，系統(tǒng)的幅度響應(yīng)特性為： 系統(tǒng)的相位頻響特性為： 說明：系統(tǒng)函數(shù)與頻率響應(yīng)有著緊密的聯(lián)系，適當(dāng)?shù)乜刂葡到y(tǒng)函數(shù)的極點(diǎn)、

38、零點(diǎn)的分布，就可以改變離散系統(tǒng)的頻率響應(yīng)特性。二、離散系統(tǒng)的頻率響應(yīng)子函數(shù) freqz() 功能：用于求離散時(shí)間系統(tǒng)的頻率響應(yīng)函數(shù)。 調(diào)用格式：1）H, w=freqz(b ,a ,n)。可以得到數(shù)字濾波器的n點(diǎn)復(fù)頻響應(yīng)值，這n個(gè)點(diǎn)均勻地分布之0, 上，并將這n個(gè)頻率抽樣點(diǎn)的頻率記錄在w中，對(duì)應(yīng)的幅度值記錄在H中，n缺省時(shí)取512點(diǎn)。 2）h f=freqz(b, a, n, Fs)；用于對(duì)在0, Fs/2上等間隔采樣n點(diǎn)，采樣點(diǎn)頻率及相應(yīng)頻響值分別記錄在f和h中，由用戶指定Fs（以Hz為單位）的值。3）h=freqz(b, a, w)；用于對(duì)在0, 上進(jìn)行采樣，采樣頻率點(diǎn)由矢量w指定。4）h

39、=freqz(b, a, f, Fs)；用于對(duì)在0, Fs上采樣，采樣頻率點(diǎn)由矢量f指定。2、angle( )功能：求相角。調(diào)用格式：p=angle(H)；用于求取復(fù)矢量或復(fù)矩陣H的相角（以弧度為單位），相角介于和之間。3、grid功能：在指定的圖形坐標(biāo)上繪制分格線。調(diào)用格式：grid緊跟在要繪制分格線的繪圖指令后面。例如：plot(t, y)；gridGrid on 繪制分格線。Grid off 不繪制分格線。4、hold功能：在當(dāng)前軸或圖形上多次重疊繪制多條曲線（x軸要一致）。調(diào)用格式：hold 使當(dāng)前圖形具備刷新性質(zhì)的雙向開關(guān)。Hold on 使當(dāng)前軸或圖形保持而不被刷新，準(zhǔn)備接受此后將

40、繪制的新曲線。Hold off 使當(dāng)前軸或圖形不再具備不被刷新的性質(zhì)。5、text功能：在圖形上標(biāo)注文字說明。調(diào)用格式：text(xt, yt, string)；在圖面上（xt, yt）坐標(biāo)處書寫文字說明，其中文字說明字符串必須使用單引號(hào)標(biāo)注。三、舉例 1、已知離散時(shí)間系統(tǒng)的系統(tǒng)函數(shù) 求系統(tǒng)在頻率范圍內(nèi)，歸一化的絕對(duì)幅度頻率響應(yīng)，相對(duì)幅度頻率響應(yīng)，相位頻率響應(yīng)和零極點(diǎn)分部圖。解：MATLAB程序如下： b=0.2, 0.1, 0.3, 0.1, 0.2; a=1, -1.1, 1.5, -0.7, 0.3; n=(0: 500)*pi/500; %在0, 的范圍內(nèi)取501個(gè)采樣點(diǎn) H, w=f

41、reqz(b, a, n); %求系統(tǒng)的頻率響應(yīng) subplot(2, 2, 1), plot(n/pi, abs(H); grid %作系統(tǒng)的絕對(duì)幅度頻響圖 axis(0,max (n/pi ) , 1.1*min(abs(H), 1.1*max(abs(H); ylabel(幅度)；title(幅頻響應(yīng)（V）); subplot(2,2, 2), plot(n/pi, angle(H); grid %作系統(tǒng)的相位頻響圖axis(0, max (n/pi ) , 1.1*min(angle(H), 1.1*max(angle(H); ylabel(相位)；xlabel(以pi為單位的頻率)；

42、title(相頻響應(yīng)); db=20*log10(abs(H); subplot(2, 2, 3), plot(n/pi, abs(H); grid %作系統(tǒng)的相對(duì)幅度頻響圖 title(幅頻響應(yīng)（db）);subplot(2, 2, 4), zplane(b, a); grid %作零極點(diǎn)分布圖title(零極點(diǎn)分布); 四、系統(tǒng)零極點(diǎn)的位置對(duì)系統(tǒng)頻率響應(yīng)的影響系統(tǒng)的零極點(diǎn)的位置對(duì)系統(tǒng)響應(yīng)有著非常明顯的影響。觀察系統(tǒng)極點(diǎn)位置對(duì)幅頻響應(yīng)的影響。已知一階離散系統(tǒng)的傳遞函數(shù)為，假設(shè)系統(tǒng)的零點(diǎn)在原點(diǎn)，極點(diǎn)分別取0.2,0.5,0.8。比較它們的幅頻響應(yīng)曲線，從中了解系統(tǒng)極點(diǎn)的位置對(duì)幅頻響應(yīng)有何影響。

43、實(shí)驗(yàn)八、 時(shí)域抽樣與信號(hào)的重建對(duì)連續(xù)信號(hào)進(jìn)行采樣：已知一個(gè)連續(xù)時(shí)間信號(hào)f(t)=sin(2pif0t)+1/3sin(6pif0t), f0=1Hz,取最高有限帶寬頻率fm=5f0,分別顯示原連續(xù)時(shí)間信號(hào)和Fs>2f0, Fs=2f0, Fs<2f0,三種情況下抽樣信號(hào)的波形。Matlab程序：dt=0.1;f0=1;T0=1/f0;fm=5*f0;Tm=1/fm;t=-2:dt:2;f=sin(2*pi*f0*t)+1/3*sin(6*pi*f0*t); %建立原連續(xù)信號(hào)subplot(4,1,1),plot(t,f);axis(min(t) max(t) 1.1*min(f)

44、1.1*max(f);%title(原連續(xù)信號(hào)和抽樣信號(hào));for i=1:3; fs=i*fm;Ts=1/fs;n=-2:Ts:2;f=sin(2*pi*f0*n)+1/3*sin(6*pi*f0*n);subplot(4,1,i+1),stem(n,f,filled);axis(min(n) max(n) 1.1*min(f) 1.1*max(f);end實(shí)驗(yàn)九（1）： 模擬原型濾波器的設(shè)計(jì)實(shí)驗(yàn)?zāi)康模海?）加深對(duì)模擬濾波器基本類型，特點(diǎn)和主要設(shè)計(jì)指標(biāo)的了解（2）掌握模擬低通濾波器原型的設(shè)計(jì)方法（3）學(xué)習(xí)MATLAB語言有關(guān)模擬原型濾波器設(shè)計(jì)的子函數(shù)的使用方法MATLAB子函數(shù)：1、butt

45、ord功能：確定巴特沃思（Butterworth）濾波器的階數(shù)和3分貝截止頻率調(diào)用格式：n, wn=buttord(wp, ws, Rp, As)；計(jì)算巴特沃思數(shù)字濾波器的階數(shù)和3分貝截止頻率。其中，其值為1時(shí)表示0.5Fs，Rp為通帶最大衰減指標(biāo)，As為阻帶最小衰減指標(biāo)。 n, wn=buttord(wp, ws, Rp, As，s)；計(jì)算巴特沃思模擬濾波器的階數(shù)和3分貝截止頻率。wp, ws,可以是實(shí)際的頻率值或角頻率值，Wn 將取相同的量綱。Rp為通帶最大衰減指標(biāo)，As為阻帶最小衰減指標(biāo)。 Wp>ws時(shí)，為高通濾波器，當(dāng)wp, ws為二元向量時(shí)，為帶通或帶阻濾波器，此時(shí)Wn也為二元

46、向量。2、cheb1ord功能：確定切比雪夫I型濾波器的階數(shù)和通帶截止頻率調(diào)用格式：n, wn=cheblord(wp, ws, Rp, As )；計(jì)算切比雪夫I型數(shù)字濾波器的階數(shù)和通帶截止頻率。其中，其值為1時(shí)表示0.5Fs，Rp為通帶最大衰減指標(biāo)，As為阻帶最小衰減指標(biāo)。n, wn= cheblord (wp, ws, Rp, As，s)；計(jì)算切比雪夫I型模擬濾波器的階數(shù)和通帶截止頻率。wp, ws,可以是實(shí)際的頻率值或角頻率值，Wn 將取相同的量綱。Rp為通帶最大衰減指標(biāo)，As為阻帶最小衰減指標(biāo)。 Wp>ws時(shí)，為高通濾波器，當(dāng)wp, ws為二元向量時(shí)，為帶通或帶阻濾波器，此時(shí)Wn

47、也為二元向量。3、cheb2ord功能：確定切比雪夫II型濾波器的階數(shù)和阻帶截止頻率調(diào)用格式：n, wn=cheb2ord(wp, ws, Rp, As )；計(jì)算切比雪夫II型數(shù)字濾波器的階數(shù)和通帶截止頻率。其中，其值為1時(shí)表示0.5Fs，Rp為通帶最大衰減指標(biāo)，As為阻帶最小衰減指標(biāo)。n, wn= cheb2ord (wp, ws, Rp, As，s)；計(jì)算切比雪夫II型模擬濾波器的階數(shù)和通帶截止頻率。wp, ws,可以是實(shí)際的頻率值或角頻率值，Wn 將取相同的量綱。Rp為通帶最大衰減指標(biāo)，As為阻帶最小衰減指標(biāo)。 Wp>ws時(shí)，為高通濾波器，當(dāng)wp, ws為二元向量時(shí)，為帶通或帶阻濾

48、波器，此時(shí)Wn也為二元向量。4、ellipord功能：確定橢圓濾波器的階數(shù)和通帶截止頻率調(diào)用格式：n, wn=ellipord(wp, ws, Rp, As)；計(jì)算橢圓數(shù)字濾波器的階數(shù)和通帶截止頻率。其中，其值為1時(shí)表示0.5Fs，Rp為通帶最大衰減指標(biāo)，As為阻帶最小衰減指標(biāo)。n, wn=ellipord(wp, ws, Rp, As，s)；計(jì)算橢圓模擬濾波器的階數(shù)和通帶截止頻率。wp, ws,可以是實(shí)際的頻率值或角頻率值，Wn 將取相同的量綱。Rp為通帶最大衰減指標(biāo)，As為阻帶最小衰減指標(biāo)。 Wp>ws時(shí)，為高通濾波器，當(dāng)wp, ws為二元向量時(shí)，為帶通或帶阻濾波器，此時(shí)Wn也為二元

49、向量。5、buttap功能：確定巴特沃思（Butterworth）模擬低通濾波器原型調(diào)用格式： z, p, k=buttap(n)；設(shè)計(jì)巴特沃思模擬低通濾波器原型，其傳遞函數(shù)為： 此時(shí)z為控陣，巴特沃思濾波器由通帶內(nèi)最平坦，總體上單調(diào)的幅度特性來表征。6、cheb1ap功能：切比雪夫I型模擬低通濾波器原型調(diào)用格式： z, p, k=cheb1ap(n, Rp)；設(shè)計(jì)切比雪夫I型模擬低通濾波器原型，其通帶內(nèi)的波紋系數(shù)為Rp分貝，傳遞函數(shù)為： 此時(shí)z為控陣。切比雪夫I型濾波器是通帶內(nèi)等波紋，阻帶內(nèi)單調(diào)的濾波器，其極點(diǎn)均為分布在左半平面的橢圓上。7、cheb2ap功能：切比雪夫II型模擬低通濾波器原

50、型調(diào)用格式： z, p, k=cheb2ap(n, As)；設(shè)計(jì)切比雪夫II型模擬低通濾波器原型，其阻帶內(nèi)的波紋系數(shù)小于As分貝，傳遞函數(shù)為： 切比雪夫II型濾波器是通帶內(nèi)單調(diào)，阻帶內(nèi)等波紋的濾波器，其極點(diǎn)位置為cheb1ap極點(diǎn)位置的倒數(shù)。8、ellipap功能：橢圓模擬低通濾波器原型調(diào)用格式： z, p, k=ellipap(n, Rp, As)；設(shè)計(jì)橢圓模擬低通濾波器原型，其通帶內(nèi)的波紋系數(shù)為Rp分貝，其阻帶內(nèi)的波紋系數(shù)小于As分貝。其傳遞函數(shù)為： 橢圓濾波器是通帶和阻帶內(nèi)均是等波紋的濾波器，它是有比巴特沃斯和切比雪夫更陡的下降斜率，但會(huì)損失通帶和阻帶的波紋指標(biāo)。9、poly功能：求某根

51、向量所對(duì)應(yīng)的特征多項(xiàng)式調(diào)用格式：P=poly()；求根向量為的特征多項(xiàng)式，產(chǎn)生多項(xiàng)式系數(shù)向量P。例如：降冪多項(xiàng)式P(x)=其系數(shù)行向量表達(dá)式為： 若要表示()()-()=，可建立，再利用指令：。多項(xiàng)式P是一個(gè)特征多項(xiàng)式，的元素被認(rèn)為是多項(xiàng)式P的根。10、poly2str功能：以習(xí)慣方式顯示多項(xiàng)式。調(diào)用格式： Pa=poly2str(a,s)；以習(xí)慣方式顯示s的多項(xiàng)式。例：輸入程序：A=2 ; 6;7 PA=poly(A) PPA=poly2str(PA,s) 得到：PA=1 -15 68 -84 PPA=s3 - 15 s2 + 68 s 8411、pzmap功能：顯示連續(xù)系統(tǒng)的零極點(diǎn)分布圖。

52、調(diào)用格式： pzmap(b, a)；繪制由行向量b和a構(gòu)成的系統(tǒng)的系統(tǒng)函數(shù)確定的零極點(diǎn)分布圖。例6-4：通過模擬濾波器原型設(shè)計(jì)一個(gè)巴特沃思模擬低通濾波器的系統(tǒng)函數(shù)，要求通帶截止頻率，通帶最大衰減，阻帶截止頻率，阻帶最小衰減解：程序如下：>> fp=4000;omgp=2*pi*fp;>> fs=8000;omgs=2*pi*fs;>> Rp=3;As=20;>> n,omgc=buttord(omgp,omgs,Rp,As,'s');>> %計(jì)算n階模擬低通濾波器原型，得到左半平面零極點(diǎn)>> z0,p0,k0=buttap(n);>> b0=k0*real(poly(z0);>> a0=real(poly(p0);%歸一化4階巴特沃斯低通濾波器，其分母多項(xiàng)式為>> ppa=poly2str(a0,'s');ppb=poly2str(b0,'s');例6-6 給定模擬低通濾波器的性能指標(biāo)為，在通帶內(nèi)，即在內(nèi)，幅度函數(shù)的波紋（起伏）,在阻帶內(nèi)，時(shí)，幅度函數(shù)衰減。試求用切貝雪夫?yàn)V波器