試驗(yàn)三周期信號(hào)頻域分析

1、實(shí)驗(yàn)三周期信號(hào)頻域分析一、目的(1)掌握周期信號(hào)傅立葉級(jí)數(shù)分解與合成的計(jì)算公式(2)掌握利用MATLAB實(shí)現(xiàn)周期信號(hào)傅立葉級(jí)數(shù)分解與綜合方法(3)理解并掌握周期信號(hào)頻譜特點(diǎn)二、周期信號(hào)傅立葉級(jí)數(shù)周期信號(hào)是定義在(*,g)區(qū)間內(nèi)，按一定時(shí)間問(wèn)隔(周期 T)不斷重復(fù)的信號(hào)?？杀硇閒(t) u f (t mT)式中m為任意整數(shù)，T為周期，周期的倒數(shù)成為該信號(hào)頻率。1、連續(xù)時(shí)間周期信號(hào)的分解設(shè)有周期信號(hào)f(t),周期為T,角頻率。1 = 2肛=至，且滿足狄里赫利條件，則該Ti周期信號(hào)可以展開(kāi)成傅立葉級(jí)數(shù)，即可表示為一系列不同頻率的正弦或復(fù)指數(shù)信號(hào)之和。傅立葉級(jí)數(shù)有三角形式和指數(shù)形式兩種：1)三角形式

2、傅立葉級(jí)數(shù)三角形式傅立葉級(jí)數(shù)為f (t) = a。 a1 c o s h s i n it a2cos2t b? s i nt . ： a c o st bn s i n nt .0000二a。一： a” c o si( J) 一: bn si m(4)n七式中系數(shù)按如下公式求解n1 T1=-21 f (t)dtT1 2 T1=T1 f (t)cosn ；tdt T2 T1=：i f (t) sin n、tdtT12)指數(shù)形式傅立葉級(jí)數(shù)指數(shù)形式的傅立葉級(jí)數(shù)表達(dá)式為：f(t)=EFne陽(yáng) n=0,1 苫津3,n -：式中E稱為傅立葉復(fù)系數(shù)，可由下式求得T1Fn = I f (t)exdtT1(3

3、-1)(3-2)(3-3)(3-4)(3-5)(3-6)傅立葉級(jí)數(shù)的指數(shù)形式和三角形式是等價(jià)的，具系數(shù)可相互轉(zhuǎn)換，具體參見(jiàn)信號(hào)與系 統(tǒng)教科書(shū)。2、連續(xù)時(shí)間周期信號(hào)的傅立葉綜合任何滿足狄里赫利條件的周期信號(hào)，可以表示成式(3-1)或(3-5)的和式形式, 式(3-1)或(3-5)成為連續(xù)時(shí)間周期信號(hào)(CTFS)綜合公式。一般說(shuō)來(lái)，傅立葉級(jí)數(shù)系數(shù)有無(wú)限個(gè)非零值，即任何具有有限個(gè)間斷點(diǎn)的周期信號(hào) 都一定有一個(gè)無(wú)限項(xiàng)非零系數(shù)的傅立葉級(jí)數(shù)表示。但對(duì)于數(shù)值計(jì)算來(lái)說(shuō)，這是無(wú)法實(shí)現(xiàn) 的。在實(shí)際應(yīng)用中，可以用有限項(xiàng)傅立葉級(jí)數(shù)求和來(lái)逼近。即：15f (t) = Fd8 =3 + an cos(nt) + bn s

4、in(n缸t)(3-7)n _Nn 1_n 1當(dāng)N值取得較大時(shí)，上去就是原周期看號(hào)f(t)的一不很好近似。式(3-7)常稱做f(t) 的截?cái)喔盗⑷~級(jí)數(shù)表示。MATLAB的符號(hào)積分函數(shù)int()可以用來(lái)求解連續(xù)時(shí)間周期信號(hào)的截?cái)喔盗⑷~級(jí)數(shù) 及傅立葉表示。求積函數(shù)int()的具體使用格式如下：intf=int(f,v);給出符號(hào)表達(dá)式f對(duì)指定變量v的(不帶積分常數(shù))不定積分；intf=int(f,v,a,b);給出符號(hào)表達(dá)式f對(duì)指定變量v的定積分；3、利用MATLAB實(shí)現(xiàn)周期信號(hào)的傅立葉級(jí)數(shù)分解與綜合(1)利用MATLAB求解周期矩形脈沖傅立葉級(jí)數(shù)，并繪制出各次諧波疊加的傅立葉綜 合波形圖。周期矩

5、形脈沖為f(t)=京6 -nT),式中t=1,T=5。 n =二二采用三角形式傅立葉級(jí)數(shù)分解與綜合形式，用式(3-2)(3-4)求出傅立葉級(jí)數(shù)分 解系數(shù)，運(yùn)用MATLAB的符號(hào)運(yùn)算功能，用式(3-7)實(shí)現(xiàn)信號(hào)的綜合，諧波的階數(shù)Nf =6。(a)實(shí)現(xiàn)流程利用MATLAB實(shí)現(xiàn)上述分析過(guò)程的流程如下：編寫(xiě)子函數(shù)x=time_fun_x(t),用符號(hào)表達(dá)式表示出周期信號(hào)在第一個(gè)周期內(nèi)的 符號(hào)表達(dá)式，并賦值返回給符號(hào)變量 x;編寫(xiě)子函數(shù)y=time_fun_e(t),求出該周期信號(hào)在繪圖區(qū)間內(nèi)的信號(hào)樣值，并賦 值給返回變量y;編寫(xiě)求解信號(hào)傅立葉系數(shù)及繪制合成波形圖的通用CTFShchsym.m,該函數(shù)流

6、程如下： 調(diào)用函數(shù)time_fun_x(t),獲取周期信號(hào)的符號(hào)表達(dá)式；求出信號(hào)的傅立葉索數(shù)；求出各次諧波；繪制各次諧波疊加波形圖； 調(diào)用函數(shù)time_fun_e(t),繪制原信號(hào)波形圖。(b) MATLAB算法提示及近明采用符號(hào)積分int求一個(gè)周期內(nèi)時(shí)間函數(shù)的三角級(jí)數(shù)展開(kāi)系數(shù)：& = Ao ,& =As, h =Bs,即計(jì)算式(3-2)(3-4)的值；用循環(huán)語(yǔ)句foL-end求出三角級(jí)數(shù)展開(kāi)系數(shù)a, bn的數(shù)值，分別為 A_sym,B_sym;用disp()語(yǔ)句輸出三角級(jí)數(shù)展開(kāi)系數(shù)A_sym, B_sym；用傅立葉三角級(jí)數(shù)展開(kāi)式(3-7)合成自續(xù)時(shí)間信號(hào)；化簡(jiǎn)表達(dá)式，據(jù)函數(shù)奇偶性可知，若 f

7、(t)為奇函數(shù)，則 a=0;若f(t)為偶函 數(shù)，M bn =0。(c)源程序編寫(xiě)函數(shù)文件CTFShchsym.m,這是一個(gè)計(jì)算連續(xù)時(shí)間周期信號(hào)的三角級(jí)數(shù)前6次展開(kāi)系數(shù)，再用這6次諧波合成原連續(xù)時(shí)間周期信號(hào)的程序，如下所示。function A_sym,B_sym=CTFShchsym%采用符號(hào)計(jì)算柬二個(gè)周期內(nèi)連續(xù)時(shí)間函數(shù) f的三角級(jí)數(shù)展開(kāi)系數(shù)，再用這些% 展開(kāi)系數(shù)合成連續(xù)時(shí)間函數(shù)f.傅立葉級(jí)數(shù)16% 函數(shù)的輸入輸出都是數(shù)值量% Nf=6諧波的階數(shù)% Nn輸出數(shù)據(jù)的準(zhǔn)確位數(shù)% A_sym第1元素是直流項(xiàng)，其后元素依次是 1,2,3次諧波cos項(xiàng)展開(kāi)系數(shù)% B_sym第2,3,4,.元素依次是1

8、,2,3次諧波sin項(xiàng)展開(kāi)系數(shù)% tao=1 tao/T=0.2 syms t n k x T=5;tao=0.2*T;a=0.5;if nargin4;Nf=6;endif nargin5;Nn=32;end x=time_fun_x(t);A0=int(x,t,-a,T-a)/T;%求出三角函數(shù)展開(kāi)系數(shù)A0As=int(2*x*cos(2*pi*n*t/T)/T,t,-a,T-a);%求出三角函數(shù)展開(kāi)系數(shù) AsBs=int(2*x*sin(2*pi*n*t/T)/T,t,-a,T-a);% 求出三角函數(shù)展開(kāi)系數(shù) BsA_sym(1)=double(vpa(A0,Nn);%獲取用數(shù)組 A0

9、所對(duì)應(yīng)的 ASC2桎值數(shù)組for k=1:NfA_sym(k+1)=double(vpa(subs(As,n,k),Nn);%獲取用數(shù)組 A 所對(duì)應(yīng)的 ASC2碼數(shù)他數(shù)組B_sym(k+1)=double(vpa(subs(Bs,n,k),Nn);%獲取用數(shù)組 B 所對(duì)應(yīng)的 ASC2碼數(shù)R數(shù)組 end if nargout=0c=A_sym;disp(c)%輸出c為三角級(jí)數(shù)展開(kāi)系數(shù)：第1元素是直流項(xiàng)，其后元素依次是1,2,a.次諧波cos項(xiàng)展開(kāi)系數(shù)d=B_sym;disp(d) %輸出d為三角級(jí)數(shù)展開(kāi)系數(shù)：第2,3,4,.元素依次是1,2,3次諧波sin項(xiàng)展開(kāi)系數(shù)t=-8*a:0.01:T-a

10、;f1=c(1)+c(2).*cos(2*pi*1*t/5)+0.*sin(2*pi*1*t/5); ;% 基波f2=c(3).*cos(2*pi*2*t/5)+0.*sin(2*pi*2*t/5); ;% 2 次t皆波f3=c(4).*cos(2*pi*3*t/5)+0.*sin(2*pi*3*t/5);% 3 次 t皆波f4=c(5).*cos(2*pi*4*t/5)+0.*sin(2*pi*4*t/5); ;% 4 次t皆波f5=c(6).*cos(2*pi*6*t/5)+0.*sin(2*pi*6*t/5);% 6 次 t皆波f6=f1+f2;% 基波+2次諧波f7=f6+f3;%基波

11、+2次諧波+3次諧波f8=f7+f4+f5;%基波+2次諧波+3次諧波+4次諧波+6次諧波subplot(2,2,1)plot(t,f1),hold ony=time_fun_e(t)%調(diào)用連續(xù)時(shí)間函數(shù)-周期矩形脈沖plot(t,y,r:)title(周期矩形波的形成一基波+直流)axis(-4,4.5,-0.5,1.5)subplot(2,2,2)plot(t,f6),hold on y=time_fun_e(t)17plot(t,y,r:)t田e(周期矩形波的形成一基波+2次諧波)axis(-4,4.5,-0.5,1,5)subplot(2,2,3)plot(t,f7),hold ony=

12、time_fun_e(t)plot(t,y,r:)title(基波+2次諧波+3次諧波)axis(-4,4.5,-0.5,1.5)subplot(2,2,4)plot(t,f8),hold ony=time_fun_e(t)plot(t,y,r:)title(基波+2次諧波+3次諧波+4次諧波+6次諧波)axis(-4,4.5,-0.5,1.5)end%function x=time_fun_x(t)%該函數(shù)是CTFShchsym.m的子函數(shù)。它由符號(hào)變量和表達(dá)式寫(xiě)成。h=1;x1=sym(Heaviside(t+0.5)*h;x=x1-sym(Heaviside(t-0.5)*h;%func

13、tion y=time_fun_e(t)%該函數(shù)是CTFShchsym.m的子函它由符號(hào)函數(shù)和表達(dá)式寫(xiě)成a=0.5;T=5;h=1;tao=0.2*T;t=-8*a:0.01:T-a;e1=1/2+1/2.*sign(t+tao/2);e2=1/2+1/2.*sign(t-tao/2);y=h.*(e1-e2);%連續(xù)時(shí)間函數(shù)-周期矩形脈沖(d)程序運(yùn)行結(jié)果在MATLAB命令窗口鍵入CTFShchsym,并按回車鍵，即可繪制出周期矩形波形 各次諧波合成波形，如圖3-1所示。18執(zhí)行如下命令就可以得到三角級(jí)數(shù)展開(kāi)系數(shù)：A_sym,B_sym=CTFShchsym結(jié)殺為一A_sym =- 0.20

14、000.37420.30270.20180.09350.0000-0.0624B_sym =0000000(2)實(shí)驗(yàn)內(nèi)容, f(t)n(a)(b)首先，推導(dǎo)出求解a。，an, b 利用MATLAB求解a , a ,要求：仿照例程，實(shí)現(xiàn)下述周期信號(hào)的傅立葉級(jí)數(shù)分解與合成：的公式，計(jì)算出前10次系數(shù)；bn的值，其中bn求解前10次系數(shù)，并給出利用這些系數(shù)合成的信號(hào)波形。三、周期信號(hào)頻譜分析1、周期信號(hào)頻譜對(duì)周期為Ti的信號(hào)f(t)進(jìn)行傅立葉級(jí)數(shù)展開(kāi)可得到f(t) =Fnen 二：二00cO二 a0 八 a” cos(nt)八 bn sin(nt)其中2 二 i = 一Ti_ 曲 1Fn = Fn

15、e 4=2 (an - jbn)(3-8)如果求出a和bn ,根據(jù)以下兩式可以畫(huà)出周期信號(hào)的幅度譜|Fn切和相位譜中n6。1,一一Fn =一；耳+幅頻2Q = -arctan2相頻an2、MATLAB 實(shí)現(xiàn)由上述分析可知，只要求出了周期信號(hào)的傅立葉級(jí)數(shù)幅度和相位Q ,就可以根據(jù)Fn和中n隨頻率8二門町的變化關(guān)系畫(huà)出幅度譜和相位譜?，F(xiàn)仍以周期矩形脈沖為f。)=京0騏-nT),(7=1,T=5)為例，來(lái)說(shuō)明如何用MATLAB繪制周期信號(hào)的頻譜圖，并對(duì)周期信號(hào)的頻譜特性進(jìn)行分析。19由于繪制頻譜圖的前提是必須先求出周期信號(hào)的傅立葉級(jí)數(shù)系數(shù)，因此只需對(duì)上述給出的求周期信號(hào)傅立葉級(jí)數(shù)的函數(shù) CTFShc

16、hsym.m進(jìn)行適當(dāng)修改，即可編寫(xiě)出繪制周 期信號(hào)頻譜的通用函數(shù)。需注意的是，由于周期信號(hào)的頻譜是離散的，故在繪制頻譜時(shí)，采用的是stem命令而不是plot命令，下面給出實(shí)現(xiàn)上述通用程序 CTFStpshsym.m的過(guò)程。(1)實(shí)現(xiàn)流程此處采用三角形式傅立葉級(jí)數(shù)分解，求出分解系數(shù)an和bn,再利用式(3-8)求出傅立葉復(fù)指數(shù)系數(shù)Fn,畫(huà)出Fn的幅度譜Fn 8。諧波的階數(shù)Nf可任意指定，此處指定Nf =60。(a)實(shí)現(xiàn)流程如下：編寫(xiě)子函數(shù)y=time_fun_s(t),用符號(hào)表達(dá)式表示出周期信號(hào)在第一個(gè)周期內(nèi)的符號(hào)表達(dá)式，并賦值返回給符號(hào)變量 V；編寫(xiě)子函數(shù)x=time_fun_e(t),求出該

17、周期信號(hào)在繪圖區(qū)間內(nèi)的信號(hào)樣值，并賦 值給返回變量x;編寫(xiě)求解信號(hào)傅立葉復(fù)指數(shù)系數(shù) Fn的通用函數(shù)，該函數(shù)流程如下： 調(diào)用函數(shù)time_fun_s(t),獲取周期信號(hào)的符號(hào)表達(dá)式；求出信號(hào)的三帚形氐的傅立葉系數(shù)an和bn;求出信號(hào)的傅立復(fù)指數(shù)系數(shù)Fn；繪制Fn的幅度頻譜圖； 調(diào)用函數(shù)time_fun_e(t)繪制信號(hào)波形。(b)算法提示采用符號(hào)積分int求0,T內(nèi)時(shí)間信號(hào)的三角級(jí)數(shù)展開(kāi)系數(shù)：4 = A0 , a = As,bn = Bs;用循環(huán)語(yǔ)句foL-end求出三角級(jí)數(shù)展開(kāi)系數(shù)an和bn的數(shù)值，并賦值給變量A_sym(k+1), B_sym(k+1)；心角級(jí)數(shù)展開(kāi)家數(shù)a和bn得到復(fù)指數(shù)展開(kāi)

18、系數(shù)Fn： Fn可以根據(jù)式(38)求出, 需要注意的是an、bn和Fn的自變量取值情況， 耳、R的變量h的取值范圍為 n=0,1,2,3，N，而匕的變量n的取值范圍為n = 0,1,2,3,，土N ,為了從出和 bn得到E,需要用到MATLAB的反折函數(shù)fliplr來(lái)實(shí)現(xiàn)頻譜的反折；例如已知0, 1, 2, 3, 4。要組成員序列反折后再-4, -3, -2, -1, 0, 1, 2, 3, 4)。程序一個(gè)序列a=5 4 3 2 1,采樣位置為: 疊加原序列的一個(gè)新序列(采樣位置： 如下：n=0:4;a=5 4 3 2 1;subplot(2,1,1);stem(n,a);b=fliplr(a

19、);k=-4:4;c=b,a(2:end);subplot(2,1,2);stem(k,c);結(jié)果如右圖所示。20(c)源程序代碼編寫(xiě)函數(shù)文件CTFStpshsym.m,如下所示：CTFStpshsym.mfunction A_sym,B_sym=CTFSshbpsym(T,Nf)%采用符彳計(jì)算親0,T內(nèi)時(shí)間函數(shù)的三角級(jí)數(shù)展開(kāi)系數(shù)。%函數(shù)的輸入輸出都是數(shù)值量% Nn輸出數(shù)據(jù)的準(zhǔn)確位數(shù)% A_sym第1元素是直流項(xiàng)，其后元素依次是1,2,3次諧波cos項(xiàng)展開(kāi)系數(shù)% B_sym第2,3,4,.元素依次是1,2,3.次諧波sin項(xiàng)展開(kāi)系數(shù)% T T=m*tao, 信號(hào)周期%Nf諧波的階數(shù)%Nn輸出數(shù)

20、據(jù)的準(zhǔn)確位數(shù)% m (m=T/tao)周期與脈沖寬度之比，如m=4,8,16,100等%tao脈寬:tao=T/msyms t n yNf=input(please Input所需展開(kāi)的最高諧波次數(shù)：Nf=);T=input(please Input 信號(hào)的周期 T=);Nn=32;y=time_fun_s(t);A0=2*int(y,t,0,T)/T;As=int(2*y*cos(2*pi*n*t/T)/T,t,0,T);Bs=int(2*y*sin(2*pi*n*t/T)/T,t,0,T);A_sym(1)=double(vpa(A0,Nn); for k=1:NfA_sym(k+1)=d

21、ouble(vpa(subs(As,n,k),Nn);B_sym(k+1)=double(vpa(subs(Bs,n,k),Nn); endS1=fliplr(A_sym);%對(duì)A_sym陣左右對(duì)稱交換S1(1,k+1)=A_sym(1); %A_sym 的 1*k 陣擴(kuò)展為 1*(k+1)陣S2=fliplr(1/2*S1);%對(duì)擴(kuò)展后的S1陣左右對(duì)稱交換回原位置S3=fliplr(1/2*B_sym); %對(duì) B_sym 陣左右對(duì)稱交換S3(1,k+1)=0;%B_sym 的 1*k 陣擴(kuò)展為 1*(k+1)陣S4=fliplr(S3);%3擴(kuò)展后的S3陣左右對(duì)稱交換回原位置S5=S2-i

22、*S4;%用三角函數(shù)展開(kāi)系數(shù) A、B值合成付里葉指數(shù)系數(shù)S6=fliplr(S5);N=Nf*2*pi/T;k2=-N:2*pi/T:N;S7=S6,S5(2:end);subplot(3,3,3)x=time_fun_e1(t) %調(diào)用連續(xù)時(shí)間函數(shù)-周期矩形脈沖stem(k2,abs(S7); %畫(huà)出周期矩形脈沖的頻譜(T=M*tao)title(連續(xù)時(shí)間函數(shù)周期矩形脈沖的雙邊幅度譜)axis(-80,80,0,0.12)line(-80,80,0,0)line(0,0,0,0.12)21%function y=time_fun_s(t,T)%該函數(shù)是CTFSshbpsym.m的子函數(shù)。它由

23、符號(hào)變量和表達(dá)式寫(xiě)成。syms a a1T=input(please Input 信號(hào)的周期 T=);A=1;tao=T/M;a=tao/2;y1=sym(Heaviside(t+a1)*A;y=y1-sym(Heaviside(t-a1)*A;y=subs(y,a1,a);y=simple(y);%function x=time_fun_e1(t)%該函數(shù)是CTFSshbpsym.m的子函數(shù)。它由符號(hào)變量和表達(dá)式寫(xiě)成。% t是時(shí)間數(shù)組% T是周期 duty=tao/T=0.2T=5;t=-2*T:0.01:2*T;tao=T/5;x=rectpuls(t,1); %產(chǎn)生一個(gè)寬度tao=1的矩

24、形脈沖plot(t,x)hold onx=rectpuls(t-5,1); %產(chǎn)生一個(gè)寬度tao=1的矩形脈，中心位置在t=5處 plot(t,x) hold onx=rectpuls(t+5,1); %產(chǎn)生一個(gè)寬度tao=1的矩形脈，中心位置在t=-5處 plot(t,x)title(周期為T=5,脈寬tao=1的矩形脈沖) axis(-10,10,0,1.2)(d)運(yùn)行結(jié)果及分析程序運(yùn)行結(jié)果調(diào)用CTFStpshsym.m函數(shù)文件，即可繪出周期矩形脈沖波形信號(hào)的雙邊幅度頻譜。指令如下：在 MATLAB 命令窗口鍵入 CTFStpshsym,回車，命令窗口將會(huì)出現(xiàn)：please Input所需展開(kāi)的最高