線性時不變離散時間系統(tǒng)的頻域分析

1、云南大學(xué)信息學(xué)院數(shù)字信號處理實驗報告實驗四 線性時不變離散時間系統(tǒng)的頻域分析實驗室名稱：格物樓 1204 實驗時間：2015年11月6日姓名： 成績教師簽名： 年 月 日一、實驗?zāi)康模?）學(xué)習(xí)MATLAB軟件及其在信號處理中的應(yīng)用，加深對常用連續(xù)時間信號的理解。（2）連續(xù)時間信號在時域和頻域中的抽樣效果。 （3）熟悉MATLAB模擬濾波器的。二、實驗內(nèi)容Q4.1 修改程序P3.1中，取三個不同的M值，當(dāng)時計算并畫出式(2.13)所示滑動平均濾波器的幅度和相位譜。證明由幅度和相位譜表現(xiàn)出的對稱類型。它表示了那種類型的濾波器？你現(xiàn)在能解釋習(xí)題Q2.1的結(jié)果嗎？修改后的程序P3.1也可以用于計算并畫

2、出線性時不變離散時間系統(tǒng)的頻率響應(yīng)，該系統(tǒng)用形如式(4.12)的傳輸函數(shù)描述。Q4.2 使用修改后的程序P3.1，計算并畫出當(dāng)時傳輸函數(shù) 的因果線性時不變離散時間系統(tǒng)的頻率響應(yīng)。它表示那種類型的濾波器？Q4.3 對下面的傳輸函數(shù)重做習(xí)題Q4.2 4.2和4.3給出的兩個濾波器之間的區(qū)別是什么？你將選擇哪一個濾波器來濾波，為什么？ Q4.6 使用zplane分別生成4.2和4.2所確定的兩個濾波器的零極點圖。討論你的結(jié)果。Q4.7 用程序P4.1計算并畫出近似理想低通濾波器的沖激響應(yīng)。低通有限沖激響應(yīng)濾波器的長度是多少？在程序P4.1中，那個語句確定濾波器的長度？那個參數(shù)控制截止頻率？Q4.8

3、修改程序P4.1 ，計算并畫出式(4.39)所示的長度為20，截止角頻率為的有限沖激響應(yīng)低通濾波器的沖激響應(yīng)。Q4.9 修改程序P4.1 ，計算并畫出式(4.39)所示的長度為15，截止角頻率為的有限沖激響應(yīng)低通濾波器的沖激響應(yīng)。Q4.10 編寫一個MATLAB程序，計算并畫出式(4.39)所示有限沖激響應(yīng)低通濾波器的振幅響應(yīng)。使用這個程序，選取幾個不同的N值，畫出振幅響應(yīng)并討論你的結(jié)果。Q4.11運行程序P4.2，計算并畫出一個長度為2的滑動平均濾波器的增益響應(yīng)。從圖中驗證3dB截止頻率在/2處。Q4.23 用MATLAB產(chǎn)生如下兩個因果系統(tǒng)傳輸函數(shù)的零極點圖 研究生成的零極點圖，你可以推斷

4、它們的穩(wěn)定性嗎？三、實驗器材及軟件1. 微型計算機1臺2. MATLAB 7.0軟件四、實驗原理熟悉MATLAB中產(chǎn)生信號和繪制信號的基本命令，MATLAB已被開發(fā)成能對數(shù)據(jù)向量或矩陣運算的工具。序列以向量的形式儲存，并且所有的信號被限定為因果的和有限長的，采用何種步驟執(zhí)行程序，要根據(jù)MATLAB所運行的平臺來決定。(1) 若hn表示一個線性時不變離散時間系統(tǒng)的沖激響應(yīng)，對hn做離散時間傅里葉變換得到其頻率響應(yīng)，即。(2) 通常，是一個周期為2的的復(fù)值函數(shù)，可以根據(jù)實部，虛部或者幅度相位來表示該函數(shù)。因此， 其中和分別是的實部和虛部，并且又| 稱為幅度響應(yīng)，而稱為線性時不變離散時間系統(tǒng)的相位響

5、應(yīng)。(3) 線性時不變系統(tǒng)的增益函數(shù)g()定義為 dB，增益函數(shù)的相反數(shù)a()=-g(),稱為衰減或損益函數(shù)。(4) 對于用實沖激響應(yīng)hn描述的離散時間系統(tǒng)，幅度函數(shù)是的偶函數(shù)，即| =| ；而相位函數(shù)是的奇函數(shù)，即。同樣，是的偶函數(shù)，是的奇函數(shù)。(5) 線性時不變離散時間系統(tǒng)的頻率響應(yīng)可以由輸出序列yn的傅里葉變換與輸入序列xn的傅里葉變換相比得到，即/。(6) 線性時不變離散時間系統(tǒng)的沖激響應(yīng)hn的z變換H(z)，稱為傳輸函數(shù)或系統(tǒng)函數(shù)。H(z)可由輸出序列yn的z變換Y(z)與輸入序列xn的z變換X(z)相比得到，即H(z)=Y(z)/X(z)。(7) 穩(wěn)定因果系統(tǒng)函數(shù)的傳輸函數(shù)H(z)

6、的所有極點都必須嚴(yán)格在單位圓內(nèi)。(8) 線性常系數(shù)差分方程描述的線性時不變系統(tǒng)，傳輸函數(shù)H(z)可以表示為 。五、實驗步驟(1)打開MATLAB (2)新建M文件 (3)編寫代碼 (4)運行代碼 (5)得到并分析結(jié)果六、實驗記錄（數(shù)據(jù)、圖表、波形、程序等）Q4.1n = 0:100;s1 = cos(2*pi*0.05*n); % A low-frequency sinusoids2 = cos(2*pi*0.47*n); % A high frequency sinusoidx = s1+s2;% Implementation of the moving average filterM =

7、input('Desired length of the filter = ');num = ones(1,M);den = filter(num,1,x)/M;clf;% Compute the frequency samples of the DTFTw = 0:2*pi;h = freqz(num, den, w);% Plot the DTFTsubplot(2,2,1)plot(w/pi,real(h);gridtitle('Real part of H(ejomega)')xlabel('omega /pi');ylabel('

8、;Amplitude');subplot(2,2,2)plot(w/pi,imag(h);gridtitle('Imaginary part of H(ejomega)')xlabel('omega /pi');ylabel('Amplitude');subplot(2,2,3)plot(w/pi,abs(h);gridtitle('Magnitude Spectrum |H(ejomega)|')xlabel('omega /pi');ylabel('Amplitude');subplot

9、(2,2,4)plot(w/pi,angle(h);gridtitle('Phase Spectrum argH(ejomega)')xlabel('omega /pi');ylabel('Phase in radians');M=3M=5M=10由圖可看出為低通濾波器。Q4.2w = 0:pi/511:pi;num = 0.15 0 -0.15;den=1 -0.5 0.7h = freqz(num, den, w);subplot(2,1,1)plot(w/pi,abs(h);gridtitle('Magnitude Spectrum

10、 |H(ejomega)|')xlabel('omega /pi');ylabel('Amplitude');subplot(2,1,2)plot(w/pi,angle(h);gridtitle('Phase Spectrum argH(ejomega)')xlabel('omega /pi');ylabel('Phase in radians');Q4.3修改4.2程序num = 0.15 0 -0.15;den=0.7 -0.5 1Q4.2和Q4.3的兩個濾波器，幅度譜是一樣的，相位譜Q4.3中的出現(xiàn)躍

11、變，我會選擇Q4.3 的濾波器。Q4.6式4.36的零極點圖。w = 0:pi/511:pi;num = 0.15 0 -0.15;den=1 -0.5 0.7h = zplane(num, den);式4.37的零極點圖。w = 0:pi/511:pi;num1 = 0.15 0 -0.15;den1=0.7 -0.5 1h1 = zplane(num1, den1);Q4.7clf; fc=0.25; n=-6.5:1:6.5; y=2*fc*sinc(2*fc*n);k=n+6.5; stem(k,y);title('N=13');

12、axis(0 13 -0.2 0.6); xlabel('時間序號 n');ylabel('振幅');grid n = Columns 1 through 13 -6.5000 -5.5000 -4.5000 -3.5000 -2.5000 -1.5000 -0.5000 0.5000 1.5000 2.5000 3.5000 4.5000 5.5000 Column 14 6.5000k = 0 1 2 3 4 5 6 7 8 9 10 11 12 13低通濾波器的長度為13，n=-6.5:1:6.5

13、決定了濾波器的長度。fc=0.25;控制截止頻率。Q4.8%程序P4.1 %截短的理想低通濾波器 clf; wc=0.45; fc=wc/2*pi; n=-9.5:1:9.5; y=2*fc*sinc(2*fc*n);k=n+9.5; stem(k,y);title('N=20');axis(0 20 -0.2 0.6); xlabel('時間序號n');ylabel('振幅');grid; Q4.10% Progra

14、m Q4_10  clear;  N = input('Enter the filter time shift N: '); No2 = N/2; fc = 0.25; n = -No2:1:No2; y = 2*fc*sinc(2*fc*n); w = 0:pi/511:pi; h 

15、;= freqz(y, 1, w); plot(w/pi,abs(h); grid; title(strcat('|H(ejomega)|, N=',num2str(N); xlabel('omega /pi'); ylabel('Amplitude');低通濾波器的幅度相應(yīng)（若干個n值）:從圖像可以得到觀察 隨著濾波器長度的增加，從通過到不通過變得更加陡峭，我們也可以看到吉布斯現(xiàn)象：當(dāng)濾波器增加時，幅度相應(yīng)更加趨向一個理想的低通特征。然

16、而隨著w增長，峰值是增加而不是降低。Q4.11function g,w=gain(num,den)  -gain函數(shù) w=0:pi/255:pi; h=freqz(num,den,w); g=20*log10(abs(h); M=2; -滑動平均低通濾波器的增益響應(yīng)程序 num=ones(1,M)/M; g,w=gain(num,1); plot(w/pi,g);grid; axis(0 1 -50 0.5) xlabel('omega/pi');ylabel('單位為db的增益'); title('M= ',num2str(M) 從圖中可以看出，在w=pi/2處增益對應(yīng)著3dB。Q4.23b=1; a=1,-1.848 0.85; zplane(b,a); title('H1(z)') b=1; a=1,-1.851 0.85; zplane(b,a); title('H2(