DFT頻譜分析作業(yè)

1、DFT頻譜分析作業(yè)通信工程 2013040150作業(yè)一設(shè)周期序列xn=cos0.48n+cos0.52n截取 N 點(diǎn)長(zhǎng)得到 x(n)=xnRN(n)（1）N=10，做10點(diǎn)DFT，得到X1(k);（2）N=10，做100點(diǎn)補(bǔ)零DFT,得到X2(k);（3）N=100，做100點(diǎn)DFT,得到X3(k)。要求：針對(duì)以上三種情況，分別輸出| X1(k)|、| X2(k)|、| X3(k)|的圖形，并進(jìn)行比較、分析和討論。解：代碼如下%雙頻周期信號(hào)的DFT分析n=0:1000;x=cos(0.48*pi*n)+cos(0.52*pi*n);%時(shí)域的離散時(shí)間信號(hào)%-N=10 10點(diǎn)DFT-%N=10;%

2、截取的序列長(zhǎng)度N0=10;%DFT點(diǎn)數(shù) N0Nxn=x(1:N);%取一個(gè)周期的長(zhǎng)度Xk=fft(xn,N0);subplot(3,1,1)stem(0:1:N0-1,abs(Xk),fill,MarkerSize,3.75)%取模畫圖title(N=10 10點(diǎn)DFT)xlabel(k);ylabel(|X1(k)|)grid%-N=10 100點(diǎn)補(bǔ)零DFT-%N=10;%截取的序列長(zhǎng)度N0=100;%DFT點(diǎn)數(shù) N0Nxn=x(1:N);%取一個(gè)周期的長(zhǎng)度Xk=fft(xn,N0);subplot(3,1,2)stem(0:1:N0-1,abs(Xk),fill,MarkerSize,3.

3、75)%取模畫圖gridtitle(N=10 100點(diǎn)補(bǔ)零DFT)xlabel(k);ylabel(|X2(k)|)%-N=100 100點(diǎn)DFT-%N=100;%截取的序列長(zhǎng)度N0=100;%DFT點(diǎn)數(shù) N0Nxn=x(1:N);%取一個(gè)周期的長(zhǎng)度Xk=fft(xn,N0);subplot(3,1,3)stem(0:1:N0-1,abs(Xk),fill,MarkerSize,3.75)%取模畫圖gridtitle(N=100 100點(diǎn)DFT)xlabel(k);ylabel(|X3(k)|)運(yùn)行,得結(jié)果如圖1所示：圖 1 DFT頻譜觀察對(duì)比：1) 10點(diǎn)DFT的| X1(k)|出現(xiàn)明顯的柵

4、欄效應(yīng)，相鄰譜線間隔很大，而100點(diǎn)DFT的間隔較小，說(shuō)明通過(guò)增大DFT點(diǎn)數(shù)N0即頻域抽樣點(diǎn)數(shù)，可以減小柵欄效應(yīng)；2) | X2(k)|與| X3(k)|相比出現(xiàn)了明顯的頻譜泄漏，即截?cái)嘈?yīng)。由于DFT需要對(duì)時(shí)域序列進(jìn)行截?cái)啵瑢?duì)x(n)乘一個(gè)N點(diǎn)長(zhǎng)的矩形窗函數(shù)RN(n) (主瓣寬度為4/N)， N越大主瓣越小引起的泄漏就越小。作業(yè)二設(shè)xat =cos200t+sin100t+cos(50t)。用FFT分析其頻譜結(jié)構(gòu)，選擇不同的截取長(zhǎng)度Tp，觀察截?cái)嘈?yīng)，試用加窗的方法降低旁瓣。解：加矩形窗的FFT代碼如下%單頻周期信號(hào)的DFT分析 cos(pi * t) 0.5Hzfs=400;% 采樣頻率H

5、Zf1=100;f2=50;f3=25;%模擬頻率HZn=0:1000;x=cos(2*pi*f1*n/fs)+sin(2*pi*f2*n/fs)+cos(2*pi*f3*n/fs);%時(shí)域的離散時(shí)間信號(hào)%-Tp=0.04Tp=0.04;N=fs*Tp;%截取的序列長(zhǎng)度N0=2048;%DFT點(diǎn)數(shù) N0Nxn=x(1:N);%取一個(gè)周期的長(zhǎng)度Xk=fft(xn,N0);Xk=abs(Xk)./max(abs(Xk);subplot(3,1,1)stem(0:1:N0-1,Xk,fill,MarkerSize,3.75)%取模畫圖gridtitle(Tp=0.04s);xlabel(k);yla

6、bel(|X(k)|)axis(0 2200 0 1.2)%-Tp=4*0.04Tp=4*0.04;N=fs*Tp;%截取的序列長(zhǎng)度N0=2048;%DFT點(diǎn)數(shù) N0Nxn=x(1:N);%取一個(gè)周期的長(zhǎng)度Xk=fft(xn,N0);Xk=abs(Xk)./max(abs(Xk);subplot(3,1,2)stem(0:1:N0-1,Xk,fill,MarkerSize,3.75)%取模畫圖gridtitle(Tp=4*0.04s);xlabel(k);ylabel(|X(k)|)axis(0 2200 0 1.2)%-Tp=16*0.04Tp=16*0.04;N=fs*Tp;%截取的序列長(zhǎng)

7、度N0=2048;%DFT點(diǎn)數(shù) N0Nxn=x(1:N);%取一個(gè)周期的長(zhǎng)度Xk=fft(xn,N0);Xk=abs(Xk)./max(abs(Xk);subplot(3,1,3)stem(0:1:N0-1,Xk,fill,MarkerSize,3.75)%取模畫圖gridtitle(Tp=16*0.04s);xlabel(k);ylabel(|X(k)|)axis(0 2200 0 1.2)運(yùn)行,得圖 2 DFT頻譜圖 3 DFT頻譜（模擬頻率）由模擬頻率圖可以看出，在做到Tp=160.04S 也即N取256時(shí)，出現(xiàn)了三個(gè)頻點(diǎn)，就是原信號(hào)的三個(gè)頻點(diǎn)100、50、25。當(dāng)Tp=0.04S、40

8、.04S也即N取16、64時(shí)頻率分辨率較低。加Hamming窗進(jìn)行分析：%單頻周期信號(hào)的DFT分析 cos(pi * t) 0.5Hzfs=400;% 采樣頻率HZf1=100;f2=50;f3=25;%模擬頻率HZn=0:1000;x=cos(2*pi*f1*n/fs)+sin(2*pi*f2*n/fs)+cos(2*pi*f3*n/fs);%時(shí)域的離散時(shí)間信號(hào)%-Tp=0.04Tp=0.04;N=fs*Tp;%截取的序列長(zhǎng)度N0=2048;%DFT點(diǎn)數(shù) N0Nxn=x(1:N).*hamming(N);%取一個(gè)周期的長(zhǎng)度Xk=fft(xn,N0);Xk=abs(Xk)./max(abs(X

9、k);subplot(3,1,1)stem(0:1:N0-1,Xk,fill,MarkerSize,3.75)%取模畫圖gridtitle(Tp=0.04s);xlabel(k);ylabel(|X(k)|)axis(0 2200 0 1.2)%-Tp=4*0.04Tp=4*0.04;N=fs*Tp;%截取的序列長(zhǎng)度N0=2048;%DFT點(diǎn)數(shù) N0Nxn=x(1:N).*hamming(N);%取一個(gè)周期的長(zhǎng)度Xk=fft(xn,N0);Xk=abs(Xk)./max(abs(Xk);subplot(3,1,2)stem(0:1:N0-1,Xk,fill,MarkerSize,3.75)%取

10、模畫圖gridtitle(Tp=4*0.04s);xlabel(k);ylabel(|X(k)|)axis(0 2200 0 1.2)%-Tp=16*0.04Tp=16*0.04;N=fs*Tp;%截取的序列長(zhǎng)度N0=2048;%DFT點(diǎn)數(shù) N0Nxn=x(1:N).*hamming(N);%取一個(gè)周期的長(zhǎng)度Xk=fft(xn,N0);Xk=abs(Xk)./max(abs(Xk);subplot(3,1,3)stem(0:1:N0-1,Xk,fill,MarkerSize,3.75)%取模畫圖gridtitle(Tp=16*0.04s);xlabel(k);ylabel(|X(k)|)axi

11、s(0 2200 0 1.2)圖 4 加Hamming窗的DFT頻譜對(duì)比可以明顯看出，加窗后得到的頻譜圖上非信號(hào)頻點(diǎn)處幅度很低。這說(shuō)明旁瓣確實(shí)低了。下圖對(duì)比更明顯： 右側(cè)的為加矩形窗，左側(cè)為加海明窗。也可以看出雖然海明窗降低了旁瓣幅度、但是主瓣寬度變寬，主瓣泄漏厲害，本來(lái)能夠分開的兩個(gè)頻率分量幅值很相近。這顯示了加窗的矛盾。作業(yè)三讀取一個(gè)音樂文件并播放，取其中一段數(shù)據(jù)顯示其時(shí)域波形，利用DFT進(jìn)行頻譜分析并顯示其幅度譜。解：clear all%聲音文件的DFT分析%data,fs,bits=wavread(warning);data,fs,bits=wavread(clock);sound(d

12、ata,fs)subplot(2,1,1)plot(data);axis(0 length(data) -1.2 1.2)title(聲音信號(hào)波形圖)N=65536;%截取的序列長(zhǎng)度N0=65536*2;%DFT點(diǎn)數(shù) N0Nxn=data(1:N);%取一個(gè)周期的長(zhǎng)度Xk=fft(xn,N0);subplot(2,1,2)plot(0:1:N0-1,abs(Xk)%fill,MarkerSize,3.75)%取模畫圖title(N=,num2str(N), ,num2str(N0), 點(diǎn)DFT)xlabel(k);ylabel(|X1(k)|)grid% subplot(3,1,3)% Xa=abs(Xk);% Xkdb=20*log10(Xa/max(Xa)% stem(0:1:N0-1,Xkdb,fill,MarkerSize,3.75)%取模畫圖% title(N=,num2str(N), ,num2str(N0), 點(diǎn)DFT)% axis(0 N0+5 -50 0)% xlabel(k);ylabe