1、實驗三、 離散系統的Z域分析（一）實驗要求1） 學習和掌握離散系統的頻率特性及其幅度特性、相位特性的物理意義；2） 深入理解離散系統頻率特性的對稱性和周期性；3） 認識離散系統頻率特性與系統參數之間的關系；4） 通過閱讀、修改并調試本實驗系統所給源程序，加強計算機編程能力；（二）實驗內容1、計算差分方程 （1）用MATLAB計算差分方程 當輸入序列為 時的輸出結果 。 MATLAB程序如下： N=41;a=0.8 -0.44 0.36 0.22;b=1 0.7 -0.45 -0.6;x=1 zeros(1,N-1);k=0:1:N-1;h=filter(a,b,x);stem(k,h)xlab

2、el('n');ylabel('h(n)')    請給出了該差分方程的前41個樣點的輸出，即該系統的單位脈沖響應。(說明：y=filter(a,b,x)，計算系統對輸入信號向量x的零狀態響應輸出信號向量y,x與y長度相等，其中a和b是所給差分方程的相量。詳見教材P25-27)2、用MATLAB計算差分方程 所對應的系統函數的FT。 差分方程所對應的系統函數為： 其FT為用MATLAB計算的程序如下：k=256;num=0.8 -0.44 0.36 0.02;den=1 0.7 -0.45 -0.6;w=0:pi/k:pi

3、;h=freqz(num,den,w);subplot(2,2,1);plot(w/pi,real(h);gridtitle('實部')xlabel('omega/pi');ylabel('幅度')subplot(2,2,2);plot(w/pi,imag(h);gridtitle('虛部')xlabel('omega/pi');ylabel('Amplitude')subplot(2,2,3);plot(w/pi,abs(h);gridtitle('幅度譜')xlabel('

4、;omega/pi');ylabel('幅值')subplot(2,2,4);plot(w/pi,angle(h);gridtitle('相位譜') xlabel('omega/pi');ylabel('弧度') (說明：freqz為計算數字濾波器H(z)的頻率響應函數。h=freqz(num,den,w)為計算由向量w指定的數字頻率點上數字濾波器H(z)的頻率響應，結果存于h向量中。Num和den為H(z)分子和分母多項式向量。詳見教材P65)*練習:、P43中的例3、例5。、用MATLAB編程，畫出的頻率特性圖。3、

5、求解的Z反變換。參考程序：b=1;a=poly(-0.2 0.5 0.5 -0.6 -0.6);r,p,k=residuez(b,a)(說明：例程序：b = -4 8;a = 1 6 8;r,p,k = residuez(b,a)運行結果：r = -12 8p = -4 -2k = 則表示：那么：)（三）實驗報告要求1、簡述實驗目的和實驗原理，用幾何確定法分析實驗中選定的系統的頻率特性，并與計算機計算結果相對照，根據實驗結果，對系統頻率特性進行討論和總結。2、用MATLAB編程，畫出P66中的例的頻率特性圖。3、根據MATLAB求解的結果，寫出序列的表達式。實驗四、離散傅里葉變換及其快速算法

6、（一）實驗要求1） 通過離散傅立葉變換（即DFT）的報表表示進一步了解其計算方法及意義；2） 掌握實數序列的DFT系數的對稱特點；3） 學習利用DFT計算程序計算IDFT的方法；4)學習時間抽選奇偶分解FFT算法；5)深入掌握時間抽選奇偶分解FFT程序的編制方法；(二)實驗內容1、對連續的單一頻率周期信號 按采樣頻率 采樣，截取長度N分別選N =20和N =16，觀察其DFT結果的幅度譜。 此時離散序列 ，即k=8。用MATLAB計算并作圖，函數fft用于計算離散傅里葉變換DFT，程序如下： k=8; n1=0:1:19; xa1=sin(2*pi*n1/k); subplot(2,2,1)

7、plot(n1,xa1) xlabel('t/T');ylabel('x(n)'); xk1=fft(xa1);xk1=abs(xk1); subplot(2,2,2) stem(n1,xk1) xlabel('k');ylabel('X(k)'); n2=0:1:15; xa2=sin(2*pi*n2/k); subplot(2,2,3) plot(n2,xa2) xlabel('t/T');ylabel('x(n)'); xk2=fft(xa2);xk2=abs(xk2); subplot(2,

8、2,4) stem(n2,xk2) xlabel('k');ylabel('X(k)');         計算結果示于圖2.1，(a)和(b)分別是N=20時的截取信號和DFT結果，由于截取了兩個半周期，頻譜出現泄漏；(c) 和(d) 分別是N=16時的截取信號和DFT結果，由于截取了兩個整周期，得到單一譜線的頻譜。上述頻譜的誤差主要是由于時域中對信號的非整周期截斷產生的頻譜泄漏。 2、 對如下各序列進行譜分析，繪制出其幅頻特性曲線。（1）（2）對取64點FFT。由于為周期序列，周

9、期為16，所以取周期為16。參考程序如下：%用FFT對序列進行譜分析%；%X1:存放的向量，Y1: 存放的向量%X2:存放的向量，Y2: 存放的向量x1=1,1,1,1,0,0,0,0;n=0:15;x2=cos(pi*n/4)+cos(pi*n/8);i=0:7;subplot(3,2,1);stem(i,x1,'.');axis(0 7 0 1);%規定x軸和y軸的標值范圍xlabel('n');ylabel('x1(n)');y1=fft(x1,8);% 的8點FFTsubplot(3,2,3);stem(i,abs(y1),'.&

10、#39;);xlabel('(N=8 wk=2pik/N)k');ylabel('X1(k)');y1=fft(x1,64);% 的64點FFTi=0:63;subplot(3,2,5);stem(i,abs(y1),'.');axis(0 63 0 4);%規定x軸和y軸的標值范圍xlabel('(N=8 wk=2pik/N)k');ylabel('X1(k)');y2=fft(x2); %的16點FFTfigure;%另一幅圖subplot(2,2,1);stem(n,x2);title('x2(n)的

11、時域序列'); xlabel('n');ylabel('x2(n)');subplot(2,2,3);stem(n,abs(y2);title('x2(n)的幅頻特性'); xlabel('(N=16 wk=2pik/N)k');ylabel('X2(k)');（三）實驗報告要求1、設。分別計算在頻率區間 0,2上的16點和32點等間隔采樣，并繪制采樣的幅頻特性圖和相頻特性圖。如：%DFT的MATLAB計算xn=1 1 1 1; %輸入時域序列向量Xk16=fft(xn,16); %計算xn16點DFTXk

12、32=fft(xn,32); %計算xn32點DFT繪圖程序略。2、編寫序列的DFT運算程序。2解：參考程序如下：N=16;n=0:1:N-1;%時域采樣xn=cos(n*pi/4)+ cos(n*pi/8);k=0:1:N-1; %頻域采樣WN=exp(-j*2*pi/N);nk=n'*k;WNnk=WN.nk; %點乘方Xk=xn*WNnk;Subplot(2,1,1)Stem(n,xn) title('xn的時域序列'); xlabel('n');ylabel('x (n)');Subplot(2,1,2)Stem(k,abs(Xk

13、);title('xn的幅頻特性'); xlabel('(N=16 wk=2pik/N)k');ylabel('X (k)');（Fft與Ifft的說明：y=fft(x)是利用fft函數求解x的離散傅里葉變換；y=fft(x，N)，N表示離散傅里葉變換x的數據長度；函數Ifft的參數與函數Fft完全相同。例：fft在信號分析中的應用。使用頻譜分析方法從受噪聲污染的信號x(t)中鑒別出有用信號。如程序：t=0:0.001:1;%采樣周期為0.001s，即采樣頻率為1000Hz;%產生受噪聲污染的正弦波信號；x=sin(2*pi*100*t)+ sin(2*pi*200*t)+rand(size(t);subplot(2,1,1)plot(x(1:50);%畫出時域內的信號