12月6日信號與系統(tǒng)實驗指導(dǎo)書

1、信信 號號 與與 系系 統(tǒng)統(tǒng) 實實驗驗指指導(dǎo)導(dǎo) 書書 長春工程學(xué)院電氣與信息學(xué)院長春工程學(xué)院電氣與信息學(xué)院 電工電子實驗教學(xué)中心電工電子實驗教學(xué)中心 二二八年二月十日八年二月十日 目目 錄錄 matlab 軟件仿真實驗部分軟件仿真實驗部分.1 實驗一 信號的基本運(yùn)算.1 實驗二 周期信號的傅里葉級數(shù)及 gibbs現(xiàn)象.12 實驗三 信號抽樣及信號重建 .19 實驗四 信號與系統(tǒng)復(fù)頻域分析 .26 el-ss-iii 實驗箱硬件實驗部分實驗箱硬件實驗部分.31 一、硬件資源.31 二、軟件安裝及使用 .34 三、實驗系統(tǒng)部分 .40 實驗一 濾波器.41 實驗二 一階電路的瞬態(tài)響應(yīng) .45 實驗

2、三 一階電路的零輸入響應(yīng)、零狀態(tài)響應(yīng)及完全響應(yīng).50 實驗四 二階電路的瞬態(tài)響應(yīng) .52 實驗五 二階網(wǎng)絡(luò)函數(shù)的模擬 .55 實驗六 方波信號的分解.59 實驗七 方波信號的合成.62 實驗八 抽樣定理.64 實驗九 數(shù)據(jù)采集.69 附錄一 實驗結(jié)果參考.70 附錄二 ad/da 卡調(diào)試說明.83 matlab 軟件仿真實驗部分軟件仿真實驗部分 實驗一實驗一 信號的基本運(yùn)算信號的基本運(yùn)算 一、實驗?zāi)康囊弧嶒災(zāi)康?1、熟悉掌握常用的用于信號與系統(tǒng)時域仿真分析的 matlab 函數(shù)。 2、掌握用 matlab 描述連續(xù)時間信號和離散時間信號的方法，能夠編寫 matlab 程序進(jìn)行 仿真。 3、熟

3、悉實現(xiàn)各種信號的時域變換和運(yùn)算的原理和方法，并在 matlab 環(huán)境下仿真。 4、利用延拓的方法將時限信號變成一個周期函數(shù)。 5、利用 matlab 的卷積工具實現(xiàn)兩個信號的卷積運(yùn)算。 二、實驗原理二、實驗原理 1、在信號與系統(tǒng)課程中，單位階躍信號 u(t) 和單位沖激信號 (t) 是二個非常有用的 信號。它們的定義如下 1.1(a) 0, 0)( 1)( tt dtt t 1.1(b) 0, 0 0, 1 )( t t tu 這里分別給出相應(yīng)的簡單的產(chǎn)生單位沖激信號和單位階躍信號的擴(kuò)展函數(shù)。產(chǎn)生單位沖激 信號的擴(kuò)展函數(shù)為： function y = delta(t) dt = 0.01; y

4、 = (u(t)-u(t-dt)/dt; 產(chǎn)生單位階躍信號的擴(kuò)展函數(shù)為： % unit step function function y = u(t) y = (t=0); % y = 1 for t 0, else y = 0 請將這二個 matlab 函數(shù)分別以 delta 和 u 為文件名保存在 work 文件夾中，以后，就可 以像教材中的方法使用單位沖激信號 (t) 和單位階躍信號 u(t)。 2、離散時間單位階躍信號 un定義為 1.2 0, 0 0, 1 n n nu 離散時間單位階躍信號 un除了也可以直接用前面給出的擴(kuò)展函數(shù)來產(chǎn)生，還可以利用 matlab 內(nèi)部函數(shù) ones(

5、1,n) 來實現(xiàn)。這個函數(shù)類似于 zeros(1,n)，所不同的是它產(chǎn)生的矩陣 的所有元素都為 1。 值得注意的是，利用 ones(1,n) 來實現(xiàn)的單位階躍序列并不是真正的單位階躍序列，而是 一個長度為 n 單位門(gate)序列，也就是 un-un-n。但是在一個有限的圖形窗口中，我們看 到的還是一個單位階躍序列。 3、信號的基本加法和乘法運(yùn)算 信號 f1 與 f2 之和（瞬時和）是指同一瞬時兩個信號之值對應(yīng)相加所構(gòu)成的“和信號”即 f3=f1+f2;信號 f1 與 f2 之積是指同一瞬時兩信號之值對應(yīng)相乘所構(gòu)成的“積信號”即 f3= f1*f2；離散序列相加（或相乘）可采用對應(yīng)樣點(diǎn)的值分

6、別相加（或相乘）的方法來計算。 用 matlab 程序仿真下面運(yùn)算：f1=sin(t)，f2=sin(t)，f3=f1+f2,f4=f1*f2;x=0 1 1 1 1 1, h=2 1 3 4 1 1,y=x+h,g=x.*h; 連續(xù)信號加法乘法實現(xiàn)程序 % program t=0:0.01:4*pi; f1=sin(t); f2= sin(t); f3=f1+f2; f4=f1.*f2; subplot(221); plot(t,f1); title(f1 signal); subplot(222); plot(t,f2); title(f2 signal); subplot(223); p

7、lot(t,f3); title(f1+f2 signal); subplot(224); plot(t,f4); title(f1*f2 signal); 運(yùn)行后的結(jié)果： 圖 1-1 程序運(yùn)行結(jié)果圖 離散序列加法乘法實現(xiàn)程序 x=0 1 1 1 1 1; h=2 1 3 4 1 1; y=x+h,g=x.*h; subplot(221); stem(x); title(x signal); subplot(222); stem(h); title(h signal); subplot(223); stem(y); title(x+h signal); subplot(224); stem(g

8、); title(x.*h signal); 運(yùn)行后的結(jié)果： 圖 1-2 程序運(yùn)行結(jié)果圖 4、信號的時移 信號的時移可用下面的數(shù)學(xué)表達(dá)式來描述： 設(shè)一個連續(xù)時間信號為 x(t)，它的時移 y(t) 表示為： y(t) = x(t - t0) 1.3 其中，t0為位移量。若 t0為正數(shù)，則 y(t)等于將 x(t)右移 t0秒之后的結(jié)果。反之，若 t0為 負(fù)數(shù)，則 y(t)等于將 x(t)左移 t0秒之后的結(jié)果。 在 matlab 中，時移運(yùn)算與數(shù)學(xué)上習(xí)慣表達(dá)方法完全相同。 對給定一個連續(xù)時間信號 x(t) = e-0.5tu(t)，對它分別左移 2 秒鐘和右移 2 秒鐘得到信號 x1(t) =

9、 e-0.5(t+2)u(t+2)和 x2(t) = e-0.5(t-2)u(t-2)。 實現(xiàn)程序： % program clear,close all, t = -5:0.01:5; x = exp(-0.5*t).*u(t); % generate the original signal x(t) x1 = exp(-0.5*(t+2).*u(t+2); % shift x(t) to the left by 2 second to get x1(t) x2 = exp(-0.5*(t-2).*u(t-2); % shift x(t) to the right by 2 second to

10、 get x2(t) subplot(311) plot(t,x) % plot x(t) grid on, title ( x = exp(-0.5*t).*u(t) subplot (312) plot (t,x1) % plot x1(t) grid on, title ( x1 = exp(-0.5*(t+2).*u(t+2) ) subplot (313) plot (t,x2) % plot x2(t) grid on, title ( x2 = exp(-0.5*(t-2).*u(t-2) xlabel (time t (sec) 程序運(yùn)行結(jié)果： 圖 1-3 程序運(yùn)行結(jié)果圖 注意

11、：在運(yùn)行上面的程序時，一定在所在的路徑下創(chuàng)建 u(t)的子函數(shù)： function y = u(t) y = (t=0); % y = 1 for t 0, else y = 0 保存名為 u.m 5、信號的時域反褶 對一個信號 xn的反褶運(yùn)算在數(shù)學(xué)上表示為 yn = x-n 1.4 這種反褶運(yùn)算，用 matlab 實現(xiàn)起來也是非常簡單的。有多種方法可以實現(xiàn)信號的反褶 運(yùn)算。 方法一方法一，修改繪圖函數(shù) plot(t,x)和 stem(n,x)中的時間變量 t 和 n，即用-t 和-n 替代原來的 t 和 n，這樣繪制出來的圖形，看起來就是原信號經(jīng)時域反褶后的版本。 方法二方法二，直接利用原信

12、號與其反褶信號的數(shù)學(xué)關(guān)系式來實現(xiàn)。這種方法最符合信號反褶運(yùn) 算的實際意義。 方法三方法三，使用 matlab 內(nèi)部函數(shù) fliplr()來實現(xiàn)信號的反褶運(yùn)算。其用法如下： y = fliplr(x)：其中 x 為原信號 x(t)或 xn，而 y 則為 x 的時域反褶。需要說明的是，函數(shù) fliplr()對信號作時域反褶，僅僅將信號中各個元素的次序作了一個反轉(zhuǎn)，這種反轉(zhuǎn)處理是獨(dú)立 于時間變量 t 和 n 的。因此，如果信號與其時間變量能夠用一個數(shù)學(xué)函數(shù)來表達(dá)的話，那么建 議將時間變量 t 和 n 的范圍指定在一個正負(fù)對稱正負(fù)對稱的時間區(qū)間即可。 分別編寫程序?qū)崿F(xiàn) m=sin(t);n=sin(-

13、t);xn=1 2 3 4;x-n,分析所畫信號波形， 程序如下： t=0:0.01:4*pi; n=0:1:3; m=sin(t); x=1 2 3 4; subplot(222); plot(t,m); title(sin(t) signal); subplot(221); plot(-t,m); title(sin(-t) signal); subplot(224); stem(n,x); title(xn signal); subplot(223); stem(-n,x); title(x-n signal); 程序運(yùn)行結(jié)果： 圖 1-4 程序運(yùn)行結(jié)果圖 6、信號的時域尺度變換 信號

14、x(t)的時域尺度變換在數(shù)學(xué)描述為 y(t) = x(at)， 1.5 其中 a 為任意常數(shù)。根據(jù) a 的不同取值，這種時域尺度變換對信號 x(t)具有非常不同的影 響。 當(dāng) a = 1 時，y(t) = x(t)； 當(dāng) a = -1 時，y(t) = x(-t)，即 y(t)可以通過將 x(t)反褶運(yùn)算而得到； 當(dāng) a 1 時，y(t) = x(at)，y(t)是將 x(t)在時間軸上的壓縮而得到； 當(dāng) 0 a 1 時，y(t) = x(at)，y(t)是將 x(t)在時間軸上的擴(kuò)展而得到； 當(dāng) -1 a 0 時，y(t) = x(at)，y(t)是將 x(t)在時間軸上的擴(kuò)展同時翻轉(zhuǎn)而得到

15、； 當(dāng) a 0 時 圖形右移， t x = 1 2 3 4; h = 4 3 2 1; y=conv(x,h) 在屏幕上得到顯示結(jié)果： y = 4 11 20 30 20 11 4 這表明，多項式 p1和 p2的乘積為： 411203020114 23456 3 ssssssp 正如前所述，用 matlab 處理連續(xù)時間信號時，獨(dú)立時間變量 t 的變化步長應(yīng)該是 很小的，假定用符號 dt 表示時間變化步長，那么，用函數(shù) conv()作兩個信號的卷積積分時， 應(yīng)該在這個函數(shù)之前乘以時間步長方能得到正確的結(jié)果。也就是說，正確的語句形式應(yīng)為：y = dt*conv(x,h)。 對于定義在不同時間段的

16、兩個時限信號 x(t)，t0 t t1，和 h(t)，t2 t t3。 如果用 y(t)來 表示它們的卷積結(jié)果，則 y(t)的持續(xù)時間范圍要比 x(t)或 h(t)要長，其時間范圍為 t0+t2 t t1+t3。這個特點(diǎn)很重要，利用這個特點(diǎn)，在處理信號在時間上的位置時，可以很容易地將信號 的函數(shù)值與時間軸的位置和長度關(guān)系保持一致性。 根據(jù)給定的兩個連續(xù)時間信號 x(t) = tu(t)-u(t-1)和 h(t) = u(t)-u(t-1)，編寫程序，完成這兩 個信號的卷積運(yùn)算，并繪制它們的波形圖。范例程序如下： % program t0 = -2; t1 = 4; dt = 0.01; t =

17、 t0:dt:t1; x = u(t)-u(t-1); h = t.*(u(t)-u(t-1); y = dt*conv(x,h); % compute the convolution of x(t) and h(t) subplot(221) plot(t,x), grid on, title(signal x(t), axis(t0,t1,-0.2,1.2) subplot(222) plot(t,h), grid on, title(signal h(t), axis(t0,t1,-0.2,1.2) subplot(212) t = 2*t0:dt:2*t1; % again speci

18、fy the time range to be suitable to the % convolution of x and h. plot(t,y), grid on, title(the convolution of x(t) and h(t), axis(2*t0,2*t1,-0.1,0.6), xlabel(time t sec) 程序運(yùn)行結(jié)果 圖 1-7 程序運(yùn)行結(jié)果圖 在有些時候，做卷積和運(yùn)算的兩個序列中，可能有一個序列或者兩個序列都非常長，甚至 是無限長，matlab 處理這樣的序列時，總是把它看作是一個有限長序列，具體長度由編程 者確定。實際上，在信號與系統(tǒng)分析中所遇到的無限長

19、序列，通常都是滿足絕對可和或絕對可 積條件的信號。因此，對信號采取這種截短處理盡管存在誤差，但是通過選擇合理的信號長度， 這種誤差是能夠減小到可以接受的程度的。若這樣的一個無限長序列可以用一個數(shù)學(xué)表達(dá)式表 示的話，那么，它的長度可以由編程者通過指定時間變量 n 的范圍來確定。 例如，對于一個單邊實指數(shù)序列 xn = 0.5nun，通過指定 n 的范圍為 0 n 100，則對應(yīng) 的 xn的長度為 101 點(diǎn)，雖然指定更寬的 n 的范圍，xn將與實際情況更相符合，但是，注意 到，當(dāng) n 大于某一數(shù)時，xn之值已經(jīng)非常接近于 0 了。對于序列 xn = 0.5nun，當(dāng) n = 7 時， x7 =

20、0.0078，這已經(jīng)是非常小了。所以，對于這個單邊實指數(shù)序列，指定更長的 n 的范圍是 沒有必要的。當(dāng)然，不同的無限長序列具有不同的特殊性，在指定 n 的范圍時，只要能夠反映 序列的主要特征就可以了。 9、關(guān)于 matlab 工具在信號處理中應(yīng)用的補(bǔ)充 在繪制信號的波形圖時，有時我們需要將若干個圖形繪制在圖一個圖形窗口中，這就需要 使用 matlab 的圖形分割函數(shù) subplot()，其用法是在繪圖函數(shù) stem 或 plot 之前，使用圖形 分割函數(shù) subplot(n1,n2,n3)，其中的參數(shù) n1，n2 和 n3 的含義是，該函數(shù)將把一個圖形窗口分 割成 n1xn2 個子圖，即將繪制

21、的圖形將繪制在第 n3 個子圖中。 常用的圖形控制函數(shù) axis(xmin,xmax,ymin,ymax)：圖型顯示區(qū)域控制函數(shù)，其中 xmin 為橫軸的顯示起點(diǎn)，xmax 為橫軸的顯示終點(diǎn)，ymin 為縱軸的顯示起點(diǎn)，ymax 為縱軸的顯示 終點(diǎn)。 有時，為了使圖形具有可讀性，需要在所繪制的圖形中，加上一些網(wǎng)格線來反映信號的幅 度大小。matlab 中的 grid on/grid off 可以實現(xiàn)在你的圖形中加網(wǎng)格線。 grid on：在圖形中加網(wǎng)格線。 grid off：取消圖形中的網(wǎng)格線。 x = input(type in signal x(t) in closed form:) 三、

22、實驗內(nèi)容及步驟三、實驗內(nèi)容及步驟 實驗前，必須首先閱讀本實驗原理，讀懂所給出的全部范例程序。實驗開始時，先在計算 機(jī)上運(yùn)行這些范例程序，觀察所得到的信號的波形圖。并結(jié)合范例程序應(yīng)該完成的工作，進(jìn)一 步分析程序中各個語句的作用，從而真正理解這些程序。 實驗前，一定要針對下面的實驗項目做好相應(yīng)的實驗準(zhǔn)備工作，包括事先編寫好相應(yīng)的實 驗程序等事項。 練習(xí)練習(xí) 1：結(jié)合編寫的階躍函數(shù)編寫一門函數(shù)(門寬為 4 幅度為 1)寫出程序并會出門函 數(shù)信號的波形； 練習(xí)練習(xí) 2：結(jié)合實驗原理的信號的基本運(yùn)算的程序，編寫程序繪制 m=sin(t);g=sin(2t-pi/2); xn=1 2 5 6 3 ;x(1

23、/2)n-1;的四個信號的波形，并分析圖形的變換過程。 練習(xí)練習(xí) 3：根據(jù)周期函數(shù)的定義以及實驗原理中延拓的方法實現(xiàn)脈沖函數(shù)（脈沖寬度 1 周期 2 幅度 1）并繪制其圖形； 練習(xí)練習(xí) 4：利用 matlab 內(nèi)部所帶的卷積工具對兩個門函數(shù)進(jìn)行卷積，分別繪制出兩 個門函數(shù)和卷積后的波形，并分析門函數(shù)卷積的規(guī)律， （一個門函數(shù)門寬為 1 幅度 1 另一個門函數(shù)門寬為 2 幅度為 1，進(jìn)行兩者卷積） 四、實驗報告要求四、實驗報告要求 1、按要求完整書寫你所編寫的全部 matlab 程序 2、詳細(xì)記錄實驗過程中的有關(guān)信號波形圖（存于自帶的 u 盤中） ，圖形要有明確的標(biāo)題。 全部的 matlab 圖

24、形應(yīng)該用打印機(jī)打印，然后貼在本實驗報告中的相應(yīng)位置，禁止復(fù)印件。 3、實事求是地回答相關(guān)問題，嚴(yán)禁抄襲。 實驗二實驗二 周期信號的傅里葉級數(shù)及周期信號的傅里葉級數(shù)及 gibbs 現(xiàn)象現(xiàn)象 一、實驗?zāi)康囊?、實驗?zāi)康?1、掌握連續(xù)時間周期信號的傅里葉級數(shù)的物理意義和分析方法； 2、觀察截短傅里葉級數(shù)而產(chǎn)生的“gibbs 現(xiàn)象” ，了解其特點(diǎn)以及產(chǎn)生的原因； 3、掌握連續(xù)周期時間信號傅里葉變換的分析方法及其物理意義； 二、實驗原理及方法二、實驗原理及方法 1、連續(xù)時間周期信號的傅里葉級數(shù) ctfs 分析 任何一個周期為 t1的正弦周期信號，只要滿足狄利克利狄利克利條件，就可以展開成傅里葉級數(shù)。 其中

25、三角傅里葉級數(shù)三角傅里葉級數(shù)為： 2.1 1 000 )sin()cos()( k kk tkbtkaatx 或： 2.2 1 00 )cos()( k kk tkcatx 其中，稱為信號的基本頻率基本頻率（fundamental frequency） ，分別是 1 0 2 t kk baa，和， 0 信號的直流分量直流分量、余弦分量幅度余弦分量幅度和正弦分量幅度正弦分量幅度，為合并同頻率項之后各正弦諧)(tx kk c、 波分量的幅度和初相位，它們都是頻率的函數(shù)，繪制出它們與之間的圖像，稱為信號 0 k 0 k 的頻譜圖頻譜圖（簡稱“頻譜頻譜” ） ，圖像為幅度譜幅度譜，圖像為相位譜相位譜。

26、 k c 0 k k 0 k 三角形式傅里葉級數(shù)表明，如果一個周期信號 x(t)，滿足狄里克利條件，那么，它就可以 被看作是由很多不同頻率的互為諧波關(guān)系諧波關(guān)系（harmonically related）的正弦信號所組成，其中每 一個不同頻率的正弦信號稱為正弦諧波分量正弦諧波分量 (sinusoid component)，其幅度幅度（amplitude）為 。也可以反過來理解三角傅里葉級數(shù)：用無限多個正弦諧波分量可以合成一個任意的非正 k c 弦周期信號。 指數(shù)形式的傅里葉級數(shù)指數(shù)形式的傅里葉級數(shù)為： 2.3 k tjk ke atx 0 )( 其中，為指數(shù)形式的傅里葉級數(shù)的系數(shù)，按如下公式計

27、算： k a 2.4 2/ 2/1 1 1 0 )( 1 t t tjk k dtetx t a 指數(shù)形式的傅里葉級數(shù)告訴我們，如果一個周期信號 x(t)，滿足狄里克利條件，那么，它 就可以被看作是由很多不同頻率的互為諧波關(guān)系（harmonically related）的周期復(fù)指數(shù)信號所 組成，其中每一個不同頻率的周期復(fù)指數(shù)信號稱為基本頻率分量，其復(fù)幅度（complex amplitude）為。這里“復(fù)幅度（complex amplitude） ”指的是通常是復(fù)數(shù)。 k a k a 上面的傅里葉級數(shù)的合成式說明，我們可以用無窮多個不同頻率的周期復(fù)指數(shù)信號來合成 任意一個周期信號。然而，用計算機(jī)

28、（或任何其它設(shè)備）合成一個周期信號，顯然不可能做到 用無限多個諧波來合成，只能取這些有限個諧波分量來近似合成。 假設(shè)諧波項數(shù)為 n，則上面的和成式為： 2.5 n nk tjk ke atx 0 )( 顯然，n 越大，所選項數(shù)越多，有限項級數(shù)合成的結(jié)果越逼近原信號 x(t)。本實驗可以比 較直觀地了解傅里葉級數(shù)的物理意義，并觀察到級數(shù)中各頻率分量對波形的影響包括 “gibbs”現(xiàn)象：即信號在不連續(xù)點(diǎn)附近存在一個幅度大約為 9%的過沖，且所選諧波次數(shù)越多， 過沖點(diǎn)越向不連續(xù)點(diǎn)靠近。這一現(xiàn)象在觀察周期矩形波信號和周期鋸齒波信號時可以看得很清 楚。 例題 1、給定一個周期為 t1 = 2s 的連續(xù)時

29、間周期方波信號，如圖所示，其一個周期內(nèi)的 數(shù)學(xué)表達(dá)式為： 21, 0 10, 1 )( 1 t t tx 解：首先，我們根據(jù)前面所給出的公式，計 算該信號的傅里葉級數(shù)的系數(shù)。 2/ 2/ 1 1 1 1 0 )( 1 t t tjk k dtetx t a 1 0 0 2 1 dte tjk 1 0 0 0 )( 2 1 0 tjkde kj tjk 0 1 0 2 0 kj e tjk 0 2 1 0 kj e jk 0 22 2 2 00 0 kj ee e k j k j k j 0 2 0 0) 2 sin( k j e k k 因為：0 = 2/t1 = ，代入上式得到： k k j

30、a k k ) 2 sin( )( 在 matlab 命令窗口，依次鍵入： k = -10:10; 1 12 -1-2 x(t) t 0 圖圖 2.1 周期方波信號周期方波信號 ak = (-j).k).* (sin(k+eps)*pi/2)./(k+eps)*pi) % the expression of ak ak = columns 1 through 4 -0.0000 0 + 0.0354i -0.0000 0 + 0.0455i columns 5 through 8 -0.0000 0 + 0.0637i -0.0000 0 + 0.1061i columns 9 through

31、 12 -0.0000 0 + 0.3183i 0.5000 0 - 0.3183i columns 13 through 16 -0.0000 0 - 0.1061i -0.0000 0 - 0.0637i columns 17 through 20 -0.0000 0 - 0.0455i -0.0000 0 - 0.0354i column 21 -0.0000 從 matlab 命令窗口，我們得到了該周期信號從到共 21 個系數(shù)。 10 a 10 a 緊接著再鍵入以下命令： subplot(221) stem(k,abs(ak),k.) title(the fourier series

32、coefficients) xlabel(frequency index k) 就得到一幅如右圖所示的描述與 k 之間的關(guān)系的 k a 圖形。 以上是我們通過手工計算得到的這個周期信號的傅 里葉級數(shù)表達(dá)式及其頻譜圖，下面給出完成傅里葉級數(shù) 系數(shù)計算的相應(yīng) matlab 范例程序。 % program2_1 % this program is used to evaluate the fourier series coefficients ak of a periodic square wave clear, close all t = 2; dt = 0.00001; t = -2:dt:2;

33、 x1 = u(t) - u(t-1-dt); x = 0; for m = -1:1 % periodically extend x1(t) to form a periodic signal x = x + u(t-m*t) - u(t-1-m*t-dt); end w0 = 2*pi/t; n = 10; % the number of the harmonic components l = 2*n+1; for k = -n: n; % evaluate the fourier series coefficients ak ak(n+1+k) = (1/t)*x1*exp(-j*k*w

34、0*t)*dt; end phi = angle(ak); % evaluate the phase of ak 執(zhí)行程序 program2_1 后，就完成了信號的傅里葉級數(shù)的系數(shù)的計算，在命令窗口鍵入 ak 命令窗口就可以顯示傅里葉級數(shù)的 21 個系數(shù)： ak = columns 1 through 4 0.0000 + 0.0000i 0.0000 + 0.0354i 0.0000 - 0.0000i 0.0000 + 0.0455i columns 5 through 8 0.0000 - 0.0000i 0.0000 + 0.0637i 0.0000 - 0.0000i 0.0000

35、+ 0.1061i columns 9 through 12 0.0000 - 0.0000i 0.0000 + 0.3183i 0.5000 0.0000 - 0.3183i columns 13 through 16 0.0000 + 0.0000i 0.0000 - 0.1061i 0.0000 + 0.0000i 0.0000 - 0.0637i columns 17 through 20 0.0000 + 0.0000i 0.0000 - 0.0455i 0.0000 + 0.0000i 0.0000 - 0.0354i column 21 0.0000 - 0.0000i 將這里的

36、 ak 之值同前面手工計算得到的 ak 比較，可見兩者是完全相同的。 再次特別提示：程序中，時間變量的變化步長 dt 的大小對傅里葉級數(shù)系數(shù)的計算精度的 影響非常大，dt 越小，精度越高，本程序中的 dt 之所以選擇 0.00001 就是為了提高計算精度。 但是，計算機(jī)所花的計算時間越長。 在程序 program2_1 中添加相應(yīng)的計算| ak |和繪圖語句，就可以繪制出信號的幅度譜和相 位譜的譜線圖。 2、周期信號的合成以及 gibbs 現(xiàn)象 從傅里葉級數(shù)的合成式合成式（synthesis equation） k tjk ke atx 0 )( 可以看出，用無窮多個不同頻率和不同振幅的周期復(fù)

37、指數(shù)信號可以合成一個周期信號。然 而，我們無法用計算機(jī)實現(xiàn)對無窮多個周期復(fù)指數(shù)信號的合成。但是，用有限項來合成卻是可 行的，在實際應(yīng)用中，多半也就是這么做的。然而，這樣做的一個必然結(jié)果，就是引入了誤差。 如果一個周期信號在一個周期有內(nèi)斷點(diǎn)存在，那么，引入的誤差將除了產(chǎn)生紋波之外，還 將在斷點(diǎn)處產(chǎn)生幅度大約為 9%的過沖過沖（overshot） ，這種現(xiàn)象被稱為吉伯斯現(xiàn)象吉伯斯現(xiàn)象（gibbs phenomenon） 。 為了能夠觀察到合成信號與原信號的不同以及 gibbs 現(xiàn)象，我們可以利用前面已經(jīng)計算出 的傅里葉級數(shù)的系數(shù)，計算出截短的傅里葉級數(shù)： n nk tjk ke atx 0 )(

38、這個計算可用 l = 2n+1 次循環(huán)來完成： tnrj k eraxx 0 )1( )(22 其中 r 作為循環(huán)次數(shù)，x2 在循環(huán)之前應(yīng)先清零。完成這一計算的 matlab 程序為： x2 = 0; l = 2*n+1; for r = 1:l; x2 = x2+ak(r)*exp(j*(r-1-n)*w0*t); end; 完成了所有的計算之后，就可以用繪圖函數(shù)：plot()和 stem()將計算結(jié)果包括 x1, x2, abs(ak)和 angle(ak)以圖形的形式給出，便于我們觀察。 觀察吉伯斯現(xiàn)象的最好的周期信號就是圖 2-1 所示的周期方波信號，這種信號在一個周期 內(nèi)有兩個斷點(diǎn)，

39、用有限項級數(shù)合成這個信號時，吉伯斯現(xiàn)象的特征非常明顯，便于觀察。 例題 2：修改程序 program2_1，使之能夠用有限項級數(shù)合成例題 2-1 所給的周期方波信號， 并繪制出原始周期信號、合成的周期信號、信號的幅度譜和相位譜。 為此，只要將前述的 for 循環(huán)程序段和繪圖程序段添加到程序 program2_1 中即可，范例 程序如下： % program2_2 % this program is used to compute the fourier series coefficients ak of a periodic square wave clear,close all t = 2;

40、 dt = 0.00001; t = -2:dt:2; x1 = u(t)-u(t-1-dt); x = 0; for m = -1:1 x = x + u(t-m*t) - u(t-1-m*t-dt); % periodically extend x1(t) to form a periodic signal end w0 = 2*pi/t; n = input(type in the number of the harmonic components n = :); l = 2*n+1; for k = -n:1:n; ak(n+1+k) = (1/t)*x1*exp(-j*k*w0*t)

41、*dt; end phi = angle(ak); y=0; for q = 1:l; % synthesiz the periodic signal y(t) from the finite fourier series y = y+ak(q)*exp(j*(-(l-1)/2+q-1)*2*pi*t/t); end; subplot(221), plot(t,x), title(the original signal x(t), axis(-2,2,-0.2,1.2), subplot(223), plot(t,y), title(the synthesis signal y(t), axi

42、s(-2,2,-0.2,1.2), xlabel(time t), subplot(222) k=-n:n; stem(k,abs(ak),k.), title(the amplitude |ak| of x(t), axis(-n,n,-0.1,0.6) subplot(224) stem(k,phi,r.), title(the phase phi(k) of x(t), axis(-n,n,-2,2), xlabel(index k) n=12 時，程序運(yùn)行結(jié)果： 圖 2-1 程序運(yùn)行結(jié)果 在用這個程序觀察吉伯斯現(xiàn)象時，可以反復(fù)執(zhí)行該程序，每次執(zhí)行時，輸入不同之 n 值， 比較所的圖形的

43、區(qū)別，由此可以觀察到吉伯斯現(xiàn)象的特征。 三、實驗內(nèi)容和要求三、實驗內(nèi)容和要求 實驗前，必須首先閱讀本實驗原理，讀懂所給出的全部范例程序。實驗開始時，先在計算 機(jī)上運(yùn)行這些范例程序，觀察所得到的信號的波形圖。并結(jié)合范例程序應(yīng)該完成的工作，進(jìn)一 步分析程序中各個語句的作用，從而真正理解這些程序。 實驗前，一定要針對下面的實驗項目做好相應(yīng)的實驗準(zhǔn)備工作，包括事先編寫好相應(yīng)的實 驗程序等事項。 練習(xí)練習(xí) 1、周期信號的傅里葉級數(shù) 給定如下兩個周期信號： t 11 1 22 )( 1 tx （1） 、手工計算 x1(t)傅里葉級數(shù)的系數(shù)： 信號 x1(t) 在其主周期內(nèi)的數(shù)學(xué)表達(dá)式為： 計算 x1(t)

44、 的傅里葉級數(shù)的系數(shù)的計算過程如下： 僅供參考 x1(t)=(t+1)*(u(t+1)-u(t)+(1-t)*(u(t)-u(t-1) 2/ 2/ 1 1 1 1 0 )( 1 t t tjk k dtetx t a 計算小程序： k = -10:10; ak = 的表達(dá)式 % the expression of ak k a 通過計算得到的 x1(t)的傅里葉級數(shù)的系數(shù)的數(shù)學(xué)表達(dá)式是： （2） 、用 matlab 計算的傅里葉級數(shù)的系數(shù) ak 從-10 到 10 共 21 個系數(shù)。仿照程序 program2_1，編寫程序以計算 x1(t)的傅里葉級數(shù)的系數(shù)。程序如下： 執(zhí)行程序后，就完成了信

45、號的傅里葉級數(shù)的系數(shù)的計算，在命令窗口鍵入 ak 命令窗口就可以顯示傅里葉級數(shù)的 21 個系數(shù)： （3） 、通過執(zhí)行程序所得到的 x1(t)的傅里葉級數(shù)的 ak 從-10 到 10 共 21 個系數(shù) 與你手工計算的 ak 相比較，是否相同，如有不同，是何原因造成的？ 答答： 練習(xí)練習(xí) 2、反復(fù)執(zhí)行程序 program2_2，每次執(zhí)行該程序時，輸入不同的 n 值，并觀察所合成 的周期方波信號。分析吉伯斯現(xiàn)象的特點(diǎn)，觀察合成的信號波形中，是否會產(chǎn)生 gibbs 現(xiàn)象？ 為什么？； 答答： 四、實驗報告要求四、實驗報告要求 1、按要求完整書寫你所編寫的全部 matlab 程序 2、詳細(xì)記錄實驗過程中

46、的有關(guān)信號波形圖（存于自帶的 u 盤中） ，圖形要有明確的標(biāo)題。 全部的 matlab 圖形應(yīng)該用打印機(jī)打印，然后貼在本實驗報告中的相應(yīng)位置，禁止復(fù)印件。 3、實事求是地回答相關(guān)問題，嚴(yán)禁抄襲。 本實驗完成時間： 年 月 日 )( 1 tx 實驗三實驗三 信號抽樣及信號重建信號抽樣及信號重建 一、實驗?zāi)康囊?、實驗?zāi)康?1、進(jìn)一步理解信號的抽樣及抽樣定理； 2、進(jìn)一步掌握抽樣信號的頻譜分析； 3、掌握和理解信號抽樣以及信號重建的原理； 二、實驗原理及方法二、實驗原理及方法 1、信號的抽樣及抽樣定理 抽樣（抽樣（sampling） ，就是從連續(xù)時間信號中抽取一系列的信號樣本，從而，得到一個離散 時

47、間序列（discrete-time sequence） ， 圖 3-1 給出了信號理想抽樣的原理圖： )(tx )(tp )(txs )(jx m m 圖圖 3-1 (a) 抽樣原理圖，抽樣原理圖，(b) 帶限信號的頻譜帶限信號的頻譜 (a)(b) 上圖中，假設(shè)連續(xù)時間信號是一個帶限信號（bandlimited signal） ，其頻率范圍為 ，抽樣脈沖為理想單位沖激串（單位沖激串（unit impulse train） ，其數(shù)學(xué)表達(dá)式為： mm 3.1 )()( s ntttp 由圖可見，模擬信號 x(t)經(jīng)抽樣后，得到已抽樣信號（已抽樣信號（sampled signal）xs(t)，且：

48、3.2)()()(tptxtxs 將 p(t)的數(shù)學(xué)表達(dá)式代入上式得到： 3.3 )()()( sss nttntxtx 顯然，已抽樣信號 xs(t) 也是一個沖激串，只是這個沖激串的沖激強(qiáng)度被 x(nts) 加權(quán)了。 從頻域上來看，p(t) 的頻譜也是沖激序列，且為： 3.4 )()( ss ntpf 根據(jù)傅里葉變換的頻域卷積定理，時域兩個信號相乘，對應(yīng)的積的傅里葉變換等于這兩 個信號的傅里葉變換之間的卷積。所以，已抽樣信號 xs(t)的傅里葉變換為： 3.5 n s s s njx t jx)( 1 )( 表達(dá)式 4.5 告訴我們，如果信號 x(t)的傅里葉變換為 x(j)，則已抽樣信號

49、xs(t) 的傅里葉 變換 xs(j)等于無窮多個加權(quán)的移位的 x(j)之和，或者說，已抽樣信號的頻譜等于原連續(xù)時 間信號的頻譜以抽樣頻率s為周期進(jìn)行周期復(fù)制的結(jié)果。如圖 3-2 所示： )(tx )(tp )(txs )(jx )(jp )(jxs s s s m m s t t t t s t/1 圖圖 3-2 信號抽樣及其頻譜圖信號抽樣及其頻譜圖 3.6 在 matlab 中，對信號抽樣的仿真， 例題 設(shè)連續(xù)時間信號為一個正弦信號 x(t) = cos(0.5t)，抽樣周期為 ts = 1/4 秒，編程序 繪制信號 x(t)和已抽樣信號 xn的波形圖。 范例程序 sampling 如下：

50、 % sampling clear, close all, t = 0:0.01:10; ts = 1/4; % sampling period n = 0:ts:10; % make the time variable to be discrete x = cos(0.5*pi*t); xn = cos(0.5*pi*n); % sampling subplot(221) plot(t,x), title(a continuous-time signal x(t), xlabel(time t) subplot(222) stem(n,xn,.), title(the sampled ver

51、sion xn of x(t), xlabel(time index n) 執(zhí)行該程序后，得到的波形圖如圖 3-3 所示。 圖圖 3-3 連續(xù)時間信號及其抽樣后的離散時間序列連續(xù)時間信號及其抽樣后的離散時間序列 在這個范例程序中，先將連續(xù)時間 t 進(jìn)行離散化，使之成為以 ts = 1/4 秒的離散時間 n， 然后，將 n 代入到信號 x(t) 的數(shù)學(xué)表達(dá)式中計算，就完成了抽樣過程，且得到了抽樣后的離 散時間序列 xn。 2、 信號抽樣過程中的頻譜混疊信號抽樣過程中的頻譜混疊 為了能夠觀察到已抽樣信號的頻譜是否會存在混疊現(xiàn)象，或者混疊程度有多么嚴(yán)重，有 必要計算并繪制出已抽樣信號的傅里葉變換。

52、根據(jù)式 3.5 可計算出已抽樣信號的頻譜。其中，主要利用了一個 for 循環(huán)程序完成周期延 拓運(yùn)算。 % program clear, close all, tmax = 4; dt = 0.01; t = 0:dt:tmax; ts = 1/10; ws = 2*pi/ts; w0 = 20*pi; dw = 0.1; w = -w0:dw:w0; n = 0:1:tmax/ts; x = exp(-4*t).*u(t); xn = exp(-4*n*ts); subplot(221) plot(t,x), title(a continuous-time signal x(t), xlabe

53、l(time t), axis(0,tmax,0,1), grid on subplot(223) stem(n,xn,.), title(the sampled version xn of x(t), xlabel(time index n), axis(0,tmax/ts,0,1), grid on xa = x*exp(-j*t*w)*dt; x = 0; for k = -8:8; x = x + x*exp(-j*t*(w-k*ws)*dt; end subplot(222) plot(w,abs(xa) title(magnitude spectrum of x(t), grid

54、on axis(-60,60,0,1.8*max(abs(xa) subplot(224) plot(w,abs(x) title(magnitude spectrum of xn), xlabel(frequency in radians/s),grid on axis(-60,60,0,1.8*max(abs(xa) 程序運(yùn)行結(jié)果： 圖 3-4 程序運(yùn)行結(jié)果圖 本程序可以用來觀察在不同的抽樣頻率條件下，已抽樣信號的頻譜的混疊程度，從而更 加牢固地理解抽樣定理。但是，提請注意的是，在 for 循環(huán)程序段中，計算已抽樣信號的頻譜 x 時，沒有乘以系數(shù) 1/ts，是為了便于比較 x 與 xa 之

55、間的區(qū)別，從而方便觀察頻譜的混疊程 度。另外，程序中的時間步長 dt 的選擇應(yīng)該與抽樣周期 ts 保持一定的比例關(guān)系，建議 ts 不 應(yīng)小于 10dt，否則，計算得到的已抽樣信號的頻譜將出現(xiàn)錯誤。 3、 信號重建信號重建 如果滿足抽樣定理，那么，我們就可以唯一地由已抽樣信號 xn 恢復(fù)出原連續(xù)時間信號 x(t)。在理想情況下，可以將離散時間序列通過一個理想低通濾波器，圖 3.5 給出了理想情況 下信號重建的原理示意圖。 )(tx )(tp )(txp 圖圖 3-5 信號重建原理圖信號重建原理圖 ideal lowpass filter )(txr 理想低通濾波器也稱重建濾波器，它的單位沖激響應(yīng)

56、 3.7 t tt th c cc )sin( )( 已抽樣信號 xp(t)的數(shù)學(xué)表達(dá)式為：，根據(jù)系統(tǒng)輸入輸出的卷 )()()(nttntxtxp 積表達(dá)式，我們有 3.8)()()(thtxtx pr 將 xp(t)代入式 4.8，有 3.9 n c cc r ntt nttt ntxtx )( )(sin( )()( 這個公式稱為內(nèi)插公式（內(nèi)插公式（interpolation formula） ，該公式的推導(dǎo)詳見教材，請注意復(fù)習(xí) 有關(guān)內(nèi)容。須提請注意的是，這里的內(nèi)插公式是基于重建濾波器為理想低通濾波器的。若重建 濾波器不是理想低通濾波器，則不能用這個內(nèi)插公式。理想低通濾波器的頻率響應(yīng)特性曲

57、線和 其單位沖激響應(yīng)曲線如圖 3.6 所示。 )(jh t c c )(th c t t 圖圖 3-6 理想低通濾波器的幅度頻率響應(yīng)和單位沖激響應(yīng)理想低通濾波器的幅度頻率響應(yīng)和單位沖激響應(yīng) 范例程序就是根據(jù)這個內(nèi)插公式來重構(gòu)原始信號。本程序已經(jīng)做了較為詳細(xì)的注釋，請 結(jié)合教材中的內(nèi)插公式仔細(xì)閱讀本程序，然后執(zhí)行，以掌握和理解信號重建的基本原理。范例 程序如下。 % program % signal sampling and reconstruction % the original signal is the raised cosin signal: x(t) = 1+cos(pi*t).*u

58、(t+1)-u(t-1). clear; close all, wm = 2*pi; % the highest frequency of x(t) a = input(type in the frequency rate ws/wm=:); % ws is the sampling frequency wc = wm; % the cutoff frequency of the ideal lowpass filter t0 = 2; t = -t0:0.01:t0; x = (1+cos(pi*t).*(u(t+1)-u(t-1); subplot(221); % plot the ori

59、ginal signal x(t) plot(t,x); grid on, axis(-2,2,-0.5,2.5); title(original signal x(t);xlabel(time t); ws = a*wm; % sampling frequency ts = 2*pi/ws; % sampling period n = fix(t0/ts); % determine the number of samplers n = -n:n; nts = n*ts; % the discrete time variable xs = (1+cos(pi*nts).*(u(nts+1)-u

60、(nts-1); % the sampled version of x(t) subplot(2,2,2) % plot xs stem(n,xs,.); xlabel(time index n); grid on, title(sampled version xn); xr = zeros(1,length(t); % specify a memory to save the reconstructed signal l = length(-n:n); xa = xr; figure(2); % open a new figure window to see the demo of sign