實驗一時域離散信號與系統(tǒng)變換域分析(2015)(2)

1、實驗一 時域離散信號與系統(tǒng)變換域分析一、實驗目的1了解時域離散信號的產(chǎn)生及基本運算實現(xiàn)。2掌握離散時間傅里葉變換實現(xiàn)及系統(tǒng)分析方法。3. 熟悉離散時間傅里葉變換性質(zhì)。 4. 掌握系統(tǒng)Z域分析方法。5. 培養(yǎng)學生運用軟件分析、處理數(shù)字信號的能力。二、實驗設(shè)備 1、計算機2、Matlab7.0以上版本三、實驗內(nèi)容 1、對于給定的時域離散信號會進行頻譜分析，即序列的傅里葉變換及其性質(zhì)分析。2、對于離散系統(tǒng)會進行頻域分析及Z域分析。包括頻譜特性、零極點畫圖、穩(wěn)定性分析。3、對于差分方程會用程序求解，包括求單位沖擊序列響應(yīng)，零輸入響應(yīng)、零狀態(tài)響應(yīng)、全響應(yīng)，求其系統(tǒng)函數(shù)，及其分析。4、信號時域采樣及其頻譜

2、分析，序列恢復。5、擴展部分主要是關(guān)于語音信號的讀取及其播放。四、實驗原理1、序列的產(chǎn)生及運算。 在Matlab中自帶了cos、sin、exp（指數(shù)）等函數(shù)，利用這些函數(shù)可以產(chǎn)生實驗所需序列。序列的運算包括序列的加法、乘法，序列的移位，翻褶等。序列的加法或乘法指同序號的序列值逐項對應(yīng)相加或相乘，但Matlab中“+”“.*”運算是對序列的值直接進行加或乘，不考慮兩序列的序號是否相同，因此編程時考慮其序號的對應(yīng)。2、序列的傅里葉變換及其性質(zhì)。序列的傅里葉變換定義：，其幅度特性為，在Matlab中采用abs函數(shù)；相位特性為，在Matlab中采用angle函數(shù)。序列傅里葉變換的性質(zhì)：（1） FT的周

3、期性，實序列傅里葉變換的對稱性。對實序列和復序列分別進行傅里葉變換，通過圖形結(jié)果觀察周期性即對稱性。（2） FT的頻移特性，對序列在時域乘以，然后進傅里葉變換，比較其結(jié)果和直接對序列進行傅里葉變換的不同。（3） 時域卷積定理：若，對序列和進行線性卷積得到，分別對它們進行傅里葉變換，應(yīng)滿足。3、離散時間系統(tǒng)的Z域分析。已知離散時間系統(tǒng)的差分方程為，對等號兩邊進行Z變換，得到其系統(tǒng)函數(shù)及系統(tǒng)零極點，對系統(tǒng)函數(shù)進行反變換得到單位取樣響應(yīng)，根據(jù)單位取樣響應(yīng)或系統(tǒng)函數(shù)的系數(shù)可以得到頻率響應(yīng)，根據(jù)極點位置判斷系統(tǒng)穩(wěn)定性。4、 信號時域采樣及恢復。給定連續(xù)信號，對其用不同的采樣頻率進行采樣，根據(jù)時域采樣定理

4、，采樣信號的頻譜是原模擬信號頻譜沿頻率軸以為周期延拓而成的，并且要不失真地還原出模擬信號時，要滿足，因此當采樣頻率滿足和不滿足采樣定理時，所得到的頻譜是不同的。根據(jù)采樣信號進行信號恢復時，采用內(nèi)插公式實現(xiàn)。五、實驗步驟1、序列的基本運算1.1、產(chǎn)生余弦信號及帶噪信號 0<=n<=50（噪聲采用randn函數(shù)）1.2、已知 ， ,求兩個序列的和、乘積、序列x1的移位序列（右移2位）,序列x2的翻褶序列，畫出原序列及運算結(jié)果圖。2、序列的傅里葉變換2.1、已知序列。試求它的傅里葉變換，并且畫出其幅度、相角、實部和虛部的波形，并分析其含有的頻率分量主要位于高頻區(qū)還是低頻區(qū)。2.2、令，求

5、其傅立葉變換。分別用和對其進行采樣，求出離散時間傅立葉變換，畫出相應(yīng)頻譜，分析結(jié)果的不同及原因。3、序列的傅里葉變換性質(zhì)分析3.1、已知序列，求其傅里葉變換，并討論其傅里葉變換的周期性和對稱性。3.2、已知序列，求其傅里葉變換，并討論其傅里葉變換的周期性和對稱性。為了方便，考慮在兩個周期，例如中2M+1個均勻頻率點上計算FT，并且觀察其周期性和對稱性。為此給出function文件如下，求解FT變換：functionX,w=ft1(x,n,k)w=(pi/abs(max(k)/2)*kX=x*(exp(-j*pi/abs(max(k)/2).(n'*k)3.3、編寫程序驗證序列傅里葉變換

6、頻移性質(zhì)，時域卷積定理（時域卷積后的頻域特性）。（所需信號自行選擇）4、時域差分方程的求解4.1求解差分方程y(n)a1y(n-1)a2y(n-2)=b0x(n)b1x(n-1)的零狀態(tài)響應(yīng)和全響應(yīng)。已知X(n)為單位取樣序列，y(-1)=1,y(-2)=2,a1=0.5,a2=0.06,b0=2,b13。5、離散系統(tǒng)的Z域分析5.1、利用系統(tǒng)函數(shù)分析系統(tǒng)的穩(wěn)定性。假設(shè)系統(tǒng)函數(shù)如下式：，試判斷系統(tǒng)是否穩(wěn)定。5.2、已知線性時不變系統(tǒng)的系統(tǒng)函數(shù)，編寫程序求其單位取樣響應(yīng)，頻率響應(yīng)及系統(tǒng)零極點，并畫出相應(yīng)圖形。6、創(chuàng)新訓練拓展內(nèi)容6.1、利用Matlab自帶的錄音功能，或利用Goldwave等音頻

7、編輯軟件，對語音或其他音頻信號進行采集并保存為*.wav文件。要求：（1）采用不同的采樣頻率（2000Hz，4000Hz，8000Hz，16000Hz等）。（2）對采集得到的信號進行播放，并畫圖。（3）分析在不同采樣頻率下得到的信號有何不同。6.2、設(shè)定一個連續(xù)時間信號，進行抽樣和恢復，要求分析不同采樣頻率對恢復結(jié)果的影響，給出實驗程序及各關(guān)鍵步驟圖形結(jié)果。6.3、設(shè)計內(nèi)容：設(shè)計一個離散系統(tǒng)，給定系統(tǒng)函數(shù)或差分方程，設(shè)定激勵及初始條件。要求：（1）繪制系統(tǒng)函數(shù)零極點圖，判斷穩(wěn)定性；（2）求單位脈沖響應(yīng)h（n）；（3）求系統(tǒng)零輸入響應(yīng)及零狀態(tài)響應(yīng)，要求零狀態(tài)響應(yīng)采樣三種方法求解（卷積的方法、迭代

8、解法、變換域求解方法），激勵自定；（4）分析系統(tǒng)頻響特性，畫出頻響函數(shù)幅頻曲線和相頻曲線。六、實驗要求第一部分：驗證實驗內(nèi)容根據(jù)給定的實驗內(nèi)容，部分實驗給出了參考程序段，見下面各段程序。請基于Matlab環(huán)境進行驗證實驗。第二部分：編程實驗內(nèi)容對于給定的實驗內(nèi)容中，沒有參考程序段的部分，進行編程，給出實驗結(jié)果，并進行相應(yīng)的分析。第三部分：創(chuàng)新訓練拓展內(nèi)容此部分內(nèi)容，要求給出程序設(shè)計流程圖（畫法見附錄3），給出程序內(nèi)容的解釋，并對結(jié)果進行分析。7、 思考題 下面四個二階網(wǎng)絡(luò)的系統(tǒng)函數(shù)具有一樣的極點分布：1） 2）3） 4）請分析研究零點分布對于單位脈沖響應(yīng)的影響。 要求： （1） 分別畫出各系統(tǒng)

9、的零、 極點分布圖；（2） 分別求出各系統(tǒng)的單位脈沖響應(yīng)，并畫出其波形；（3） 分析零點分布對于單位脈沖響應(yīng)的影響。八、實驗參考資料1、高西全,丁玉美.數(shù)字信號處理M.西安：西安電子科技大學出版社,20082、張德豐.詳解MATLAB 數(shù)字信號處理M.北京：電子工業(yè)出版社,20103、王月明,張寶華.MATLAB基礎(chǔ)與應(yīng)用教程M.北京:北京大學出版社,2012附錄1：實驗報告要求：實驗名稱：-班級： 組號： 姓名1（學號）、姓名2（學號）、姓名3（學號）一、實驗目的二、實驗主要內(nèi)容(要根據(jù)自己組所做內(nèi)容寫，做了的寫，沒做的不要寫)例如：1. 對序列的產(chǎn)生和運算方法進行實現(xiàn)2. 序列的傅里葉變換

10、實現(xiàn)、性質(zhì)及分析等等三、實驗主要儀器、設(shè)備及軟件四、實驗步驟、結(jié)果與分析例如：1. 序列的運算序列為，進行加法、乘法、運算 運算結(jié)果為2. 序列的傅里葉變換實現(xiàn)及分析（1） 已知序列。試求它的傅里葉變換，并且畫出其幅度、相角、實部和虛部的波形，并分析其含有的頻率分量主要位于高頻區(qū)還是低頻區(qū)。 程序 結(jié)果 分析（2）序列的傅里葉變換性質(zhì)分析 注1：（包括程序框圖及代碼、圖形、數(shù)據(jù)等），其中程序框圖、代碼、圖形可以直接打印，結(jié)果分析手寫。注2：對已給出（程序、結(jié)果及分析）的驗證性實驗，自己運行即可，可以不用寫在報告中。 對已給出（程序）的驗證性實驗，程序可以不用寫在實驗報告中，只寫出結(jié)果和分析。五

11、、實驗結(jié)論與總結(jié)（手寫）六、思考題（分析手寫）七、實驗參考資料附錄2：實驗所需部分函數(shù)及驗證性程序:1、序列的基本運算%1.單位取樣序列 x(n)=delta(n-n0) 要求n1<=n0<=n2 functionx,n=impseq(n0,n1,n2)n=n1:n2; x=(n-n0)=0; = 是邏輯判斷%2.單位階躍序列 x(n)=u(n-n0) 要求n1<=n0<=n2functionx,n=stepseq(n0,n1,n2)n=n1:n2; x=(n-n0)>=0;%3.信號加 y(n)=x1(n)+x2(n)%find函數(shù)：找出非零元素的索引號%x1:

12、第一個序列的值，n1:序列x1的索引號%x2:第二個序列的值，n2:序列x2的索引號functiony,n=sigadd(x1,n1,x2,n2)n=min(min(n1),min(n2):max(max(n1),max(n2);y1=zeros(1,length(n); y2=y1;y1(find(n>=min(n1)&(n<=max(n1)=1)=x1; y2(find(n>=min(n2)&(n<=max(n2)=1)=x2;y=y1+y2;%4.信號乘 y(n)=x1(n)*x2(n)functiony,n=sigmult(x1,n1,x2,n2

13、)n=min(min(n1),min(n2):max(max(n1),max(n2);y1=zeros(1,length(n); y2=y1;y1(find(n>=min(n1)&(n<=max(n1)=1)=x1;y2(find(n>=min(n2)&(n<=max(n2)=1)=x2;y=y1.*y2;%5.移位 y(n)=x(n-n0)functiony,n=sigshift(x,m,n0)n=m+n0; y=x;%6.翻褶 y(n)=x(-n)functiony,n=sigfold(x,n)y=fliplr(x); n=-fliplr(n);2、

14、序列的傅里葉變換%7. 求序列的傅里葉變換w=0:1:500*pi/500X=exp(j*w)./(exp(j*w)-0.5*ones(1,501)magX=abs(X)angX=angle(X)realX=real(X)imagX=imag(X)subplot(2,2,1)plot(w/pi,magX)gridxlabel('frequency in pi units')title('Magnitude Part')ylabel('Magnitude')subplot(2,2,3)plot(w/pi,angX)gridxlabel('f

15、requency in pi units')title('Angle Part')ylabel('Radians')subplot(2,2,2)plot(w/pi,realX)gridxlabel('frequency in pi units')title('Real Part')ylabel('Real')subplot(2,2,4)plot(w/pi,imagX)gridxlabel('frequency in pi units')title('Imaginary Part

16、9;)ylabel('Imaginary') 程序執(zhí)行結(jié)果：%8 令，繪制其傅立葉變換。用不同頻率對其進行采樣，分別畫出。Dt=0.00005; %步長為0.00005st=-0.005:Dt:0.005; xa=exp(-1000*abs(t); %取時間從-0.005s到0.005s這段模擬信號Wmax=2*pi*2000; %信號最高頻率為2*2000K=500; %頻域正半軸取500個點進行計算k=0:1:K;W=k*Wmax/K; % 求模擬角頻率Xa=xa*exp(-j*t'*W)*Dt; %計算連續(xù)時間傅立葉變換（利用矩陣運算實現(xiàn)） Xa=real(Xa)

17、; %取實部W=-fliplr(W),W(2:501); %將角頻率范圍擴展為從-到+Xa=fliplr(Xa),Xa(2:501); subplot(2,2,1);plot(t*1000,xa); %畫出模擬信號，橫坐標為時間（毫秒），縱坐標為幅度xlabel('time(millisecond)');ylabel('xa(t)'); title('anolog signal');subplot(2,2,2);plot(W/(2*pi*1000),Xa*1000); %畫出連續(xù)時間傅立葉變換 xlabel('frequency(kHZ)

18、'); %橫坐標為頻率（kHz）ylabel('xa(jw)'); %縱坐標為幅度title('FT');%下面為采樣頻率5kHz時的程序T=0.0002; %采樣間隔為n=-25:1:25;x=exp(-1000*abs(n*T); %離散時間信號K=500;k=0:1:K;w=pi*k/K; %w為數(shù)字頻率X=x*exp(-j*n'*w); %計算離散時間傅立葉變換（序列的傅立葉變換）X=real(X); w=-fliplr(w),w(2:K+1);X=fliplr(X),X(2:K+1);subplot(2,2,3);stem(n*T*10

19、00,x); %畫出采樣信號（離散時間信號）xlabel('time(millisecond)');ylabel('x1(n)');title('discrete signal');subplot(2,2,4);plot(w/pi,X); %畫出離散時間傅立葉變換xlabel('frequency(radian)'); %橫坐標為弧度ylabel('x1(jw)');title('DTFT');3、序列的傅里葉變換性質(zhì)分析%9 已知序列，求其傅里葉變換，并討論其傅里葉變換的周期性和對稱性。n=0:1

20、0x=(0.9*exp(j*pi/3).nk=-200:200X,w=ft1(x,n,k)magX=abs(X)angX=angle(X)subplot(2,1,1)plot(w/pi,magX)gridxlabel('frequency in pi units')ylabel('/X/')title('Magnitude Part')subplot(2,1,2)plot(w/pi,angX/pi)gridxlabel('frequency in pi units')ylabel('Radians/pi')titl

21、e('Angle Part')由圖可見，序列的傅里葉變換對是周期的，但不是共軛對稱的。%10、已知序列，求其傅里葉變換，并討論其傅里葉變換的周期性和對稱性。n=-5:5x=(-0.9).nk=-200:200X,w=ft1(x,n,k)magX=abs(X)angX=angle(X)subplot(2,1,1)plot(w/pi,magX)gridxlabel('frequency in pi units')ylabel('/X/')title('Magnitude Part')subplot(2,1,2)plot(w/pi,an

22、gX/pi)gridxlabel('frequency in pi units')ylabel('Radians/pi')title('Angle Part')由圖可見，序列的傅里葉變換對是周期的，是共軛對稱的。4、時域差分方程的求解采用filter函數(shù)實現(xiàn)線性常系數(shù)差分方程的遞推求解，函數(shù)調(diào)用格式如下：l yn=filter(B,A,xn) 計算輸入信號xn的零狀態(tài)響應(yīng)ynl yn=filter(B,A,xn,xi) 計算輸入信號xn的全響應(yīng)yn，xi為等效初始條件的輸入序列l(wèi) xi=filtic(B,A,ys,xs) 由初始條件計算xi的函數(shù)

23、4.1求解差分方程y(n)a1y(n-1)a2y(n-2)=b0x(n)b1x(n-1)的零狀態(tài)響應(yīng)和全響應(yīng)。已知X(n)為單位取樣序列，y(-1)=1,y(-2)=2,a1=0.5,a2=0.06,b0=2,b13。程序：xn=1 zeros(1,20) B=2,3 A=1,0.5,0.06ys=1,2xi=filtic(B,A,ys) yn1=filter(B,A,xn) yn2=filter(B,A,xn,xi) subplot(2,1,1)n1=0:length(yn1)-1stem(n1,yn1,'.')axis(0,21,-3,3)subplot(2,1,2)n2=