哈爾濱工程大學(xué) 數(shù)字信號(hào)處理實(shí)驗(yàn)二

1、 實(shí) 驗(yàn) 報(bào) 告課程名稱實(shí)驗(yàn)項(xiàng)目名稱實(shí)驗(yàn)類型實(shí)驗(yàn)學(xué)時(shí)班級(jí)學(xué)號(hào)姓名指導(dǎo)教師實(shí)驗(yàn)室名稱實(shí)驗(yàn)時(shí)間實(shí)驗(yàn)成績預(yù)習(xí)部分實(shí)驗(yàn)過程表現(xiàn)實(shí)驗(yàn)報(bào)告部分總成績教師簽字日期實(shí)驗(yàn)二 離散時(shí)間傅立葉變換 實(shí)驗(yàn)報(bào)告1、 實(shí)驗(yàn)原理 1、經(jīng)由正、逆離散時(shí)間傅里葉變換表達(dá)的信號(hào)傅里葉表達(dá)式是信號(hào)分析的一個(gè)關(guān)鍵部分，下面方程分別是分析方程和綜合方程。X()=類似地，當(dāng)LTI系統(tǒng)用于濾波時(shí)，作為沖擊響應(yīng)離散時(shí)間傅里葉變換的頻率響應(yīng)，提供了LTI系統(tǒng)簡介的描述。離散時(shí)間傅里葉變換X()是的周期復(fù)值函數(shù)，周期總是2，并且基周期通常選在區(qū)間-，）上。對離散時(shí)間傅里葉變換DTFT來說有兩個(gè)問題：第一個(gè)問題是：DTFT的定義對無限長信號(hào)是有效

2、的第二個(gè)問題是：DTFT是連續(xù)變量的函數(shù)在MATLAB中，任何信號(hào)必須是有限長度的，這就使第一點(diǎn)成為問題。因此，不可能使用MATLAB計(jì)算無限長信號(hào)的DTFT。有一個(gè)值得注意的例外情形，當(dāng)能從變換定義式推導(dǎo)出解析式并只是計(jì)算它時(shí)，可以使用MATLAB計(jì)算無限長信號(hào)的DTFT。例如在xn=un，xn具有有理的DTFT的情形下。 2、第二個(gè)問題是頻率抽樣問題。MATLAB擅長在有限網(wǎng)格點(diǎn)上計(jì)算DTFT。通常選擇足夠多的頻率以使繪出的圖平滑，逼近真實(shí)的DTFT。對計(jì)算有利的最好選擇是在（-，）區(qū)間上一組均勻地隔開的頻率，或者對共軛對稱變換選擇0,區(qū)間。采用上述抽樣辦法，DTFT式變成X()= （3.

3、11）DTFT的周期性意味著在-<0區(qū)間上的數(shù)值是那些對k>N/2的數(shù)值。因?yàn)樯鲜绞窃谟邢迶?shù)量的頻率點(diǎn)=2k/N處計(jì)算，并在有限范圍內(nèi)求和，因此它是可計(jì)算的。由于信號(hào)長度必須是有限的（0n<L），這個(gè)求和式不適用于xn=un的情形。在對DTFT進(jìn)行抽樣時(shí)，并不要求N=L，盡管通常經(jīng)由DFT進(jìn)行計(jì)算（使用FFT算法）時(shí)，如果N=L很便利。的確，在使用N點(diǎn)DFT時(shí)，如果N>L，只須想像xn是增補(bǔ)零的即可；如果N<L，會(huì)較復(fù)雜。這種情形下，在正確應(yīng)用FFT計(jì)算N點(diǎn)DFT前，需要對xn進(jìn)行時(shí)間混疊，因此此時(shí)應(yīng)確保總是比遠(yuǎn)處是時(shí)間序列（即N>=L）點(diǎn)數(shù)多得多的頻率處來

4、計(jì)算DTFT。3、 計(jì)算DTFT需要兩個(gè)函數(shù)，MATLAB的freqz函數(shù)計(jì)算無限長信號(hào)，dtft（h，H）函數(shù)計(jì)算有限長信號(hào)的DTFT。二：實(shí)驗(yàn)內(nèi)容 實(shí)驗(yàn)內(nèi)容包括脈沖信號(hào)的DTFT，asinc的M文件，無限長信號(hào)的DTFT，指數(shù)信號(hào)，復(fù)指數(shù)信號(hào)這五個(gè)部分。1. 脈沖信號(hào)的DTFTa、已知矩形脈沖rn的定義式為在0n<L上rn=1，n為其余值時(shí)，rn均為0。將rn代入DTFT的分析方程X()=中，計(jì)算化簡后可得：。b、建立DTFT函數(shù)（.M文件）實(shí)驗(yàn)程序清單1：functionA,B=dtft(h,N)%DTFT calculate DTFT at N equally spaced fr

5、equencies%Usage:%A,B=dtft(h,N)%h:finite-length input vector,whose length is L%N:number of frequencies for evaluation over -pi,pi)%=>constraint:N>=L%H:DTFT values(complex)%W:(2nd output)vector of freqs where DTFT is computed%N=fix(N);L=length(h); h=h(:); %<-for vectors ONLY!if(N<L) error(

6、'DTFT: # data samples cannot exceed # freq samples')endB=(2*pi/N)*0:(N-1)'mid=ceil(N/2)+1;B(mid:N)=B(mid:N)-2*pi;%movepi,2pi)to-pi,0)B=fftshift(B);A=fftshift(fft(h,N);%move negative freq components實(shí)驗(yàn)程序清單2:（當(dāng)w=80時(shí)）h=ones(1,12);A,B=dtft(h,80); %w=80時(shí)subplot(311),plot(B,real(A); %把繪圖窗口分成三行一

7、列三個(gè)子圖，然后在第一塊區(qū)域作圖xlabel('B'),ylabel('A,real part'); %以B為橫軸，A的實(shí)部為縱軸subplot(312),plot(B,imag(A); %把繪圖窗口分成三行一列三個(gè)子圖，然后在第二塊區(qū)域作圖xlabel('B'),ylabel('A,imag part'); %以B為橫軸，A的虛部為縱軸subplot(313),plot(B,abs(A); %把繪圖窗口分成三行一列三個(gè)子圖，然后在第三塊區(qū)域作圖xlabel('B'),ylabel('A amplitude

8、'); %以B為橫軸，A的幅值為縱軸實(shí)驗(yàn)仿真圖：實(shí)驗(yàn)程序清單3:（當(dāng)w=100時(shí)）h=ones(1,12);A,B=dtft(h,100); %w=100時(shí)subplot(311),plot(B,real(A);xlabel('B'),ylabel('A,real part');subplot(312),plot(B,imag(A);xlabel('B'),ylabel('A,imag part');subplot(313),plot(B,abs(A);xlabel('B'),ylabel('A a

9、mplitude');實(shí)驗(yàn)仿真圖：實(shí)驗(yàn)程序清單4:（當(dāng)w=120時(shí)）h=ones(1,12);A,B=dtft(h,120); %w=120時(shí)subplot(311),plot(B,real(A);xlabel('B'),ylabel('A,real part');subplot(312),plot(B,imag(A);xlabel('B'),ylabel('A,imag part');subplot(313),plot(B,abs(A);xlabel('B'),ylabel('A amplitud

10、e');實(shí)驗(yàn)仿真圖：實(shí)驗(yàn)結(jié)果分析：由以上三組不同頻率（w=80,w=100,w=120）的仿真圖形可知，當(dāng)w越大，取的點(diǎn)數(shù)越多時(shí)，繪制的圖形越來越平滑。c、注意asinc函數(shù)零點(diǎn)的位置是規(guī)則分布的。對奇數(shù)長脈沖，比如L=15的脈沖重復(fù)進(jìn)行DTFT計(jì)算并繪出幅度；同樣再次檢驗(yàn)零點(diǎn)位置，注意峰值高度。實(shí)驗(yàn)程序清單：n=0:14;h=ones(1,15);X,W=dtft(h,150);subplot(211),plot(W,abs(X); %把繪圖窗口分成二行一列二個(gè)子圖，然后在第一塊區(qū)域作圖grid; %圖層為網(wǎng)格xlabel('W'),ylabel('|X(w)

11、|'); %X軸為w，Y軸為X(w)幅值H,X=dtft(X,150);subplot(212),plot(X,abs(H); %把繪圖窗口分成二行一列二個(gè)子圖，然后在第二塊區(qū)域作圖grid;xlabel('X'),ylabel('|H(w)|');%X軸為X，Y軸為H(w)幅值實(shí)驗(yàn)仿真圖： 實(shí)驗(yàn)結(jié)果分析：由仿真結(jié)果可以看出函數(shù)零點(diǎn)關(guān)于x=0對稱。在對稱軸處幅值最大，離對稱軸越遠(yuǎn)，幅值逐漸減小。d、對asinc函數(shù)零點(diǎn)的間距與asinc函數(shù)的直流值，確定出通用規(guī)則。 零點(diǎn)間距為/6，由數(shù)學(xué)推導(dǎo)：在a部分中，當(dāng)L=12時(shí)，由R()=0可得sin(wL/2)

12、=0，即wL/2=k* 則w=k*/6。故得出以上結(jié)論。在實(shí)驗(yàn)a部分的仿真圖中，可觀察到幅值高度為12。則零點(diǎn)間距*幅值=(/6)*12=2*。2、 asinc的M文件編寫一個(gè)MATLAB的函數(shù)如asinc(w,L),直接從混疊函數(shù)表達(dá)式中計(jì)算在頻率格上的asinc(w,L)。該函數(shù)應(yīng)該有兩個(gè)輸入：唱的L和頻率w的向量。函數(shù)必須檢查被零除的情況，如w=0時(shí)。直接計(jì)算混疊sinc函數(shù)得到脈沖信號(hào)的DTFT。繪出幅度，保存該圖以便將其與用dtft得到的結(jié)果進(jìn)行比較。function H,W=asinc(w,L) %定義asinc函數(shù)N=fix(L);if(N<0) error('Pl

13、ease input a positive integer.(L>0)');end%if(sin(w/2)=0)%endW=(2*pi/N)*0:(N-1)'H=sin(w*N/2)/sin(w/2); 3、 指數(shù)信號(hào)對于信號(hào)xn=(0.9)n*un，使用freqz函數(shù)計(jì)算其DTFT X()。a、對w在區(qū)間-，）上繪出幅值與相位特性。這需要從freqz返回的X,W向量的移位。解釋為什么幅度特性是w的偶函數(shù)，而相位特性是w的奇函數(shù)。實(shí)驗(yàn)程序清單1：N=100;a=1,-0.9;b=1;H,W=freqz(b,a,N);W=-pi:0.1:pi;H=freqz(b,a,W);

14、subplot(211);%把繪圖窗口分成二行一列二個(gè)子圖，然后在第一塊區(qū)域作圖plot(W,abs(H);grid,title('MAGNITUDE RESPONSE')xlabel('FREQUENCY W'),ylabel('|H(w)|')subplot(212);%把繪圖窗口分成二行一列二個(gè)子圖，然后在第二塊區(qū)域作圖plot(W,angle(H);grid,title('PHASE RESPONSE')xlabel('FREQUENCY W'),ylabel('DEGREES')實(shí)驗(yàn)仿真圖

15、：實(shí)驗(yàn)結(jié)果分析：將該指數(shù)信號(hào)代入分析方程X()=中，計(jì)算化簡后可從結(jié)果X()看出：該結(jié)果的幅值關(guān)于w是偶函數(shù)，相位關(guān)于w是奇函數(shù)。b、推算一階系統(tǒng)的幅度特性與相位特性的表達(dá)式。設(shè)一階差分方程：yn=xn+ayn-1,其系統(tǒng)方程H()=其幅度相位c、直接以這些表達(dá)式來計(jì)算幅度特性與相位特性，并與用freqz函數(shù)計(jì)算出的結(jié)果相對比。 xn=un的幅度特性為 相位特性為實(shí)驗(yàn)程序清單：N=1000;a=0.9;w=(2*pi/N)*-(N-1)/2:(N-1)/2; m=1./(1-a*cos(w)+(a*sin(w).2);n=1./(1-a*cos(w).2+(a*sin(w).2);H=sqrt

16、(m.2+n.2);%H的表達(dá)式subplot(2,1,1);%把繪圖窗口分成二行一列二個(gè)子圖，然后在第一塊區(qū)域作圖plot(w,sqrt(H);grid;xlabel(' W'),ylabel('|H(w)|');phase=(a*sin(w)./(1-a*cos(w);subplot(2,1,2);plot(w,atan(-phase);%繪制連續(xù)圖形grid,title('PHASE ');xlabel(' W'),ylabel('DEGREES');實(shí)驗(yàn)仿真圖：實(shí)驗(yàn)結(jié)果分析：以上兩個(gè)圖對比，可以發(fā)現(xiàn)用表達(dá)式

17、來計(jì)算幅度特性和相位特性這種方式與用freqz函數(shù)計(jì)算這種方式，它們最終得出的結(jié)果是一樣的。3、復(fù)指數(shù)信號(hào)：a、取z=0.95ej3pi/11,對0,30，繪出xn=zn*un。用subplot指令將實(shí)部和虛部對n繪出雙子圖。實(shí)驗(yàn)程序清單：z0=0.95*exp(j*3*pi/11);n=0:30;x=z0.n;subplot(311);%用subplot指令繪出子圖stem(n,x); %以x為橫軸，n為縱軸繪制離散信號(hào)波形title('復(fù)指數(shù)信號(hào)');xlabel('n'),ylabel('x');subplot(3,1,2); stem(n

18、,real(x);title('REAL PART');xlabel('n'),ylabel('real');subplot(313);stem(n,imag(x); %以n為橫軸，x的虛部為縱軸繪制離散信號(hào)波形title('IMAG PART');xlabel('n'),ylabel('imag');實(shí)驗(yàn)仿真圖：b、再取一次z=0.95ej3pi/11，計(jì)算其DTFT，并對w繪出幅度。注意作為w的函數(shù)，幅度響應(yīng)的尖峰位于何處。把尖峰位置與z的極坐標(biāo)形式聯(lián)系起來。實(shí)驗(yàn)程序清單：b=1;a=1,-0.

19、95*exp(j*3*pi/11);w=-(8/11)*pi:0.05*pi:(14/11)*pi;H=freqz(b,a,w);stem(w,abs(H); %以w為橫軸，H的幅度為縱軸繪制離散圖形grid;xlabel('w'),ylabel('|H(w)|');實(shí)驗(yàn)仿真圖：c、如果z的角度變?yōu)?/5，粗略繪出預(yù)期的DTFT，并對w繪出幅度，并通過對freqz的計(jì)算進(jìn)行繪圖來驗(yàn)證。實(shí)驗(yàn)程序清單：w=-0.4*pi:0.05*pi:1.6*pi; %W的表達(dá)式z0=0.95*exp(j*3*pi/5);%Z0的表達(dá)式H=1./(1-z0*exp(-j.*w);

20、%H的表達(dá)式subplot(211);stem(w,abs(H);grid;xlabel('w'),ylabel('|H(w)|');b=1;a=1,-0.95*exp(j*3*pi/5);c=-0.4*pi:0.05*pi:1.6*pi;H=freqz(b,a,c); %調(diào)用freqz函數(shù)subplot(212); %繪制子圖stem(c,abs(H); %繪制離散圖形grid;xlabel('c'),ylabel('|H(c)|');實(shí)驗(yàn)仿真圖：d、試建立一個(gè)聯(lián)系帶寬與r的簡單表達(dá)式。實(shí)驗(yàn)程序清單：a1=1,-0.975*exp(j*3*pi/5);a2=1,-0.95*exp(j*3*pi/5);a3=1,-0.9*exp(j*3*pi/5);a4=1,-0.8*exp(j*3*pi/