DSP大作業(yè)--快速傅立葉變換實(shí)驗(yàn)與設(shè)計(jì)

1、DSP原理及應(yīng)用大作業(yè) 快速傅立葉變換 專業(yè)：XXXX姓名：XXX 學(xué)號(hào)：08201081XX指導(dǎo)老師：XX時(shí)間：2XXXX快速傅立葉變換（FFT）實(shí)驗(yàn)一、設(shè)計(jì)目的 1.在理論學(xué)習(xí)的基礎(chǔ)上，通過本實(shí)驗(yàn)，加深對(duì)FFT的理解，熟悉FFT子程序。 2.熟悉應(yīng)用FFT對(duì)典型信號(hào)進(jìn)行頻譜分析的方法 3.了解應(yīng)用FFT進(jìn)行信號(hào)頻譜分析過程中可能出現(xiàn)的問題以便在實(shí)際中正確應(yīng)用FFT。4.掌握用窗函數(shù)法設(shè)計(jì) FFT 快速傅里葉的原理和方法； 5熟悉 FFT 快速傅里葉特性；二、所需設(shè)備 PC兼容機(jī)一臺(tái)，操作系統(tǒng)為 Windows2000( 或 Windows98 ， WindowsXP ，以下默認(rèn)為 Wind

2、ows2000)，安裝 Code Composer Studio 2.0 軟件。三、設(shè)計(jì)內(nèi)容 本試驗(yàn)要求使用 FFT 變換求一個(gè)時(shí)域信號(hào)的頻域特性，并從這個(gè)頻域特性求出該信號(hào)的頻率值。使用c語(yǔ)言實(shí)現(xiàn)對(duì) FFT 算法的仿真，然后使用 DSP 匯編語(yǔ)言實(shí)現(xiàn)對(duì) FFT 的 DSP編程。本實(shí)驗(yàn)采用軟件仿真，不需設(shè)置硬件。四、設(shè)計(jì)原理 在各種信號(hào)序列中，有限長(zhǎng)序列信號(hào)處理占有很重要地位，對(duì)有限長(zhǎng)序列，我們可以使用離散Fouier變換(DFT)。這一變換不但可以很好的反映序列的頻譜特性，而且易于用快速算法在計(jì)算機(jī)上實(shí)現(xiàn)，當(dāng)序列x(n)的長(zhǎng)度為N時(shí)，它的DFT定義為：，反換為：有限長(zhǎng)序列的DFT是其Z變換在

3、單位圓上的等距采樣，或者是序列Fourier變換的等距采樣，因此可以用于序列的譜分析。 FFT并不是與DFT不同的另一種變換，而是為了減少DFT運(yùn)算次數(shù)的一種快速算法。它是對(duì)變換式進(jìn)行一次次分解，使其成為若干小點(diǎn)數(shù)的組合，從而減少運(yùn)算量。常用的FFT是以2為基數(shù)的，其長(zhǎng)度 N=2L，它的效率高，程序簡(jiǎn)單使用非常方便，當(dāng)要變換的序列長(zhǎng)度不等于2的整數(shù)次方時(shí)，為了使用以2為基數(shù)的FFT，可以用末位補(bǔ)零的方法，使其長(zhǎng)度延長(zhǎng)至2的整數(shù)次方。在運(yùn)用DFT進(jìn)行頻譜分析的過程中可能產(chǎn)生幾種問題：(1) 混疊    序列的頻譜時(shí)被采樣信號(hào)的周期延拓，當(dāng)采樣速率不滿足Ny

4、quist定理時(shí)，就會(huì)發(fā)生頻譜混疊，使得采樣后的信號(hào)序列頻譜不能真實(shí)的反映原信號(hào)的頻譜。 避免混疊現(xiàn)象的唯一方法是保證采樣速率足夠高，使頻譜混疊現(xiàn)象不致出現(xiàn)，即在確定采樣頻率之前，必須對(duì)頻譜的性質(zhì)有所了解，在一般情況下，為了保證高于折疊頻率的分量不會(huì)出現(xiàn)，在采樣前，先用低通模擬濾波器對(duì)信號(hào)進(jìn)行濾波。(2) 泄漏 實(shí)際中我們往往用截短的序列來近似很長(zhǎng)的甚至是無限長(zhǎng)的序列，這樣可以使用較短的DFT來對(duì)信號(hào)進(jìn)行頻譜分析，這種截短等價(jià)于給原信號(hào)序列乘以一個(gè)矩形窗函數(shù)，也相當(dāng)于在頻域?qū)⑿盘?hào)的頻譜和矩形窗函數(shù)的頻譜卷積，所得的頻譜是原序列頻譜的擴(kuò)展。 泄漏不能與混疊完全分開，因?yàn)樾孤?dǎo)致頻譜的擴(kuò)

5、展，從而造成混疊。為了減少泄漏的影響，可以選擇適當(dāng)?shù)拇昂瘮?shù)使頻譜的擴(kuò)散減至最小。 DFT是對(duì)單位圓上Z變換的均勻采樣，所以它不可能將頻譜視為一個(gè)連續(xù)函數(shù)，就一定意義上看，用DFT來觀察頻譜就好像通過一個(gè)柵欄來觀看一個(gè)圖景一樣，只能在離散點(diǎn)上看到真實(shí)的頻譜，這樣就有可能發(fā)生一些頻譜的峰點(diǎn)或谷點(diǎn)被“尖樁的柵欄”所攔住，不能別我們觀察到。    減小柵欄效應(yīng)的一個(gè)方法就是借助于在原序列的末端填補(bǔ)一些零值，從而變動(dòng)DFT的點(diǎn)數(shù)，這一方法實(shí)際上是人為地改變了對(duì)真實(shí)頻譜采樣的點(diǎn)數(shù)和位置，相當(dāng)于搬動(dòng)了每一根“尖樁柵欄”的位置，從而使得頻譜的峰點(diǎn)或谷點(diǎn)暴露出來。 用FF

6、T可以實(shí)現(xiàn)兩個(gè)序列的圓周卷積。在一定的條件下，可以使圓周卷積等于線性卷積。一般情況，設(shè)兩個(gè)序列的長(zhǎng)度分別為N1和N2，要使圓周卷積等于線性卷積的充要條件是FFT的長(zhǎng)度 NN1N2 對(duì)于長(zhǎng)度不足N的兩個(gè)序列，分別將他們補(bǔ)零延長(zhǎng)到N。 當(dāng)兩個(gè)序列中有一個(gè)序列比較長(zhǎng)的時(shí)候，我們可以采用分段卷積的方法。有兩種方法： 重疊相加法。將長(zhǎng)序列分成與短序列相仿的片段，分別用FFT對(duì)它們作線性卷積，再將分段卷積各段重疊的部分相加構(gòu)成總的卷積輸出。 重疊保留法。這種方法在長(zhǎng)序列分段時(shí)，段與段之間保留有互相重疊的部分，在構(gòu)成總的卷積輸出時(shí)只需將各段線性卷積部分直接連接起來，省掉了輸出段的直接相加。 (3) 柵欄效應(yīng)

7、 DFT是對(duì)單位圓上z變換的均勻采樣，所以它不可能將頻譜視為一個(gè)連續(xù) 函數(shù)，從某種意義上講，用DFT來觀察頻譜就如同通過一個(gè)柵欄來觀看景象一 樣，只能在離散點(diǎn)上看到真實(shí)的頻譜，這樣一些頻譜的峰點(diǎn)或谷點(diǎn)就可能被"尖 樁的柵欄"擋住，也就是正好落在兩個(gè)離散采樣點(diǎn)之間，不能被觀察到。減小柵欄效應(yīng)的一個(gè)方法是在原序列的末端填補(bǔ)一些零值，從而變動(dòng)DFT 的點(diǎn)數(shù)，這一方法實(shí)際上是人為地改變了對(duì)真實(shí)頻譜采樣的點(diǎn)數(shù)和位置，相當(dāng)于 搬動(dòng)了"尖樁柵欄"的位置，從而使得頻譜的峰點(diǎn)或谷點(diǎn)暴露出來。(4) DFT的分辨率填補(bǔ)零值可以改變對(duì)DTFT的采樣密度，人們常常有一種誤解，認(rèn)

8、為補(bǔ)零可以提高DFT的頻率分辨率，事實(shí)上，DFT的頻率分辨率通常規(guī)定為，這里的N是指信號(hào)的有效長(zhǎng)度，而不是補(bǔ)零的長(zhǎng)度。不同長(zhǎng)度的，其DTFT的結(jié)果是不同的;而相同長(zhǎng)度的盡管補(bǔ)零的長(zhǎng)度不同其DTFT的結(jié)果應(yīng)是相同的，它們的DFT只是反映了對(duì)相同的DTFT采用了不同的采樣密度?？偨Y(jié)一下: 要提高DFT分辨率只有增加信號(hào)的截取長(zhǎng)度N。五、實(shí)驗(yàn)內(nèi)容1、原來的程序下做出的圖將程序改為四個(gè)節(jié)點(diǎn)時(shí)的程序如下：#include "myapp.h"#include "csedu.h"#include "scancode.h"#include <m

9、ath.h>#define PI 3.1415926#define SAMPLENUMBER 128void InitForFFT();void MakeWave();int INPUTSAMPLENUMBER,DATASAMPLENUMBER;float fWaveRSAMPLENUMBER,fWaveISAMPLENUMBER,wSAMPLENUMBER;float sin_tabSAMPLENUMBER,cos_tabSAMPLENUMBER;main()int i;InitForFFT();MakeWave();for ( i=0;i<SAMPLENUMBER;i+ )fW

10、aveRi=INPUTi;fWaveIi=0.0f;wi=0.0f;FFT(fWaveR,fWaveI);for ( i=0;i<SAMPLENUMBER;i+ )DATAi=wi;while ( 1 );/ break pointvoid FFT(float dataRSAMPLENUMBER,float dataISAMPLENUMBER)int x0,x1,x2,x3,xx;int i,j,k,b,p,L;float TR,TI,temp;/* following code invert sequence */for ( i=0;i<SAMPLENUMBER;i+ )x0=x

11、1=x2=x3=0;x0=i&0x01; x1=(i/2)&0x01; x2=(i/4)&0x01; x3=(i/8)&0x01;xx=x0*8+x1*4+x2*2+x3;dataIxx=dataRi;for ( i=0;i<SAMPLENUMBER;i+ )dataRi=dataIi; dataIi=0; /* following code FFT */for ( L=1;L<=4;L+ ) /* for(1) */b=1; i=L-1;while ( i>0 ) b=b*2; i-; /* b= 2(L-1) */for ( j=0;j&l

12、t;=b-1;j+ ) /* for (2) */p=1; i=4-L;while ( i>0 ) /* p=pow(2,8-L)*j; */p=p*2; i-;p=p*j;for ( k=j;k<128;k=k+2*b ) /* for (3) */TR=dataRk; TI=dataIk; temp=dataRk+b;dataRk=dataRk+dataRk+b*cos_tabp+dataIk+b*sin_tabp;dataIk=dataIk-dataRk+b*sin_tabp+dataIk+b*cos_tabp;dataRk+b=TR-dataRk+b*cos_tabp-da

13、taIk+b*sin_tabp;dataIk+b=TI+temp*sin_tabp-dataIk+b*cos_tabp; /* END for (3) */ /* END for (2) */ /* END for (1) */for ( i=0;i<SAMPLENUMBER/2;i+ ) wi=sqrt(dataRi*dataRi+dataIi*dataIi); /* END FFT */void InitForFFT()int i;for ( i=0;i<SAMPLENUMBER;i+ )sin_tabi=sin(PI*2*i/SAMPLENUMBER);cos_tabi=co

14、s(PI*2*i/SAMPLENUMBER);void MakeWave()int i;for ( i=0;i<SAMPLENUMBER;i+ )INPUTi=sin(PI*2*i/SAMPLENUMBER*3)*1024;得出的圖為： 程序改為八個(gè)節(jié)點(diǎn)時(shí)：#include "myapp.h"#include "csedu.h"#include "scancode.h"#include <math.h>#define PI 3.1415926#define SAMPLENUMBER 128void InitForFF

15、T();void MakeWave();int INPUTSAMPLENUMBER,DATASAMPLENUMBER;float fWaveRSAMPLENUMBER,fWaveISAMPLENUMBER,wSAMPLENUMBER;float sin_tabSAMPLENUMBER,cos_tabSAMPLENUMBER;main()int i;InitForFFT();MakeWave();for ( i=0;i<SAMPLENUMBER;i+ )fWaveRi=INPUTi;fWaveIi=0.0f;wi=0.0f;FFT(fWaveR,fWaveI);for ( i=0;i<

16、;SAMPLENUMBER;i+ )DATAi=wi;while ( 1 );/ break pointvoid FFT(float dataRSAMPLENUMBER,float dataISAMPLENUMBER)int x0,x1,x2,x3,x4,x5,x6,x7,xx;int i,j,k,b,p,L;float TR,TI,temp;/* following code invert sequence */for ( i=0;i<SAMPLENUMBER;i+ )x0=x1=x2=x3=x4=x5=x6=0;x0=i&0x01; x1=(i/2)&0x01; x2

17、=(i/4)&0x01; x3=(i/8)&0x01;x4=(i/16)&0x01; x5=(i/32)&0x01; x6=(i/64)&0x01;x7=(i/128)&x01;xx=x0*128+x1*64+x2*32+x3*16+x4*8+x5*4+x6*2+x7;dataIxx=dataRi;for ( i=0;i<SAMPLENUMBER;i+ )dataRi=dataIi; dataIi=0; /* following code FFT */for ( L=1;L<=8;L+ ) /* for(1) */b=1; i=L-1

18、;while ( i>0 ) b=b*2; i-; /* b= 2(L-1) */for ( j=0;j<=b-1;j+ ) /* for (2) */p=1; i=8-L;while ( i>0 ) /* p=pow(2,8-L)*j; */p=p*2; i-;p=p*j;for ( k=j;k<128;k=k+2*b ) /* for (3) */TR=dataRk; TI=dataIk; temp=dataRk+b;dataRk=dataRk+dataRk+b*cos_tabp+dataIk+b*sin_tabp;dataIk=dataIk-dataRk+b*sin_tabp+dataIk+b*cos_tabp;data