1、實驗報告2014-2015學年第二學期開 課 單 位 適用年級、專業 13級所有專業 課 程 名 稱 數字信號處理實驗 主 講 教 師 課 程 序 號 課 程 代 碼 實 驗 名 稱 實 驗 學 時 學 號 姓 名 姓名： 學實驗1 用MATLAB產生時域離散信號一、.實驗目的：1、了解常用時域離散信號及其特點。2、掌握用MATLAB產生時域離散信號的方法。二、.實驗原理：1、時域離散信號的概念在時間軸的離散點上取值的信號，稱為離散時間信號。通常，離散時間信號用x(n)表示，其幅度可以在某一范圍內連續取值。由于信號處理設備或裝置（如計算機、專用的信號處理芯片等）均以有限位的二進制數來表示信號的

2、幅度，因此，信號的幅度也必須離散化。我們把時間和幅度均取離散值的信號稱為時域離散信號或數字信號。在MATLAB語言中，時域離散信號可以通過編寫程序直接產生。三、.實驗內容：1、閱讀并上機驗證實驗原理部分的例題程序，理解每一條語句的含義。改變例題中的有關參數（如信號的頻率、周期、幅度、顯示時間的取值范圍、采樣點數等），觀察對信號波形的影響。2、編寫程序，產生以下離散序列：（1）f(n)=(n) (-3n4)編寫程序：n1=-3;n2=4;n0=0;n=n1:n2;x=n=n0; stem(n,x, filled);axis(n1,n2,0,1.1*max(x);xlabel(時間(n);ylab

3、el(幅度x(n);title(單位脈沖序列);程序運行結果如圖：（2）f(n)=u(n) (-5n5)編寫程序：n1=-5;n2=5;n0=0;n=n1:n2;x=n=n0;stem(n,x,filled);axis(n1,n2,0,1.1*max(x);xlabel(時間(n);ylabel(幅度 (n);title(單位階躍序列);box程序運行結果如圖：（3）f(n)= e(0.1+j1.6)n (0n16)編寫程序：n1=16;a=0.1;w=1.6;n=0:n1;x=exp(a+j*w)*n);subplot(2,2,1);plot(n,real(x);title(復指數信號的實部

4、);subplot(2,2,3);stem(n,real(x),filled);title(復指數序列的實部);subplot(2,2,2);plot(n,imag(x);title(復指數信號的虛部);subplot(2,2,4);stem(n,imag(x),filled);title(復指數序列的虛部);box程序運行結果如圖：（4）f(n)=3sin(n/4) (0n20)編寫程序：f=1/8;Um=3;nt=2.5;T=1/f;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.

5、1*max(x);ylabel(x(t);ylabel(x(n);box 程序運行結果如圖：3、一個連續的周期性方波信號頻率為200Hz，信號幅度在-1+1V之間，要求在圖形窗口上顯示其兩個周期的波形。以4kHz的頻率對連續信號進行采樣，編寫程序生成連續信號和其采樣獲得的離散信號波形。編寫程序：f=200;nt=2;N=20;T=1/f;dt=T/N;n=0:nt*N-1;tn=n*dt;x=square(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(t);subplot(2,1,

6、2);stem(tn,x);axis(0,nt*T,1.1*min(x),1.1*max(x);ylabel(x(n);box程序運行結果如圖：四、實驗總結1、通過本次實驗，我學會了MATLAB產生時域離散信號波形。使用MATLAB這個軟件和編寫程序的一些基本知識，同時還加深了我對各種時域離散信號公式與波形的記憶和理解。思考題：1、通過例題程序，你發現采樣頻率Fs、采樣點數N、采樣時間間隔dt在程序編寫中有怎樣的聯系？使用時需注意什么問題？答： dt=T/N N=Fs/f T=1/f 實驗2 離散LSI系統的時域分析.實驗目的：1、加深對離散系統的差分方程、單位脈沖響應、單位階躍響應和卷積分析

7、方法的理解。2、初步了解用MATLAB語言進行離散時間系統時域分析的基本方法。3、掌握求解離散時間系統的單位脈沖響應、單位階躍響應、線性卷積以及差分方程的程序的編寫方法，了解常用子函數的調用格式。二、實驗原理：1、離散LSI系統的響應與激勵由離散時間系統的時域分析方法可知，一個離散LSI系統的響應與激勵可以用如下框圖表示：其輸入、輸出關系可用以下差分方程描述：三、實驗內容：1、輸入并運行例題程序，理解每一條語句的含義。2、已知描述某離散LSI系統的差分方程為2y(n)-3y(n-1)+y(n-2)=x(n-1)，分別用impz和dstep函數、filtic和filter函數兩種方法求解系統的單

8、位序列響應和單位階躍響應。用impz和dstep函數求解程序清單如下：a=1,-3/2,1/2;b=0,1/2,0;N=32;n=0:N-1;hn=impz(b,a,n);gn=dstep(b,a,n);subplot(1,2,1);stem(n,hn,k);title(系統的單位序列響應);ylabel(h(n);xlabel(n);axis(0,N,1.1*min(hn),1.1*max(hn);subplot(1,2,2);stem(n,gn,k);title(系統的單位階躍響應);ylabel(g(n);xlabel(n);axis(0,N,1.1*min(gn),1.1*max(gn

9、);filtic和filter函數求解程序清單如下：a=1,-3/2,1/2;b=0,1/2,0;N=32;n=0:N-1;xi=filtic(b,a,0);x1=n=0;hn=filter(b,a,x1,xi);x2=n=0;gn=filter(b,a,x2,xi);subplot(1,2,1);stem(n,hn,k);title(系統的單位序列響應);ylabel(h(n);xlabel(n);axis(0,N,1.1*min(hn),1.1*max(hn);subplot(1,2,2);stem(n,gn,k);title(系統的單位階躍響應);ylabel(g(n);xlabel(n

10、);axis(0,N,1.1*min(gn),1.1*max(gn);3、編寫程序描繪下列序列的卷積波形：（1）f1(n=u(n),f2(n)=u(n-2), (0n10) 程序清單如下：n1=0:10;N1=length(n1);f1=ones(1,N1);subplot(2,2,1);stem(n1,f1,filled);title(f1(n);n2=2:12;N2=length(n2);f2=ones(1,N2);subplot(2,2,2);stem(n2,f2,filled);title(f2(n);y=conv(f1,f2);subplot(2,1,2);stem(y,filled

11、);x(n)=sin(n/2),h(n)=(0.5)n (-3n4) 程序清單如下：n1=-3:4*pi;f1=sin(n1/2);n2=-3:4*pi;f2=0.5.n2;nys=n1(1)+n2(1);nyf=n1(end)+n2(end);y=conv(f1,f2);ny=nys:nyf;subplot(2,2,1);stem(n1,f1);subplot(2,2,2);stem(n2,f2);subplot(2,1,2);stem(ny,y);4、已知某離散LSI系統的單位序列響應為h(n)=3(n-3)+0.5(n-4)+0.2(n-5)+0.7(n-6)-0.8(n-7)求輸入為x

12、(n)=e-0.5nu(n)時的系統響應。程序清單如下:N=16;n=0:N-1;x=exp(-0.5*n);subplot(2,2,1);stem(n,x);a=1;b=0,0,0,3,0.5,0.2,0.7,-0.8;hn=impz(b,a,n);subplot(2,2,2);stem(n,hn);y=conv(x,hn);subplot(2,1,2);stem(y);程序運行結果如下圖：已知描述某離散LSI系統的差分方程為y(n)=0.7y(n-1)+2x(n)-x(n-2)，求輸入為x(n)=u(n-3)時的系統響應。 程序清單如下:N=16;n1=3:N+2;f1=zeros(1,3

13、),ones(1,(N-3);subplot(2,2,1);stem(n1,f1);a=1,-0.7,0;b=2,0,-1;f2=impz(b,a,n1);subplot(2,2,2);stem(n1,f2);nys=2*n1(1);nyf=2*n1(end);y=conv(f1,f2);ny=nys:nyf;subplot(2,1,2);stem(ny,y);程序運行結果如下圖： 實驗3 離散LSI系統的頻域分析一、實驗目的：1、加深對離散系統變換域分析z變換的理解，掌握使用MATLAB進行z變換和逆z變換的常用函數的用法。2、了解離散系統的零極點與系統因果性和穩定性的關系，熟悉使用MATL

14、AB進行離散系統的零極點分析的常用函數的用法。3、加深對離散系統的頻率響應特性基本概念的理解，掌握使用MATLAB進行離散系統幅頻響應和相頻響應特性分析的常用方法。二、實驗原理：1、z變換和逆z變換（1）用ztrans函數求無限長序列的z變換。該函數只給出z變換的表達式，而沒有給出收斂域。另外，由于這一函數還不盡完善，有的序列的z變換還不能求出，逆z變換也存在同樣的問題。例3-1 求以下各序列的z變換 x1(n)=an x2(n)=n x3(n)=n(n-1)/2 x4(n)=ejon x5(n)=1/n(n-1)程序清單如下：syms w0 n z a;x1=an;X1=ztrans(x1)

15、x2=n;X2=ztrans(x2)x3=(n*(n-1)/2;X3=ztrans(x3)x4=exp(j*w0*n);X4=ztrans(x4)x5=1/n*(n-1);X5=ztrans(x5)程序運行結果如下：X1 =z/a/(z/a-1) X2 =z/(z-1)2X3 =1/2*z*(z+1)/(z-1)3-1/2*z/(z-1)2X4 =z/exp(i*w0)/(z/exp(i*w0)-1)X5 =z/(z-1)-ztrans(1/n,n,z)（2）用iztrans函數求無限長序列的逆z變換。例3-2 求下列函數的逆z變換。 程序清單如下：syms n z a;X1=z/(z-1);

16、x1=iztrans(X1)X2=a*z/(a-z)2;x2=iztrans(X2)X3=z/(z-1)3;x3=iztrans(X3)X4=(1-z-n)/(1-z-1);x4=iztrans(X4)程序運行結果如下：x1 =1x2 =an*nx3 =1/2*n2-1/2*nx4 =iztrans(1-z(-n)/(1-1/z),z,n)2、離散系統的零極點分析（系統極點位置對系統響應的影響）例3-3 研究z右半平面的實數極點對系統的影響。已知系統的零極點增益模型分別為：求這些系統的零極點分布圖以及系統的單位序列響應，判斷系統的穩定性。程序清單如下：z1=0;p1=0.85;k=1;b1,a

17、1=zp2tf(z1,p1,k);subplot(3,2,1);zplane(z1,p1);title(極點在單位圓內);subplot(3,2,2);impz(b1,a1,20);z2=0;p2=1;b2,a2=zp2tf(z2,p2,k);subplot(3,2,3);zplane(z2,p2);title(極點在單位圓上);subplot(3,2,4);impz(b2,a2,20);z3=0;p3=1.5;b3,a3=zp2tf(z3,p3,k);subplot(3,2,5);zplane(z3,p3);title(極點在單位圓外);subplot(3,2,6);impz(b3,a3,2

18、0);程序運行結果如圖3-1所示。由圖可見，這三個系統的極點均為實數且處于z平面的右半平面。由圖可知，當極點位于單位圓內，系統的單位序列響應隨著頻率的增大而收斂；當極點位于單位圓上，系統的單位序列響應為等幅振蕩；當極點位于單位圓外，系統的單位序列響應隨著頻率的增大而發散。由此可知系統1、2為穩定系統。圖3-1例3-4 研究z左半平面的實數極點對系統的影響。已知系統的零極點增益模型分別為：求這些系統的零極點分布圖以及系統的單位序列響應，判斷系統的穩定性。程序清單如下：z1=0;p1=-0.85;k=1;b1,a1=zp2tf(z1,p1,k);subplot(3,2,1);zplane(z1,p

19、1);title(極點在單位圓內);subplot(3,2,2);impz(b1,a1,20);z2=0;p2=-1;b2,a2=zp2tf(z2,p2,k);subplot(3,2,3);zplane(z2,p2);title(極點在單位圓上);subplot(3,2,4);impz(b2,a2,20);z3=0;p3=-1.5;b3,a3=zp2tf(z3,p3,k);subplot(3,2,5);zplane(z3,p3);title(極點在單位圓外);subplot(3,2,6);impz(b3,a3,20);程序運行結果如圖3-2所示。由圖可見，這三個系統的極點均為實數且處于z平面的

20、左半平面。由圖可知，當極點位于單位圓內，系統的單位序列響應隨著頻率的增大而收斂；當極點位于單位圓上，系統的單位序列響應為等幅振蕩；當極點位于單位圓外，系統的單位序列響應隨著頻率的增大而發散。由此可知系統1、2為穩定系統。圖3-2例3-5 研究z右半平面的復數極點對系統響應的影響已知系統的零極點增益模型分別為：求這些系統的零極點分布圖以及系統的單位序列響應，判斷系統的穩定性。程序清單如下：z1=0.3,0;p1=0.5+0.7j,0.5-0.7j;k=1;b1,a1=zp2tf(z1,p1,k);subplot(3,2,1);zplane(z1,p1);title(極點在單位圓內);subplo

21、t(3,2,2);impz(b1,a1,20);z2=0.3,0;p2=0.6+0.8j,0.6-0.8j;b2,a2=zp2tf(z2,p2,k);subplot(3,2,3);zplane(z2,p2);title(極點在單位圓上);subplot(3,2,4);impz(b2,a2,20);z3=0.3,0;p3=1+j,1-j;b3,a3=zp2tf(z3,p3,k);subplot(3,2,5);zplane(z3,p3);title(極點在單位圓外);subplot(3,2,6);impz(b3,a3,20);程序運行結果如圖3-3所示。由圖可見，這三個系統的極點均為復數且處于z平

22、面的右半平面。由圖可知，當極點位于單位圓內，系統的單位序列響應隨著頻率的增大而收斂；當極點位于單位圓上，系統的單位序列響應為等幅振蕩；當極點位于單位圓外，系統的單位序列響應隨著頻率的增大而發散。由此可知系統1、2為穩定系統。圖3-3由以上三例可得結論：系統只有在極點處于單位圓內才是穩定的。例3-6 已知某離散時間系統的系統函數為求該系統的零極點及零極點分布圖，并判斷系統的因果穩定性。程序清單如下：b=0.2,0.1,0.3,0.1,0.2;a=1,-1.1,1.5,-0.7,0.3;rz=roots(b)rp=roots(a)subplot(2,1,1);zplane(b,a);title(系

23、統的零極點分布圖);subplot(2,1,2);impz(b,a,20);title(系統的單位序列響應);xlabel(n);ylabel(h(n);程序運行結果如下：rz = -0.5000 + 0.8660i -0.5000 - 0.8660i 0.2500 + 0.9682i 0.2500 - 0.9682irp = 0.2367 + 0.8915i 0.2367 - 0.8915i 0.3133 + 0.5045i 0.3133 - 0.5045i圖3-4由零極點分布圖可見，該系統的所有極點均在單位圓內，因此該系統是一個因果穩定系統。3、離散系統的頻率響應（1）離散系統的頻率響應的

24、基本概念已知穩定系統傳遞函數的零極點增益模型為則系統的頻響函數為其中，系統的幅頻特性為系統的相頻特性為由以上各式可見，系統函數與頻率響應有著密切的聯系。適當地控制系統函數的零極點分布，可以改變離散系統的頻響特性： 在原點（z=0）處的零點或極點至單位圓的距離始終保持不變，其值|ej|=1，所以，對幅度響應不起作用； 單位圓附近的零點對系統幅度響應的谷值位置及深度有明顯影響； 單位圓內且靠近單位圓附近的極點對系統幅度的峰值位置及大小有明顯的影響。（2）系統的頻響特性分析例3-7 已知某離散時間系統的系統函數為求該系統在0頻率范圍內的相對幅頻響應與相頻響應。程序清單如下：b=0.1321,0,-0

25、.3963,0,0.3963,0,-0.1321;a=1,0,0.34319,0,0.60439,0,0.20407;freqz(b,a);程序運行結果如圖3-5所示。該系統是一個IIR數字帶通濾波器。其中幅頻特性采用歸一化的相對幅度值，以分貝（dB）為單位。圖3-5例3-8 已知某離散時間系統的系統函數為求該系統在0頻率范圍內的絕對幅頻響應與相頻響應。程序清單如下：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;h,w=freqz(b,a,n);subplot(2,1,1);plot(n/pi,abs(h);grid;

26、axis(0,1,1.1*min(abs(h),1.1*max(abs(h);xlabel(omega/pi);ylabel(幅度);subplot(2,1,2);plot(n/pi,angle(h);grid;axis(0,1,1.1*min(angle(h),1.1*max(angle(h);xlabel(omega/pi);ylabel(相位);程序運行結果如圖3-6所示。該系統為一低通濾波器。圖3-6例3-9 已知某離散時間系統的系統函數為求該系統在0頻率范圍內的絕對幅頻響應與相頻響應、相對幅頻響應與相頻響應及零極點分布圖。程序清單如下：b=0.1,-0.4,0.4,-0.1;a=1,

27、0.3,0.55,0.2;n=(0:500)*pi/500;h,w=freqz(b,a,n);db=20*log10(abs(h);subplot(2,2,1);plot(w/pi,abs(h);grid;axis(0,1,1.1*min(abs(h),1.1*max(abs(h);title(幅頻特性（V）);xlabel(omega/pi);ylabel(幅度(V);subplot(2,2,2);plot(w/pi,angle(h);grid;axis(0,1,1.1*min(angle(h),1.1*max(angle(h);xlabel(omega/pi);ylabel(相位);tit

28、le(相頻特性);subplot(2,2,3);plot(w/pi,db);gridaxis(0,1,-100,5);title(幅頻特性（dB）);subplot(2,2,4);zplane(b,a);title(零極點分布);程序運行結果如圖3-7所示：圖3-7（3）一個求解頻率響應的實用函數。在實際使用freqz進行離散系統頻響特性分析時。通常需要求解幅頻響應、相頻響應、群時延，幅頻響應又分為絕對幅頻和相對幅頻兩種表示方法。下面定義函數freqz_m，利用該函數，可方便求出上述各項。freqz_m函數定義如下：functiondb,mag,pha,grd,w=freqz_m(b,a);H

29、,w=freqz(b,a,1000,whole);H=(H(1:501);w=(w(1:501);mag=abs(H);db=20*log10(mag+eps)/max(mag);pha=angle(H);grd=grpdelay(b,a,w);例3-10 已知某離散時間系統的系統函數為求該系統在0頻率范圍內的絕對幅頻響應與相頻響應、相對幅頻響應與相頻響應及群時延。程序清單如下：b=0.1321,0,0.3963,0,0.3963,0,0.1321;a=1,0,-0.34319,0,0.60439,0,-0.20407;db,mag,pha,grd,w=freqz_m(b,a);subplot

30、(2,2,1);plot(w/pi,mag);gridaxis(0,1,1.1*min(mag),1.1*max(mag);title(幅頻特性（V）);xlabel(omega/pi);ylabel(幅度(V);subplot(2,2,2);plot(w/pi,pha);grid;axis(0,1,1.1*min(pha),1.1*max(pha);xlabel(omega/pi);ylabel(相位);title(相頻特性);subplot(2,2,3);plot(w/pi,db);gridaxis(0,1,-100,5);title(幅頻特性（dB）);subplot(2,2,4);pl

31、ot(w/pi,grd);gridaxis(0,1,0,10)title(群時延);程序運行結果如圖3-8所示：圖3-8三、實驗內容：1、輸入并運行例題程序，理解每條語句的含義。2、求下列各序列的z變換：程序清單如下：symsw0na;x1=n*an;X1=ztrans(x1)x2=sin(w0*n);X2=ztrans(x2)x3=exp(-a*n)*sin(w0*n);X3=ztrans(x3)程序運行結果如下：X1=z*a/(-z+a)2X2=-sin(w0)*z/(-z2+2*z*cos(w0)-1)X3=-sin(w0)*z/exp(-a)/(-z2/exp(-a)2+2*z/exp

32、(-a)*cos(w0)-1)3、求下列函數的逆z變換程序清單如下：symsw0nza;X1=z/(z-a);x1=iztrans(X1)X2=z/(a-z)2;x2=iztrans(X2)X3=z/z-exp(j*w0);x3=iztrans(X3)X4=(1-z-3)/(1-z-1);x4=iztrans(X4)程序運行結果如下：x1=anx2=n*an/ax3=charfcn0(n)-iztrans(exp(i*w0),w0,n)x4=charfcn2(n)+charfcn1(n)+charfcn0(n)4、求下列系統函數所描述的離散系統的零極點分布圖，并判斷系統的穩定性（1）z1=-0

33、.3;p1=-1+j,-1-j;k=1;b1,a1=zp2tf(z1,p1,k);subplot(3,2,1);zplane(z1,p1);title(極點在單位圓內);subplot(3,2,2);impz(b1,a1,20);由圖可見當極點位于單位圓內，系統的單位序列響應隨著頻率的增大而收斂；當極點位于單位圓上，系統的單位序列響應為等幅振蕩；當極點位于單位圓外，系統的單位序列響應隨著頻率的增大而發散。由此可知系統為穩定系統。（2）b=4,-1.6,-1.6,4;a=1,0.4,0.35,-0.4;rz=roots(b)rp=roots(a)subplot(2,1,1);zplane(b,a

34、);title(系統的零極點分布圖);subplot(2,1,2);impz(b,a,20);title(系統的單位序列響應);xlabel(n);ylabel(h(n);rz=-1.00000.7000+0.7141i0.7000-0.7141irp=-0.4500+0.7730i-0.4500-0.7730i0.5000由零極點分布圖可見，該系統的所有極點均在單位圓內，因此該系統是一個因果穩定系統。5、已知某離散時間系統的系統函數為求該系統在0頻率范圍內的絕對幅頻響應與相頻響應、相對幅頻響應與相頻響應及群時延。b=0.187632,0,-0.241242,0,0.241242,0,-187

35、632;a=1,0,0.602012,0,0.495684,0,0.035924;db,mag,pha,grd,w=freqz_m(b,a);subplot(2,2,1);plot(w/pi,mag);gridaxis(0,1,1.1*min(mag),1.1*max(mag);title(幅頻特性（V）);xlabel(omega/pi);ylabel(幅度(V);subplot(2,2,2);plot(w/pi,pha);grid;axis(0,1,1.1*min(pha),1.1*max(pha);xlabel(omega/pi);ylabel(相位);title(相頻特性);subpl

36、ot(2,2,3);plot(w/pi,db);gridaxis(0,1,-100,5);title(幅頻特性（dB）);subplot(2,2,4);plot(w/pi,grd);gridaxis(0,1,0,10);title(群時延);四、實驗預習：1、認真閱讀實驗原理部分，明確實驗目的，復習有關離散LSI系統頻率響應的理論知識。2、讀懂實驗原理部分的例題程序，熟悉與本實驗有關的MATLAB函數。3、根據實驗內容預先編寫實驗程序，并思考本實驗提出的有關MATLAB函數在調用時應注意哪些問題。4、預習思考題： 系統函數零極點的位置與系統單位序列響應有何關系？ 離散系統的零極點對系統幅頻響應

37、有何影響？答：系統函數零極點的位置與系統單位序列響應有何關系？當極點位于單位圓內，系統的單位序列響應隨著頻率的增大而收斂；當極點位于單位圓上，系統的單位序列響應為等幅振蕩；當極點位于單位圓外，系統的單位序列響應隨著頻率的增大而發散。離散系統的零極點對系統幅頻響應有何影響？在原點（z=0）處的零點或極點至單位圓的距離始終保持不變，其值|ej|=1，所以，對幅度響應不起作用；單位圓附近的零點對系統幅度響應的谷值位置及深度有明顯影響；單位圓內且靠近單位圓附近的極點對系統幅度的峰值位置及大小有明顯的影響。實驗4 DFS、DFT與FFT一、實驗目的：1、加深對周期序列DFS、有限長序列DFT和FFT的基

38、本概念及理論的理解。2、掌握用MATLAB語言求解DFS、DFT、FFT以及相應反變換的方法。3、觀察周期序列的重復周期數對序列頻譜特性的影響。 二、實驗原理：1、周期序列的離散傅里葉級數（DFS）（1）DFS的基本概念：離散時間序列x(n)滿足x(n)= x(n+rN)，稱為離散周期序列，用表示。其中，N為信號的周期，x(n)稱為離散周期序列的主值。周期序列可以用離散傅里葉級數（DFS）表示：其中，是周期序列DFS第k次諧波分量的系數，也稱為周期序列的頻譜，可表示為以上兩式也是周期序列的一對傅里葉級數變換對。令，以上DFS變換對又可以寫成：與連續周期信號的傅里葉級數相比，周期序列的離散傅里葉

39、級數有以下特點： 連續周期信號的傅里葉級數由無窮多個與基波頻率成整數倍的諧波分量疊加而成，而周期為N的周期序列的傅里葉級數僅有N個獨立的諧波分量。 周期序列的頻譜也是一個以N為周期的周期序列。（2）周期序列的DFS和IDFS例4-1 已知一個周期性矩形脈沖寬度占整個周期的1/4，一個周期的采樣點數為16點，編程顯示3個周期的序列波形，并： 用傅里葉級數求信號的幅度譜和相位譜。 求傅里葉級數逆變換的圖形，并與原序列進行比較。程序清單如下：N=16;xn=ones(1,N/4),zeros(1,3*N/4);xn=xn,xn,xn;n=0:3*N-1;k=0:3*N-1;Xk=xn*exp(-j*

40、2*pi/N).(n*k);x=(Xk*exp(j*2*pi/N).(n*k)/N;subplot(2,2,1);stem(n,xn);title(x(n);axis(-1,3*N,1.1*min(xn),1.1*max(xn);subplot(2,2,2);stem(n,abs(x);title(IDFS|X(k)|);axis(-1,3*N,1.1*min(x),1.1*max(x);subplot(2,2,3),stem(k,abs(Xk);title(|X(k)|);axis(-1,3*N,1.1*min(abs(Xk),1.1*max(abs(Xk);subplot(2,2,4),s

41、tem(k,angle(Xk);title(arg|X(k)|);axis(-1,3*N,1.1*min(angle(Xk),1.1*max(angle(Xk);程序運行結果如圖4-1所示：圖4-1由離散傅里葉級數逆變換圖形可見，與原序列相比，幅度擴大了32倍。這是因為周期序列為原主值序列周期的3倍，做逆變換時未做處理。可將逆變換程序改為：x=(Xk*exp(j*2*pi/N).(n*k)/3*3*N;由上例可見，周期序列的DFS和IDFS是依據變換公式編程的，無論信號序列如何變化，求解的公式總是一樣的。因此，可將其編寫成通用子程序： 離散傅里葉級數變換通用子程序dfs.mfunctionXk

42、=dfs(xn,N)n=0:N-1;k=0:N-1;WN=exp(-j*2*pi/N);nk=n*k;Xk=xn*WN.nk; 離散傅里葉級數逆變換通用子程序idfs.mfunctionxn=idfs(Xk,N)n=0:N-1;k=0:N-1;WN=exp(j*2*pi/N);nk=n*k;xn=(Xk*WN.nk)/N;例4-2 利用上述兩個子程序，重做例4-1程序清單如下：N=16;xn=ones(1,N/4),zeros(1,3*N/4);n=0:N-1;k=0:N-1;Xk=dfs(xn,N);x=idfs(Xk,N);subplot(2,2,1);stem(n,xn);title(x

43、(n);axis(-1,3*N,1.1*min(xn),1.1*max(xn);subplot(2,2,2);stem(n,abs(x);title(IDFS|X(k)|);axis(-1,3*N,1.1*min(x),1.1*max(x);subplot(2,2,3),stem(k,abs(Xk);title(|X(k)|);axis(-1,3*N,1.1*min(abs(Xk),1.1*max(abs(Xk);subplot(2,2,4),stem(k,angle(Xk);title(arg|X(k)|);axis(-1,3*N,1.1*min(angle(Xk),1.1*max(angl

44、e(Xk);程序運行結果如圖4-2所示。由于子程序僅適用于對主值區間進行變換，周期次數無法傳遞給子程序，因此程序執行結果僅顯示一個周期的變換情況。圖4-2（2）周期重復次數對序列頻譜的影響理論上講，周期序列不滿足絕對可積條件，因此不能用傅里葉級數來表示。實際處理時可先取K個周期進行處理，然后令K趨于無窮大，分析其極限情況。根據這一分析思路，可以觀察序列由非周期到周期變化時，頻譜由連續譜逐漸向離散譜過渡的過程。例4-3 已知一矩形脈沖寬度占整個周期的1/2，一個周期的采樣點數為10點，用DFS求序列的重復周期數分別為1、4、7、10時的幅度譜。程序清單如下：xn=ones(1,5),zeros(

45、1,5);Nx=length(xn);Nw=1000;dw=2*pi/Nw;k=floor(-Nw/2+0.5):(Nw/2+0.5);for r=0:3 K=3*r+1; nx=0:(K*Nx-1); x=xn(mod(nx,Nx)+1); Xk=x*(exp(-j*dw*nx*k)/K;subplot(4,2,2*r+1);stem(nx,x);axis(0,K*Nx-1,0,1.1);ylabel(x(n);subplot(4,2,2*r+2);plot(k*dw,abs(Xk);axis(-4,4,0,1.1*max(abs(Xk);ylabel(X(k);end程序運行結果如圖4-3

46、所示。由圖可見，序列的重復周期數越多，頻譜越是向幾個頻點集中。當序列的周期數趨于無窮大時，頻譜轉化為離散譜。圖4-32、離散傅里葉變換（DFT）（1）DFT與IDFT在實際中常常使用有限長序列。如果有限長序列為x(n)，則該序列的離散傅里葉變換對可表示為從離散傅里葉變換定義式可以看出，有限長序列在時域上是離散的，在頻域上也是離散的，式中，即僅在單位圓上N個等間距的點上取值，這為使用計算機進行處理帶來了方便。由有限長序列的傅里葉變換和逆變換定義可知，DFT和DFS的變換公式非常相似，因此，在程序編寫上也基本一致。例4-4 已知x(n)=0,1,2,3,4,5,6,7，求其DFT和IDFT。要求：

47、 畫出序列傅里葉變換對應的|X(k)|和argX(k)圖形。 畫出x(n)圖形，并與IDFTX(k)圖形進行比較。程序清單如下：xn=0,1,2,3,4,5,6,7;N=length(xn);n=0:N-1;k=0:N-1;Xk=xn*exp(-j*2*pi/N).(n*k);x=(Xk*exp(j*2*pi/N).(n*k)/N;subplot(2,2,1);stem(n,xn);title(x(n);axis(-1,N,1.1*min(xn),1.1*max(xn);subplot(2,2,2);stem(n,abs(x);title(IDFT|X(k)|);axis(-1,N,1.1*m

48、in(x),1.1*max(x);subplot(2,2,3),stem(k,abs(Xk);title(|X(k)|);axis(-1,N,1.1*min(abs(Xk),1.1*max(abs(Xk);subplot(2,2,4),stem(k,angle(Xk);title(arg|X(k)|);axis(-1,N,1.1*min(angle(Xk),1.1*max(angle(Xk);程序運行結果如圖4-4所示。由圖可見，與周期序列不同，有限長序列本身是僅有N點的離散序列，相當于周期序列的主值部分。因此，其頻譜也對應序列的主值部分，是長度為N的離散序列。圖4-4（2）DFT與DFS的聯

49、系將周期序列的傅里葉級數變換對和有限長序列的離散傅里葉變換對進行比較可見，兩者的區別僅僅是將周期序列換成了有限長序列x(n)，同時，由于式中的周期性，因而有限長序列的離散傅里葉變換實際上隱含著周期性。例4-5 已知周期序列的主值x(n)=0,1,2,3,4,5,6,7，求x(n)的周期重復次數為4次時的DFS。要求： 畫出主值序列周期序列的波形。 畫出周期序列傅里葉變換對應的和arg的圖形。程序清單如下：xn=0,1,2,3,4,5,6,7;N=length(xn);m=0:N-1;n=0:4*N-1;k=0:4*N-1;xn1=xn(mod(n,N)+1);Xk=xn1*exp(-j*2*p

50、i/N).(n*k);subplot(2,2,1);stem(m,xn);title(x(n);subplot(2,2,2);stem(n,xn1);title(周期序列);subplot(2,2,3),stem(k,abs(Xk);title(|X(k)|);axis(-1,4*N,1.1*min(abs(Xk),1.1*max(abs(Xk);subplot(2,2,4),stem(k,angle(Xk);title(arg|X(k)|);axis(-1,4*N,1.1*min(angle(Xk),1.1*max(angle(Xk);程序運行結果如圖4-5所示。與例4-4相比，有限長序列x

51、(n)可以看成是周期序列的一個周期；反之，周期序列可以看成是有限長序列以N為周期的周期延拓。頻域上的情況也一樣。從這個意義上說，周期序列只有有限個序列值有意義。（3）DFT與DTFT的聯系 若離散時間非周期序列為x(n)，則它的離散傅里葉變換（DTFT）對定義為其中稱為序列的頻譜。可表示為，稱為序列的幅圖4-5度譜，稱為序列的相位譜。由DTFT的定義可見，序列在時域是離散的、非周期的，在頻域是連續的、周期的。與有限長序列相比，僅在單位圓上取值，X(k)是在單位圓上N個等間距的點上取值。因此，連續譜可由離散譜X(k)經插值后得到。例4-6 求有限長序列x(n)=0,1,2,3,4,5,6,7的D

52、TFT，將（）區間分成500份。要求： 畫出序列x(n)的圖形。 畫出由DTFT求出的幅度譜和相位譜的圖形。程序清單如下：xn=0,1,2,3,4,5,6,7;N=length(xn);n=0:N-1;w=linspace(-2*pi,2*pi,500);X=xn*exp(-j*n*w);subplot(3,1,1);stem(n,xn);box on;title(x(n);subplot(3,1,2);plot(w,abs(X);title(|X(ejomega)|);axis(-2*pi,2*pi,1.1*min(abs(X),1.1*max(abs(X);subplot(3,1,3),p

53、lot(w,angle(X);title(ephi(omega);axis(-2*pi,2*pi,1.1*min(angle(X),1.1*max(angle(X);程序運行結果如圖4-6所示：圖4-6與圖4-4相比，兩者有一定差別。主要原因在于，該例進行DTFT時，在單位圓上取250個點進行分割；而圖4-4進行DFT時，X(k)是在單位圓上以N=8的等間距取值，X(k)的序列長度與相比不夠長。若將x(n)=0,1,2,3,4,5,6,7補零拓展至長度N=100，再求其DFT，如圖4-7所示。與例4-6相比，|X(k)|和argX(k)圖形接近和的圖形。注意圖4-7對應0,2區間。圖4-7程序

54、清單如下：N=100;xn=0,1,2,3,4,5,6,7,zeros(1,N-8);n=0:N-1;k=0:N-1;Xk=xn*exp(-j*2*pi/N).(n*k);x=(Xk*exp(j*2*pi/N).(n*k)/N;subplot(2,2,1);stem(n,xn);title(x(n);axis(-1,N,1.1*min(xn),1.1*max(xn);subplot(2,2,2);stem(n,abs(x);title(IDFT|X(k)|);axis(-1,N,1.1*min(x),1.1*max(x);subplot(2,2,3),stem(k,abs(Xk);title(

55、|X(k)|);axis(-1,N,1.1*min(abs(Xk),1.1*max(abs(Xk);subplot(2,2,4),stem(k,angle(Xk);title(arg|X(k)|);axis(-1,N,1.1*min(angle(Xk),1.1*max(angle(Xk);3、FFT（1）MATLAB提供的FFT函數 由理論學習可知，DFT是唯一在時域和頻域均離散的變換方法，它適用于有限長序列。盡管這種變換方法是可以用于數值計算的，但如果只是簡單地按照定義進行數據處理，當序列長度很大時，將占用很大的內存空間，且運算時間很長。 快速傅里葉變換是用于DFT運算的高效快速算法的統稱，

56、FFT只是其中的一種。FFT主要有時域抽取算法和頻域抽取算法，基本思想是將一個長度為N的序列分解成多個段序列，如基2算法、基4算法等，大大地縮短了DFT的時間。有關詳細理論可參考教材。MATLAB提供了進行FFT的函數fft和ifft分別用于計算DFT和IDFT。例4-7 已知一個長度為8的時域離散信號，n1=0，n2=7，在n0=4前為0，n0以后為1.對其進行FFT變換，作時域信號及DFT、IDFT的圖形。程序清單如下：n1=0;n2=7;n0=4;n=n1:n2;N=length(n);xn=(n-n0)=0;subplot(2,2,1);stem(n,xn);title(x(n);k=

57、0:N-1;Xk=fft(xn,N);subplot(2,1,2);stem(k,abs(Xk);title(Xk=DFT(xn);xn1=ifft(Xk,N);subplot(2,2,2);stem(n,xn1);title(x(n)=IDFT(Xk);程序運行結果如圖4-8所示：圖4-9（2）用FFT進行頻譜分析 對有限長序列進行譜分析一個序號從n1到n2的時域有限長序列x(n)，它的頻譜定義為它的離散傅里葉變換，且在Nyquist頻率范圍內有界并連續。序列的長度為N，則N=n2-n1+1。計算x(n)的DFT得到的是的N個樣本點。其中數字頻率為式中：為數字頻率的分辨率；k取對應-（N-1

58、）/2到（N-1）/2區間的整數。在實際使用中，往往要求計算出信號以模擬頻率為橫坐標的頻譜，此時對應的模擬頻率為式中：D為模擬頻率的分辨率或頻率間隔；Ts為采樣信號的周期，Ts=1/Fs；定義信號的長度L=NTs。在使用FFT進行DFT的高效運算時，一般不直接用n從n1到n2的x(n)，而是取的主值區間（n=0,1,N-1）的數據，經FFT將產生N個數據，定位在k=0,1,N-1的數字頻率點上，即對應0,2。如果要顯示-，范圍的頻譜，則可以使用fftshift(X)進行位移。例4-8 已知有限長序列x(n)=1，2，3，2，1，其采樣頻率Fs=10Hz。請使用FFT計算其頻譜。程序清單如下Fs

59、=10;xn=1,2,3,2,1;N=length(xn);D=2*pi*Fs/N;k=floor(-(N-1)/2:(N-1)/2);X=fftshift(fft(xn,N);subplot(1,2,1);plot(k*D,abs(X),o:);title(幅度頻譜);xlabel(rad/s);subplot(1,2,2);plot(k*D,angle(X),o:);title(相位頻譜);xlabel(rad/s);程序運行結果如圖4-10所示：圖4-10由圖4-10可知，當有限長序列的長度N=5時，頻譜的樣本點數也為5，頻率點之間的間距非常大，即分辨率很低。即使使用了plot命令的插值

60、功能，顯示出的曲線仍是斷續的，與真實曲線有較大誤差。改變分辨率的基本方法是給輸入序列補零，即增加頻譜的密度。這種方法只是改善了圖形的視在分辨率，并不增加頻譜的細節信息。將上述有限長序列x(n)1,2,3,2,1末尾補零到N=1000點，將程序改為：Fs=10;N=1000;xn=1,2,3,2,1;Nx=length(xn);xn=1,2,3,2,1,zeros(1,N-Nx-1);D=2*pi*Fs/N;k=floor(-(N-1)/2:(N-1)/2);X=fftshift(fft(xn,N);subplot(1,2,1);plot(k*D,abs(X);title(幅度頻譜);xlabe