現(xiàn)代測試信號分析與處理

1、現(xiàn)代測試信號分析與處理 現(xiàn)代測試信號分析與處理課程作業(yè)學生姓名： 李 亮 學 號： s201301109 導 師： 焦敬品 教授 2013年12 月 05日目 錄作業(yè)11第1題1第2題3第3題4第4題6作業(yè)二8第1題8第2題12第3題13第4題16作業(yè)320第1題20第2題23第3題24第4題27第5題32作業(yè)1第1題 用matlab編程，繪制出周期矩形脈沖信號的幅頻譜圖，其中周期矩形脈沖的周期為 t=1ms，脈沖寬度為=0.1ms，如下圖所示。 解：生成周期矩形脈沖及其頻譜繪制的matlab程序如下：35%作業(yè)1_1.mclcclear allclose all%產(chǎn)生周期矩形波syms t

2、n t tau_t;tau_t=0.1;t=1;a=1;t=-3.5:1/1e4:3.5;d=-4:4;f=pulstran(t,d,rectpuls,0.1); figure(1);plot(t,f);box off;axis(-4,4,0,1.2);xlabel(時間t/ms,fontsize,14);ylabel(幅值x(t),fontsize,14);title(周期矩形脈沖時域波形,fontsize,14);set(gca,xticklabel, ,-3t,-2t,-t,0,t,2t,3t, ,fontsize,12);%周期矩形脈沖的傅里葉級數(shù)三角函數(shù)展開式a0=1/t*int(1

3、,-tau_t/2,tau_t/2); %計算直流分量a_0=subs(a0);an=2/t*int(cos(2*pi*n*t),-tau_t/2,tau_t/2);%計算余弦分量；bn=2/t*int(sin(2*pi*n*t),-tau_t/2,tau_t/2);%計算正弦分量n0=0;n1=40;n=0:n1;a_n=subs(an);b_n=subs(bn);a_n=sqrt(a_n.2+b_n.2);p_n=atan2(-b_n,a_n)*180/pi;n=0,n;a_n=a_0,a_n;p_n=0,p_n;figure(2);stem(n,a_n,fill);box off;yla

4、bel(幅值an,fontsize,14);xlabel(nomega,fontsize,14);title(周期矩形脈沖三角函數(shù)展開式幅值譜,fontsize,14);set(gca,xticklabel, ,2/,4/,6/,8/,fontsize,14)figure(3);stem(n,p_n,fill);box off;ylabel(phi_n(度),fontsize,14);xlabel(n,fontsize,14);title(周期矩形脈沖相位譜,fontsize,14);%周期矩形脈沖的傅里葉級數(shù)雙邊復指數(shù)展開式clear all;syms t n;t=1;tau_t=0.1;a

5、=1;x=a*exp(-j*n*2*pi/t*t);xn=int(x,t,-tau_t/2,tau_t/2);xn=simple(xn);n=-40:1,realmin,1:40;xn=subs(xn,n,n);figure(4);stem(n,xn,fill);box off;ylabel(幅值xn,fontsize,14);xlabel(nomega,fontsize,14);title(傅里葉級數(shù)的雙邊復指數(shù)展開式,fontsize,14)set(gca,xticklabel,-8/,-6/,-4/,-2/,0,2/,4/,6/,8/,fontsize,12)運行上述“作業(yè)1_1.m”程

6、序，得到周期矩形脈沖時域波形如圖1.1所示，得到周期矩形脈沖傅里葉級數(shù)三角函數(shù)展開式幅頻譜如圖1.2所示，得到周期矩形信號傅里葉級數(shù)復指數(shù)展開式幅頻譜如圖1.3所示： 圖1.1 周期矩形脈沖時域波形圖圖1.2 周期矩形脈沖三角函數(shù)展開式幅值譜圖1.3 周期矩形脈沖的傅里葉級數(shù)的雙邊復指數(shù)展開式分析：從周期矩形脈沖的傅里葉級數(shù)的三角函數(shù)展開式和雙邊復指數(shù)展開式可以看出，周期信號的頻譜具有以下3個特點:1) 周期信號的頻譜是離散譜；2) 周期信號的譜線僅出現(xiàn)在基波及各次諧波頻率處；3) 周期信號的幅值譜中各頻率分量的幅值隨著頻率的升高而減小，頻率越高，幅值越小。此外，比較傅里葉級數(shù)的兩種展開式的頻

7、譜圖(圖1.2和圖1.3)可知，由三角函數(shù)表達的傅里葉級數(shù)的頻譜為單邊譜，而以復指數(shù)函數(shù)表達的傅里葉級數(shù)的頻譜為雙邊譜。兩種形式的幅值譜在幅值上的關系是xn=an/2，即雙邊譜中各諧波的幅值為單邊譜中各對應諧波幅值的一半。第2題 一周期為2的方波信號f(t)為：ft=1, &0t-1, &t2用matlab語言編程實現(xiàn)：1）畫出方波圖形；2）分別用其基波，基波和3次諧波，基波、3次諧波、5次諧波、7次諧波、9次諧波合成波形，并與原圖形相比較，觀察其逼近程度。解:生成周期方波時域波形及其不同階數(shù)的諧波疊加的matlab程序如下：%作業(yè)1_2.mclcclear allclose all%繪制周期

8、方波的時域波形圖n_max=1 3 9;n=length(n_max);t=linspace(-10,10,1000);f=square(t);subplot(2,2,1);plot(t,f,r,linewidth,2);title(方波時域波形圖);xlabel(時間/s);ylabel(幅值);%繪制不同頻次諧波疊加后的時域波形圖w0=1;for k=1:n n=; n=1:2:n_max(k); b_n=4./(pi*n); x=b_n*sin(w0*n*t); subplot(2,2,k+1),plot(t,f,r:,t,x,linewidth,2); title(最高諧波次數(shù)=,nu

9、m2str(n_max(k);xlabel(時間/s);ylabel(幅值);end運行上述“作業(yè)1_2.m”程序，得到不同諧波次數(shù)疊加后的時域波形如圖2所示：圖2 周期方波時域波形圖及不同最高諧波次數(shù)疊加合成的波形圖分析：從圖2可以看出，周期信號本身可以用傅里葉級數(shù)的三角函數(shù)展開式展開成無窮多項諧波信號的疊加，所以周期信號可以用傅里葉級數(shù)中的諧波的和來逼近，而且，所取的項數(shù)越多，亦即n越大，近似的精度越高，合成后的信號越接近原來的周期信號。第3題 帶限信號 f(t)ft=sat=sintt其最高頻率m=1, 用matlab繪制三種不同采樣頻率s的抽樣信號。1) s=102) s=0.33)

10、s=4生成帶限信號時域波形及不同采樣頻率下所得到的的信號波形的matlab程序如下所示：%作業(yè)1_3.mclcclear allclose all%帶限信號f(t)時域波形繪制t1=-4*pi,t2=4*pi,t=1/10,10/3,1/4;fs=10,0.3,4;%定義采樣頻率n=length(t);t=linspace(t1,t2,2000);ft=sin(t)./t;subplot(2,2,1)plot(t,ft);title(帶限信號時域波形圖,fontsize,12);box off;xlabel(時間/s,fontsize,12);ylabel(幅值,fontsize,12);ax

11、is(-4*pi,4*pi,-0.4,1.2);%繪制不同采樣頻率下的帶限信號的時域波形for k=1:n t=t1:t(k):t2 subplot(2,2,k+1); ft=sin(t)./t; stem(t,ft); axis(-4*pi,4*pi,-0.4,1.2);box off; title(采樣頻率為,num2str(fs(k); xlabel(時間/s,fontsize,12);ylabel(幅值,fontsize,12)end運行上述“作業(yè)1_3.m”程序，帶限信號時域波形及不同采樣頻率下所得信號的波形如圖3所示：圖3 帶限信號時域波形及不同采樣頻率下所得信號的波形分析：從圖3

12、中可以看出，采樣頻率越大，采樣間隔越小，所采集得到的數(shù)據(jù)點就越多，這些采集的數(shù)據(jù)點的包絡線越接近原連續(xù)信號的時域波形；采樣頻率越小，采樣間隔越大，采集的數(shù)據(jù)點越少，這些數(shù)據(jù)點的包絡線偏離原連續(xù)信號的波形越大。所以，在實際操作中，為了減少失真，應選擇合適的采樣頻率。第4題 矩形脈沖信號 f(t) ft=1, &-0.5t0.50, &(其他)用matlab語言編程實現(xiàn)：1) 繪制該信號的幅頻譜；2) 通過對該信號頻譜進行反傅里葉變換，重新獲得該信號波形。解:生成矩形脈沖信號的時域波形及矩形脈沖信號傅里葉變換后的幅頻譜的matlab程序如下：%作業(yè)1_4.mclcclear allclose al

13、l%生成矩形脈沖函數(shù)f=sym(heaviside(t+0.5)-heaviside(t-0.5);%用階躍信號生成矩形脈沖信號subplot(221);ezplot(f,-1,1);xlabel(時間t/s);ylabel(幅值f(t);box off;title(矩形脈沖信號時域波形圖);axis(-1,1,0,1.2);%矩形脈沖函數(shù)傅里葉變換f=fourier(f);%傅里葉變換subplot(222);ezplot(abs(f),-15*pi,15*pi);axis(-15*pi,15*pi,0,1.2);f1=simple(f);ylabel(幅值f(omega);xlabel(頻

14、率omega);title(矩形脈沖信號幅頻譜圖);box off;subplot(223);%矩形脈沖函數(shù)反傅里葉變換f2=ifourier(f);ezplot(f2,-1,1);axis(-1,1,0,1.2);xlabel(時間t/s);ylabel(幅值f2);box off;title( 反傅里葉變換后時域波形圖);運行程序“作業(yè)1_4.m”,得到矩形脈沖信號的時域圖、頻譜圖和反傅里葉變換后時域波形如圖4所示:圖4 矩形脈沖信號的時域圖、頻譜圖和反傅里葉變換后時域波形分析：從圖4中可以看出，矩形脈沖信號的頻譜是sinc函數(shù)，信號的能量主要集中在第一個主值區(qū)間內(nèi)；此外，矩形脈沖信號在時

15、域上是有限的，對應的傅里葉變換的頻譜在頻域上則是無限的。從矩形脈沖信號的傅里葉變換和反傅里葉變換圖形可以看出，矩形脈沖函數(shù)與sinc函數(shù)是一個傅里葉變換對，二者可以相互轉(zhuǎn)換，而轉(zhuǎn)換的紐帶是傅里葉變換。作業(yè)二第1題 對連續(xù)諧波信號ft=cossf0t，f0=20hz進行截斷，在采樣間隔為0.005s進行采樣處理后，進行dft分析。1）若截斷長度為0.6s，在窗函數(shù)分別選擇矩形窗函數(shù)、高斯窗函數(shù)、漢寧窗函數(shù)、三角窗函數(shù)時，繪制其截斷后信號的時域波形及頻譜圖。2）若截斷長度變?yōu)?.61s，其波形和頻譜有何變化，請進行對比分析。解：繪制原始信號f（t）時域波形圖、頻譜圖及用不同長度、不同類型的窗函數(shù)對

16、原始信號進行截斷后的時域波形圖和頻譜圖的matlab程序如下：%作業(yè)2_1.mclcclear allclose all%截斷前連續(xù)諧波信號ft時域波形f0=20,ts=0.005,fs=1/ts;%諧波信號頻率、采樣間隔、采樣頻率n1=0.6/0.005;n2=0.61/0.005;%采樣點數(shù)t=linspace(-0.5,0.5,200);f=cos(2*pi*f0*t);figure(1);subplot(121);plot(t,f);box off;title(連續(xù)諧波信號ft時域波形);xlabel(時間/s);ylabel(幅值);n1=0:n1-1;fr1=n1*fs/n1;%0

17、.6s的頻率轉(zhuǎn)換n2=0:n2-1;fr2=n2*fs/n2;%0.61s的頻率轉(zhuǎn)換%截斷前連續(xù)諧波信號的fftft=fft(f,n1);subplot(122);plot(fr1,abs(ft)*2/n1);box off;title(連續(xù)諧波信號fft);xlabel(頻率/hz);ylabel(幅值f(jomega);axis(0,100,0,1);set(gca,xtickmode,manual,xtick,0,20,50,100);%窗函數(shù)寬度為0.6t1=-0.3:1/fs:0.3;%定義窗函數(shù)寬度為0.6fm=cos(2*pi*20*t1);%定義不同的窗函數(shù)h1=window(

18、rectwin,length(t1);%矩形窗h2=window(gausswin,length(t1); %高斯窗h3=window(hanning,length(t1); %漢寧窗h4=window(triang,length(t1);%三角窗%對連續(xù)余弦諧波信號進行截斷f1=fm.*h1;f2=fm.*h2;f3=fm.*h3;f4=fm.*h4;%對截斷后的余弦諧波信號進行fftff1=fft(f1,n1);ff2=fft(f2,n1);ff3=fft(f3,n1);ff4=fft(f4,n1);%截斷長度為t2=0.61t2=-0.305:1/fs:0.305;fn=cos(2*pi

19、*20*t2);c1=window(rectwin,length(t2);c2=window(gausswin,length(t2); c3=window(hamming,length(t2); c4=window(triang,length(t2);y1=fn.*c1;y2=fn.*c2;y3=fn.*c3;y4=fn.*c4;fy1=fft(y1,n2);fy2=fft(y2,n2);fy3=fft(y3,n2);fy4=fft(y4,n2);%截斷長度為0.6s和0.61s的時域波形對比figure(2);subplot(221),plot(t1,f1,b);box off;hold o

20、nplot(t2,y1,r);box off;axis(-0.6,0.6,-1,1);title(矩形窗0.6s截斷與0.61s截斷比較);xlabel(時間/s);ylabel(幅值)legend(0.6s,0.61s)subplot(222),plot(t1,f2,b);box off;hold onplot(t2,y2,r);box off;axis(-0.6,0.6,-1,1);title(高斯窗0.6s截斷與0.61s截斷比較);xlabel(時間/s);ylabel(幅值)legend(0.6s,0.61s)subplot(223),plot(t1,f3,b);box off;ho

21、ld onplot(t2,y3,r);box off;axis(-0.6,0.6,-1,1);title(漢寧窗0.6s截斷與0.61s截斷比較);xlabel(時間/s);ylabel(幅值)legend(0.6s,0.61s)subplot(224),plot(t1,f4,b);box off;hold onplot(t2,y4,r);box off;axis(-0.6,0.6,-1,1);title(三角窗0.6s截斷與0.61s截斷比較);xlabel(時間/s);ylabel(幅值)legend(0.6s,0.61s)%頻譜對比figure(3)subplot(221),plot(f

22、r1,abs(ff1)*2/n1,b);box off;hold onplot(fr2,abs(fy1)*2/n2,r);box off;title(矩形窗0.6s截斷與0.61s截斷比較);set(gca,xtickmode,manual,xtick,0,20,50,100);axis(0,100,0,1)xlabel(頻率/hz);ylabel(幅值f(jomega)legend(0.6s,0.61s)subplot(222),plot(fr1,abs(ff2)*2/n1,b);box off;hold onplot(fr2,abs(fy2)*2/n2,r);box off;title(高

23、斯窗0.6s截斷與0.61s截斷比較);set(gca,xtickmode,manual,xtick,0,20,50,100);axis(0,100,0,1)xlabel(頻率/hz);ylabel(幅值f(jomega)legend(0.6s,0.61s)subplot(223),plot(fr1,abs(ff3)*2/n1,b);box off;hold onplot(fr2,abs(fy3)*2/n2,r);box off;set(gca,xtickmode,manual,xtick,0,20,50,100);title(漢寧窗0.6s截斷與0.61s截斷比較);axis(0,100,0

24、,1)xlabel(頻率/hz);ylabel(幅值f(jomega)legend(0.6s,0.61s)subplot(224),plot(fr1,abs(ff4)*2/n1,b);box off;hold onplot(fr2,abs(fy4)*2/n2,r);box off;title(三角窗0.6s截斷與0.61s截斷比較);set(gca,xtickmode,manual,xtick,0,20,50,100);axis(0,100,0,1)xlabel(頻率/hz);ylabel(幅值f(jomega)legend(0.6s,0.61s)運行上述matlab程序，得到原始信號的時域和

25、頻域波形圖如圖1.1所示，得到不同長度、不同類型窗函數(shù)截斷后的信號的頻譜圖如圖1.2所示。圖1.1 原始信號時域波形圖和頻譜圖圖1.2 時間長度0.6s和0.61s不同窗函數(shù)對原始信號截斷后時域波形圖圖1.3 時間長度0.6s和0.61s不同窗函數(shù)對原始信號截斷后頻譜圖分析：原始信號ft=cos2f0t，在時域上為無限長的周期信號，而計算機只能處理有限長的信號，所以在對原信號進行處理前需要對其進行截斷處理。所謂截斷，就是將無限長的信號與一定寬度的窗函數(shù)相乘。時域上的截斷會在頻域上導致泄露，從圖1.3中可以看出，截斷后的信號頻譜在f=20hz附近出現(xiàn)了泄露。為減小泄露，我們可以采用整周期采樣的方

26、法，從圖1.3中可以看出，當采用截斷長度為0.6s，即整周期截斷時泄露較小。不同類型的窗函數(shù)的主瓣和旁瓣寬度有所差異，主瓣越窄，旁瓣越小，則截斷后的能量泄露就越小。從圖1.3中可以看出，在矩形窗、高斯窗、漢寧窗和三角窗中，矩形窗的主瓣最窄，所以用矩形窗截斷后的信號的頻譜在f=20hz時最集中；但矩形窗旁瓣較高，并且還有負旁瓣。相比之下，漢寧窗和三角窗的旁瓣較小，但主瓣寬，所以由旁瓣引起的能量泄露減小，由主瓣引起的能量泄露大。綜上，在對信號進行截斷處理時，為減小能量泄露，可以采用整周期采樣，并根據(jù)不同的信號，選取合適的窗函數(shù)。第2題 對連續(xù)的單頻周期信號xat=cos(2fat)，fa=60hz

27、，按 fs=8fa采樣，截取長度分別選n=50、 n=30、 n=16和n=4，觀察其幅頻譜。解:對連續(xù)單頻周期信號xat=cos(2fat)，fa=60hz，進行采樣和不同長度截斷的頻譜分析的matlab程序如下：%作業(yè)2_2.mclcclear allclose all%基本參數(shù)設置fa=60,fs=8*fa,n=50,30,16,4;%采樣頻率，截取長度l=length(n);t=linspace(-0.5,0.5,200);%編寫dft變換函數(shù)for k1=1:l n=0:(n(k1)-1); xn=cos(2*pi*fa*n/fs);%原連續(xù)信號序列化 k=0:(n(k1)-1);

28、wn=exp(-j*2*pi/n(k1); nk=n*k; wnnk=wn.nk; xk=xn*wnnk; subplot(2,2,k1); f=(0:n(k1)-1)*fs/n(k1);%頻率轉(zhuǎn)換 plot(f,abs(xk)*2/n(k1),linewidth,2);box off; set(gca,xtickmode,manual,xtick,0,60,100,150,200); xlabel(頻率/hz,fontsize,12);ylabel(幅值,fontsize,12); title(截取長度n=,num2str(n(k1),fontsize,12); axis(0,200,0,1

29、)end運行上述matlab程序，得到不同截取長度的幅頻譜如圖2所示。圖2 截取長度分別為：n=50,30,16,4時信號的幅頻譜圖分析：從圖2中可以看出，截取長度n越大，得到的幅頻譜越接近信號的真實幅頻譜圖，當截取長度n取的很小時，信號的幅頻譜發(fā)生了嚴重失真，從n=4時的幅頻譜圖可以看出，此時已不能真實反映原信號的頻率成分，但是，如果選取很大的n值，則會引起信號處理時間的大大延長。綜上所述，在利用計算機對連續(xù)信號進行離散化和dft時，要充分考慮截取長度對信號頻譜的影響，既要能反映信號真實的頻率成分，又要提高信號處理的速度。第3題 一個連續(xù)信號含兩個頻率分量，經(jīng)采樣得xn=sin20.125n

30、+cos20.125+fn，n=0,1,2,n-1.已知n=16,f分別為1/16和1/64，觀察其頻譜；當n=128時，f不變，其結(jié)果如何。解：對經(jīng)過采樣后的信號xn=sin20.125n+cos20.125+fn分別在n=16，n=128，f=116和f=164時的頻譜分析的matlab程序如下：%作業(yè)2_3.mclcclear allclose all%基本參數(shù)設置n1=16,n2=128;%采樣點數(shù)fs=8;df1=1/16,df2=1/64;n1=0:(n1-1);n2=0:(n2-1);t1=n1/fs;t2=n2/fs;xt1=sin(2*pi*t1)+cos(2*pi*(1+8

31、*df1)*t1);%df1=1/16xt2=sin(2*pi*t2)+cos(2*pi*(1+8*df1)*t2);xt3=sin(2*pi*t1)+cos(2*pi*(1+8*df2)*t1);xt4=sin(2*pi*t2)+cos(2*pi*(1+8*df2)*t2);%df1=1/64%n=16,df=1/16,df=1/16時的時域波形和頻域波形繪制f1=n1*fs/n1;%頻率轉(zhuǎn)換f2=n2*fs/n2;fft1=fft(xt1,n1);figure(1);subplot(221);plot(t1,xt1);title(n=16,df=1/16時的時域波形圖);box off;x

32、label(時間t);ylabel(幅值);subplot(223);plot(f1,abs(fft1)*2/n1);title(n=16,df為1/16的幅頻譜圖);axis(0,4,0,1.2);box off;xlabel(頻率/hz);ylabel(幅值);fft3=fft(xt3,n1);subplot(222);plot(t1,xt3);title(n=16,df=1/64時的時域波形圖);box off;xlabel(時間t);ylabel(幅值);subplot(224);plot(f1,abs(fft3)*2/n1);title(n=16,df為1/64的幅頻譜圖);axis

33、(0,4,0,1.2);%n=128,df=1/16，df=1/64時的時域波形和頻域波形繪制figure(2);subplot(221);plot(t2,xt2);title(n=128,df=1/16時的時域波形圖);box off;xlabel(時間t);ylabel(幅值);subplot(223);fft1=fft(xt2,n2);plot(f2,abs(fft1)*2/n2);title(n=128,df為1/16的幅頻譜圖);box off;xlabel(頻率/hz);ylabel(幅值);axis(0,4,0,1.2);subplot(222);plot(t2,xt4);tit

34、le(n=128,df=1/64時的時域波形圖);box off;xlabel(時間t);ylabel(幅值);subplot(224);fft4=fft(xt4,n2);plot(f2,abs(fft4)*2/n2);title(n=128,df為1/64的幅頻譜圖);box off;xlabel(頻率/hz);ylabel(幅值);axis(0,4,0,1.2);運行上述matlab程序，得到n=16，f=116和f=164時的時域波形圖和幅頻譜圖如圖3.1所示，得到n=128，f=116和f=164時的時域波形圖和幅頻譜圖如圖3.2所示圖3.1 n=16，f=116和f=164時的時域波

35、形圖和幅頻譜圖圖3.2 n=128，f=116和f=164時的時域波形圖和幅頻譜圖分析：從圖3.1和圖3.2對比可知，當采樣點數(shù)過少時，信號的時域和頻域波形都會發(fā)生嚴重的失真，不能正確的反映信號的真實信息；當采樣點數(shù)選擇恰當時，可以很好的識別信號的波形和信號中的頻率成分。從圖3.2可知，當n=128時，可以很好的分辨出f=116，而當f=164時，此時已不能清晰的分辨出f=164，這是由于頻率分辨力f=fsn與采樣頻率和采樣點數(shù)有關，當采樣頻率一定時，采樣點數(shù)越多，頻率分辨力越高，越能識別出 最小頻率成分，當采樣點數(shù)n和采樣頻率fs都確定時，最小頻率分辨力也就確定了，所以此時，n=128時可以

36、分辨出f=116，而f=164時，兩條譜線就混在了一起，不能很好的分別。所以，在對信號處理的時候，我們可以通過增加采樣點數(shù)來增加頻率分辨力，這對識別頻譜中的細微差別有很大幫助。第4題 一個連續(xù)信號含兩個頻率分量，其表達式為xt=cos22t+sin(24.5t)在采樣間隔為0.05條件下，進行dft分析：1）矩形窗截斷，長度n=40；2）漢寧窗截斷，長度n=40；3）三角窗截斷，長度n=40；4）矩形窗截斷，長度n=100；5）漢寧窗截斷，長度n=100；6）三角窗截斷，長度n=100；解：對連續(xù)信號xt=cos22t+sin(24.5t)在采樣間隔為0.05條件下，用不同長度、不同類型的窗函

37、數(shù)進行dft分析的matlab程序如下：%作業(yè)2_4.mclcclear allclose all%基本參數(shù)設置ts=0.05;fs=1/ts;%采樣間隔和采樣頻率n1=40,n2=100;%窗長度t1=(0:n1-1)*ts;t2=(0:n2-1)*ts;f1=(0:n1-1)*fs/n1;f2=(0:n2-1)*fs/n2;y1=cos(2*pi*2*t1)+sin(2*pi*4.5*t1);y2=cos(2*pi*2*t2)+sin(2*pi*4.5*t2);%原始信號時域波形繪制figure(1);plot(t2,y2);xlabel(時間/s);ylabel(幅值);axis(0,2

38、,-2,2);title(時域波形圖);%定義長度為n1=40的不同的窗函數(shù)h1=window(rectwin,length(t1);%矩形窗h2=window(hanning,length(t1); %漢寧窗h3=window(triang,length(t1);%三角窗%定義長度為n2=100的不同的窗函數(shù)h4=window(rectwin,length(t2);%矩形窗h5=window(hanning,length(t2); %漢寧窗h6=window(triang,length(t2);%三角窗%用長度為n1=40的窗截斷g1=y1.*h1;g2=y1.*h2;g3=y1.*h3;%

39、對截斷長度n=40的信號進行fft變換fg1=fft(g1,n1);fg2=fft(g2,n1);fg3=fft(g3,n1);%用長度為n2=100的窗截斷g4=y2.*h4;g5=y2.*h5;g6=y2.*h6;%對截斷長度n2=100的信號進行fft變換fg4=fft(g4,n2);fg5=fft(g5,n2);fg6=fft(g6,n2);figure(2);subplot(231);plot(f1,abs(fg1)*2/n1);title(n1=40矩形窗截斷);xlabel(頻率/hz);ylabel(幅值);axis(0,10,0,1);box off;set(gca,xtic

40、kmode,manual,xtick,0,2,4.5,10);subplot(234);plot(f2,abs(fg4)*2/n2);title(n2=100、矩形窗截斷);set(gca,xtickmode,manual,xtick,0,2,4.5,10);xlabel(頻率/hz);ylabel(幅值);axis(0,10,0,1);box off;subplot(232);plot(f1,abs(fg2)*2/n1);title(n1=40漢寧窗截斷);set(gca,xtickmode,manual,xtick,0,2,4.5,10);xlabel(頻率/hz);ylabel(幅值);

41、axis(0,10,0,1);box off;subplot(235);plot(f2,abs(fg5)*2/n2);title(n2=100漢寧窗截斷);set(gca,xtickmode,manual,xtick,0,2,4.5,10);xlabel(頻率/hz);ylabel(幅值);axis(0,10,0,1);box off;subplot(233);plot(f1,abs(fg3)*2/n1);title(n1=40三角窗截斷);set(gca,xtickmode,manual,xtick,0,2,4.5,10);xlabel(頻率/hz);ylabel(幅值);axis(0,10

42、,0,1);box off;subplot(236);xlabel(頻率/hz);ylabel(幅值);axis(0,10,0,1);box off;plot(f2,abs(fg6)*2/n2);title(n2=100三角窗截斷);set(gca,xtickmode,manual,xtick,0,2,4.5,10);xlabel(頻率/hz);ylabel(幅值);axis(0,10,0,1);box off;運行上述matalb程序，得到原信號時域波形圖如圖4.1所示，得到不同長度、不同類型窗函數(shù)截斷后信號的頻譜圖如圖4.2所示。圖4.1 原信號時域波形圖圖4.2 不同長度、不同類型窗函數(shù)

43、截斷后信號的頻譜圖分析：由分辨力f=fsn可知，截取長度n越大，頻率分辨力越高。由圖4.2可知，對同一類型的窗函數(shù)，n2=100時的頻率分辨力比n1=40時頻率分辨力高。不同類型的窗函數(shù)主瓣和旁瓣不同。矩形窗、漢寧窗、三角窗相比，矩形窗的主瓣最窄，所以從圖4.2可以看出，用矩形窗截斷時，在原信號理想頻率點位置能量最集中，而漢寧窗和三角窗主瓣比矩形窗寬，所以主能量分布也較寬。同時，矩形窗旁瓣較大，從頻譜圖中可以看出，矩形窗截斷時由旁瓣引起的能量泄露較大。漢寧窗和三角窗相比，三角窗截斷后由于旁瓣的影響，旁瓣處的頻率幅值出現(xiàn)了一定的波動。所以，為了提高頻率分辨力，在對信號處理的過程中，我們可以選擇增

44、大采樣點數(shù)，同時，為了減少能量泄露，對不同的信號選取不同的窗函數(shù)。作業(yè)3第1題 用矩形窗設計一個線性相位低通濾波器hd(ej )=e-j，|c0，|c|1）寫出h（n）的表達式，確定 與n的關系；2）問有幾種類型，分別屬于哪一種線性相位濾波器；2）若改用漢寧窗設計，寫出 h（n）的表達式。若n=11，=0.25，用matlab程序繪制出濾波器的頻率特性。(1) 解：由理想的單位沖激響應函數(shù) hdn=12-hd(ej)ejnd，得：hdn=12-cce-jejnd=sincn-n- 積分結(jié)果表明hdn是以為中心的偶對稱無限長序列，為得到設計所要求的線性相位，?。?n-12，得理想單位沖激響應：h

45、dn=sincn-n-12n-n-12選矩形窗截斷，得矩形窗設計fir數(shù)字濾波器單位沖激響應：hn=hdnrn(n)=sincn-n-12n-n-12(2) 當n為奇數(shù)時，相頻特性()為：=-n-12，具有嚴格的線性相位關系；當n為偶數(shù)時，相頻特性()為：=-n-12，這種情況仍然保持嚴格的線性相位特性，但由于n為偶數(shù)，相位常數(shù)n-12已經(jīng)不是整數(shù)。綜上，根據(jù)n取值為奇數(shù)還是偶數(shù)，可分為型和型兩種類型fir數(shù)字濾波器。(3) 當改用漢寧窗設計時，漢寧窗函數(shù)為：wn=121-cos2nn-1rnn，此時：hn=hdnw(n)=12sincn-n-12n-n-121-cos2nn-1當n為奇數(shù)時，

46、計算其頻率響應函數(shù)為：hej=m=0n-12-12sincm-n-12m-n-12cosm-n-12121-cos(2mn-1)+ce-jn-12將n=11，c=0.25，代入上式得：hej=m=042sin0.25m-5(m-5)cos(m-5)121-cosm5+0.25e-5j 解：用matlab編寫矩形窗和漢寧窗設計的fir數(shù)字濾波器程序如下：%作業(yè)3_1.mclcclear allclose all%基本參數(shù)設置n=11;wc=0.25*pi;n=0:n-1;h1=fir1(n-1,wc/pi,boxcar(n);%矩形窗設計fir線性相位低通濾波器h2=fir1(n-1,wc/pi

47、,hanning(n);%漢寧窗設計fir線性相位低通濾波器w=0:0.001:pi;h1=freqz(h1,1,w);%矩形窗濾波器的頻率響應h2=freqz(h2,1,w);%漢寧窗濾波器的頻率響應dbh1=20*log10(abs(h1);%矩形窗濾波器的幅頻特性dbh2=20*log10(abs(h2);%漢寧窗窗濾波器的幅頻特性q1=unwrap(angle(h1);%矩形窗濾波器的相頻特性q2=unwrap(angle(h2);%漢寧窗濾波器的相頻特性figure(1);plot(w/pi,dbh1,-.r,w/pi,dbh2,b,linewidth,2);grid on;box

48、off;legend(矩形窗設計fir濾波器幅頻特性,漢寧窗設計fir濾波器幅頻特性);xlabel(頻率(單位：）);ylabel(幅頻特性/db);figure(2)plot(w/pi,q1,-.r,w/pi,q2,b,linewidth,2);grid on;box off;legend(矩形窗設計fir濾波器相頻特性,漢寧窗設計fir濾波器相頻特性);xlabel(頻率(單位：）);ylabel(相頻特性);運行上述作業(yè)3_1.m程序，得到分別用矩形窗和漢寧窗設計的fir數(shù)字濾波器的幅頻特性和相頻特性分別如圖1.1和圖1.2所示：圖1.1 矩形窗和漢寧窗設計的fir濾波器幅頻特性圖圖1

49、.2 矩形窗和漢寧窗設計的fir濾波器相頻特性圖分析：由圖1.1可知，改善窗函數(shù)的形狀可以改善濾波器的幅頻特性和相頻特性，即可以改變?yōu)V波器的濾波性能。從矩形窗和漢寧窗設計的fir數(shù)字濾波器的幅頻特性曲線可以看出，用漢寧窗設計的濾波器旁瓣的衰減比矩形窗的大，但主瓣比矩形窗的寬。在實際設計濾波器時，我們希望窗函數(shù)的主瓣盡可能的窄，以獲得較陡的過渡帶；盡量減少窗函數(shù)頻譜的最大旁瓣的相對幅度，以減少帶內(nèi)波動和加大阻帶內(nèi)的衰減。但主瓣寬度和旁瓣寬度往往是一對矛盾體，二者不可能同時達到最佳，所以，在實際應用過程中，要根據(jù)不同的信號選取不同的窗函數(shù)。第2題 已知函數(shù)： ft=sin200t自行設計數(shù)字分析參

50、數(shù)，實現(xiàn)以下功能。 1）利用matlab繪制出其波形；2）對其進行頻譜分析，做出幅頻譜圖；3）在頻譜分析基礎上，對其進行反變換，并繪制反變換后波形，并與原波形相對比。解：對函數(shù)ft=sin200t進行時域波形繪制、頻譜分析和反傅里葉變換的matlab程序如下：%作業(yè)3_2.mclcclear allclose all%基本參數(shù)設置w0=200*pi,f0=100;%諧波信號頻率fs=1000,ts=1/fs;%采樣頻率n=0.6/ts;%采樣點數(shù)n=0:n-1;%t=n/fs;%定義窗函數(shù)寬度為t=linspace(0,0.6,600)win=window(rectwin,length(t);

51、%定義矩形窗函數(shù)f=sin(w0*t);%諧波信號表達式%繪制原信號時域波形圖subplot(311)plot(t,f);box off;title(連續(xù)諧波信號ft時域波形);xlabel(時間/s);ylabel(幅值);axis(0,0.3,-1,1);%對信號進行頻譜分析freqtr=n*fs/n;%0.6s的頻率轉(zhuǎn)換fft=fft(f,n);subplot(312)plot(freqtr,abs(fft)*2/n);box off;title(連續(xù)諧波信號fft后的幅頻譜);xlabel(頻率/hz);ylabel(幅值f(jomega);axis(0,400,0,1);%對fft變

52、換后的信號進行ifft變換f1=ifft(fft,n);%反傅里葉變換magx=real(f1);subplot(313);plot(t,magx);title(ifft后的時域波形圖);box off;axis(0,0.3,-1,1);xlabel(時間/s);ylabel(幅值);運行上述作業(yè)3_2.m程序，得到連續(xù)諧波信號ft=sin200t的時域波形、幅頻譜圖和反傅里葉變換后的時域波形圖如圖2所示：圖2 連續(xù)諧波信號ft的時域波形、幅頻譜圖和反傅里葉變換后的時域波形圖分析：在對連續(xù)信號進行數(shù)字分析時，首先要對時域上無限長的諧波信號ft進行截斷和采樣，在本題中選取矩形窗對ft進行截斷，這是由于矩形窗的主瓣窄，對信號進行截斷后，信號在頻域能量在主瓣內(nèi)比較集中。在選擇采樣頻率時，根據(jù)采樣定理fs2fmax，原信號為單頻信