武漢理工大學(xué)8點基于DIFFFT

武漢理工大學(xué)8點基于DIFFFT武漢理工大學(xué)8點基于DIFFFT/武漢理工大學(xué)8點基于DIFFFT課程設(shè)計任務(wù)書學(xué)生姓名：李嘉辛專業(yè)班級：電信1206指導(dǎo)教師：黃朝兵工作單位：信息工程學(xué)院題目：8點基于DIF的FFT的實現(xiàn)初始條件：具備Matlab編程能力；熟悉基于DIF的FFT的實現(xiàn)原理；提供編程所需要的計算機(jī)一臺要求完成的主要任務(wù)：（包括課程設(shè)計工作量及其技術(shù)要求，以及說明書撰寫等具體要求）1、獨立編寫一個8點的基于DIF的FFT實現(xiàn)程序，不能使用matlab自帶的FFT實現(xiàn)函數(shù)2、并調(diào)用該函數(shù)實現(xiàn)16點的FFT運算，用matlab自帶函數(shù)對運行結(jié)果結(jié)果進(jìn)行驗證3、完成符合學(xué)校要求的設(shè)計說明書時間安排：一周，其中3天程序設(shè)計，2天程序調(diào)試指導(dǎo)教師簽名：年月日系主任（或責(zé)任教師）簽名：年月日目錄摘要 11概述 11.1數(shù)字信號處理定義 11.2數(shù)字信號處理的特點及實現(xiàn)方法 12理論分析 22.1DFT的定義 22.2直接計算DFT的問題及FFT思想 22.3基2按時間抽?。―IT）的FFT算法 22.4基2按頻率抽取（DIF）的FFT算法 42.5按頻率抽取的FFT的特點 62.5.1原位運算 62.5.2蝶形運算兩節(jié)點之間的“距離” 62.5.3旋轉(zhuǎn)因子的變化規(guī)律 63程序設(shè)計 73.1變址 73.2L級遞推計算 74結(jié)果及分析 95心得體會 12參考文獻(xiàn) 13摘要快速傅里葉變換（FFT）是離散傅里葉變換（DFT）的快速算法，F(xiàn)FT算法通過利用旋轉(zhuǎn)因子的性質(zhì)，將一個大點數(shù)DFT化成幾個小點數(shù)DFT，就可以大大減少運算量。DIF-FFT是利用頻率抽選的FFT算法，在Matlab中可以通過三重循環(huán)語句實現(xiàn)。關(guān)鍵詞：FFT，蝶形運算，倒序排列1概述1.1數(shù)字信號處理定義數(shù)字信號是用數(shù)字序列表示的信號，數(shù)字信號處理就是通過計算機(jī)或?qū)Ｓ锰幚碓O(shè)備，用數(shù)值計算等數(shù)字方式對數(shù)字序列進(jìn)行各種處理，將數(shù)字信號變換成符合要求的某種形式。數(shù)字信號處理主要包括數(shù)字濾波和數(shù)字頻譜分析兩大部分。例如，對數(shù)字信號進(jìn)行濾波，限制其頻帶或濾除噪聲和干擾，以提取和增強信號的有用分量；對信號進(jìn)行頻譜分析或功率分析，了解信號的頻譜組成，以對信號進(jìn)行識別。當(dāng)然，凡是用數(shù)字方式對信號進(jìn)行濾波，變換，增強，壓縮，估計和識別等都是數(shù)字信號處理的研究范疇。數(shù)字信號處理在理論上所涉及的范圍及其廣泛。數(shù)字領(lǐng)域中的微積分，概率統(tǒng)計，隨機(jī)過程，高等代數(shù)，數(shù)值分析，復(fù)變函數(shù)和各種變換（如傅里葉變換，Z變換，離散傅里葉變換，小波變換等）都是它的基本工具，網(wǎng)絡(luò)理論，信號及系統(tǒng)等則是它的理論基礎(chǔ)。在科學(xué)發(fā)展上，數(shù)字信號處理又和最優(yōu)控制，通信理論等緊密相連，目前已成為人工智能，模式識別，神經(jīng)網(wǎng)絡(luò)等新興學(xué)科的重要理論基礎(chǔ)，其實現(xiàn)技術(shù)又和計算機(jī)科學(xué)和微電子技術(shù)密不可分。因此，數(shù)字信號處理是把經(jīng)典的理論基礎(chǔ)體系作為自身的理論基礎(chǔ)，同時又使自己成為一系列新興學(xué)科的理論基礎(chǔ)。1.2數(shù)字信號處理的特點及實現(xiàn)方法及模擬信號處理相比，數(shù)字信號處理具有高精度、高穩(wěn)定性、靈活性好、易于大規(guī)模集成等顯著的優(yōu)點。數(shù)字信號處理的主要研究對象是數(shù)字信號，且采用數(shù)值運算的方法達(dá)到處理的目的。數(shù)字信號處理的實現(xiàn)方法基本上可以分為軟件實現(xiàn)方法、硬件實現(xiàn)方法和軟硬件想結(jié)合的實現(xiàn)方法。數(shù)字信號處理的理論、算法和實現(xiàn)方法三者是密不可分的。2理論分析2.1DFT的定義對于有限長離散數(shù)字信號{x[n]}，0nN-1，其離散譜{X[k]}可以由離散付氏變換（DFT）求得。DFT的定義為：，k=0,1,…N-1（2.1）通常令，稱為旋轉(zhuǎn)因子。2.2直接計算DFT的問題及FFT思想由DFT的定義可以看出，在x[n]為復(fù)數(shù)序列的情況下，完全直接運算N點DFT需要N-1的2次方復(fù)數(shù)乘法和N（N-1）次加法。因此，對于一些相當(dāng)大的N值（如1024）來說，直接計算它的DFT所作的計算量是很大的。FFT的基本思想在于，將原有的N點序列分成兩個較短的序列，這些序列的DFT可以很簡單的組合起來得到原序列的DFT。例如，若N為偶數(shù)，將原有的N點序列分成兩個（N/2）點序列，那么計算N點DFT將只需要約[(N/2)2?2]=N2/2次復(fù)數(shù)乘法。即比直接計算少做一半乘法。因子(N/2)2表示直接計算（N/2）點DFT所需要的乘法次數(shù)，而乘數(shù)2代表必須完成兩個DFT。上述處理方法可以反復(fù)使用，即（N/2）點的DFT計算也可以化成兩個（N/4）點的DFT（假定N/22.3基2按時間抽?。―IT）的FFT算法設(shè)序列長度為，L為整數(shù)（如果序列長度不滿足此條件，通過在后面補零使其滿足）。將長度為的序列，先按n的奇偶分成兩組：（r=0,1,…,N/2-1）（2.2）DFT化為：(2.3)上式中利用了旋轉(zhuǎn)因子的可約性，即：。又令，則上式可以寫成：(k=0,1,…,N/2-1) (2.4)可以看出，分別為從中取出的N/2點偶數(shù)點和奇數(shù)點序列的N/2點DFT值，所以，一個N點序列的DFT可以用兩個N/2點序列的DFT組合而成。但是，從上式可以看出，這樣的組合僅表示出了前N/2點的DFT值，還需要繼續(xù)利用表示的后半段本算法推導(dǎo)才完整。利用旋轉(zhuǎn)因子的周期性，有：，則后半段的DFT值表達(dá)式：(2.5)（k=0,1,…,N/2-1）(2.6)所以后半段（k=N/2,…,N-1）的DFT值可以用前半段k值表達(dá)式獲得，中間還利用到，得到后半段的值表達(dá)式為：（k=0,1,…,N/2-1）。這樣，通過計算兩個N/2點序列的N/2點DFT，可以組合得到N點序列的DFT值，其組合過程如下圖所示：-1圖2.1FFT形成過程2.4基2按頻率抽取（DIF）的FFT算法設(shè)序列長度為，L為整數(shù)（如果序列長度不滿足此條件，通過在后面補零使其滿足）。在把按k的奇偶分組之前，把輸入按n的順序分成前后兩半：(2.7)因為，則有，所以：X[k]=n=0按k的奇偶來討論，k為偶數(shù)時：X[2r]=n=0k為奇數(shù)時：

X[2r+1]=前面已經(jīng)推導(dǎo)過，所以上面的兩個等式可以寫為：X[2r]=n=0X[2r+1]=W通過上面的推導(dǎo)，X[k]的偶數(shù)點值X[2r]和奇數(shù)點值X[2r+1]分別可以由組合而成的N/2點的序列來求得，其中偶數(shù)點值X[2r]為輸入x[n]的前半段和后半段之和序列的N/2點DFT值，奇數(shù)點值X[2r+1]為輸入x[n]的前半段和后半段之差再及WN令， (2.13) (2.14)則有： (2.15)這樣，也可以用兩個N/2點DFT來組合成一個N點DFT，組合過程如下圖所示：-1圖2.2DIF-FFT算法這樣可以把一個N點DFT分解為兩個N/2點DFT的組合，兩個N/2點DFT還可以繼續(xù)分解，設(shè)N=2M，則經(jīng)過M-1次分解，最后可以分解成為N/2個兩點DFT，可以由一個蝶形運算來求解。例如8點DIF-FFT蝶形運算圖如圖2.2圖2.38點DIF-FFT蝶形圖輸出序列的排列規(guī)律不是從小到大按順序的，而是按照倒敘規(guī)則排序的，即先將0-7轉(zhuǎn)換為二進(jìn)制數(shù)，然后將二進(jìn)制數(shù)左右倒序，再轉(zhuǎn)為十進(jìn)制就可以得到新的數(shù)列，即：0，4，2，6，1，5，3，7。2.5按頻率抽取的FFT的特點2.5.1原位運算在DIF-FFT蝶形圖中，取第m級且兩輸入節(jié)點分別在第k，j行的蝶形為例，討論DIF-FFT的原位運算規(guī)律。由圖可得蝶形運算的關(guān)系式可表示為=，=[]。有上式可得的m-1級的第k行及第j行的輸出，在運算流圖中的作用是，用來計算第m級的第k行和第j行的輸出，，這樣當(dāng)計算完，后，，在運算流圖中就不在起作用，因此可以采用原位運算，把，直接存入原來存放，的存儲單元。同理可以把第m級蝶形的N個輸出值直接存放在第m-1級蝶形輸出的N個存儲單元中，這樣從第一級的輸入x（n）開始到最后一級的輸出X(k)，只需要N個存儲單元。2.5.2蝶形運算兩節(jié)點之間的“距離”第一級蝶形每個蝶形運算量節(jié)點的“距離”為4，第二級每個蝶形運算另節(jié)點的“距離”為2，第三級蝶形每個蝶形運算量節(jié)點的“距離”為1。依次類推：對于N等于2的L次方的DIF-FFT，可以得到第M級蝶形每個蝶形運算量節(jié)點的“距離”為2的L-M次方。2.5.3旋轉(zhuǎn)因子的變化規(guī)律以8點的FFT為例，第一級蝶形，r=0,1，2；第二級蝶形，r=0,1；第三級的蝶形，r=0。依次類推，對于M級蝶形，旋轉(zhuǎn)因子的指數(shù)為r=J?2M-1，J=0,1,2,3，……，這樣就可以算出每一級的旋轉(zhuǎn)因子。對于M級的任一蝶形運算所對應(yīng)的旋轉(zhuǎn)因子的指數(shù)，可以如下方法得到：1將待求的蝶形輸入節(jié)點中上面節(jié)點的行標(biāo)號值k寫成L位二進(jìn)制數(shù)；2將此二進(jìn)制數(shù)乘以2的M-1次方，即將L位二進(jìn)制數(shù)左移M-1位，右邊的空位補零，然后從低位到高位截取L位，即所得指數(shù)r所對應(yīng)的二進(jìn)制數(shù)。3程序設(shè)計FFT程序包括變址（倒位序）和L級遞推計算（N=，L為正整數(shù)）兩部分。3.1變址DIF-FFT是輸出倒位序的變址處理，設(shè)x(i)表示存放自然順序輸入數(shù)據(jù)的內(nèi)存單元，x（j）表示存放倒位序序數(shù)的內(nèi)存單元，I、J=0,1,…,N-1,當(dāng)I=J時，不用變址；當(dāng)IJ時，需要變址；但是當(dāng)I<J時，進(jìn)行變址在先，故在I>J時，就不需要再變址了，否則變址兩次等于不變。其中本程序使用的“反向進(jìn)位加法”。也可用bin2dec函數(shù)可以實現(xiàn)倒以位序。3.2L級遞推計算整個L級遞推過程由三個嵌套循環(huán)構(gòu)成。外層的一個循環(huán)控制L（L=）級的順序運算；內(nèi)層的兩個循環(huán)控制同一級（M相同）各蝶形結(jié)的運算，其中最內(nèi)層循環(huán)控制同一種（即中的r相同）蝶形結(jié)的運算。其循環(huán)變量為I，I用來控制同一種蝶形結(jié)運算。其步進(jìn)值為蝶形結(jié)的間距值LE=，同一種蝶形結(jié)中參加運算的兩節(jié)點的間距為LE1=點。第二層循環(huán)，其循環(huán)變量J用來控制計算不同種（系數(shù)不同）的碟形結(jié)，J的步進(jìn)值為1。也可以看出，最內(nèi)層循環(huán)完成每級的蝶形結(jié)運算，第二層循環(huán)則完成系數(shù)的運算。最外層循環(huán)，用循環(huán)變量M來控制運算的級數(shù)，M為1到L，步進(jìn)值為1，當(dāng)M改變時，則LE1，LE和系數(shù)U都會改變。MATLAB實現(xiàn)的代碼：function[Xk]=DIF_FFT(xn,N);%實現(xiàn)對輸入序列的DIF-FFT的基2算法，點數(shù)小于等于2的次冪n=nextpow2(length(xn));%取最小的滿足2次冪的nN=2^n;iflength(xn)<Nxn=[xn,zeros(1,N-length(xn))];end%蝶形運算算法M=log2(N);%M=3或4%第一層循環(huán)分解步驟form=0:M-1Num_Group=2^m;%每一次分解中組的個數(shù)IntV_Group=N/2^m;%每一次分解中組之間的間距IntV_Unit=N/2^(m+1);%每一組運算單元的間距Count=IntV_Unit-1;%運算單元循環(huán)次數(shù)第三次循環(huán)次數(shù)Wn=exp(-j*2*pi/IntV_Group);%旋轉(zhuǎn)因子%第二層循環(huán)組的循環(huán)forg=1:Num_Groupdx1=(g-1)*IntV_Group;%第g組中運算量1的偏移量dx2=(g-1)*IntV_Group+IntV_Unit;%第g組中運算量2的偏移量%第三層循環(huán)組內(nèi)運算的循環(huán)forr=0:Countk=r+1;%“組內(nèi)”序列的下標(biāo)b=xn(k+dx1);xn(k+dx1)=b+xn(k+dx2);%第m級，第g組的蝶形運算式1xn(k+dx2)=[b-xn(k+dx2)]*Wn^r;%第m級，第g組的蝶形運算式2endendend%序列排序開始n1=fliplr(dec2bin([0:N-1]));%碼位倒置十進(jìn)制轉(zhuǎn)二進(jìn)制并倒序n2=bin2dec(n1);%倒序后二進(jìn)制轉(zhuǎn)十進(jìn)制%倒序循環(huán)輸出fori=1:NXk(i)=xn(n2(i)+1);end4結(jié)果及分析編寫主函數(shù)，在主函數(shù)中輸入一個16點的序列分別調(diào)用自己編寫的FFT函數(shù)，和MATLAB本身系統(tǒng)的FFT函數(shù)并比較兩個結(jié)果是否相等，以判斷自己編寫的FFT程序是否正確xn=[0:15];m=1:16；N=16x1=DIF_FFT(xn,N)x2=fft(xn)x3=abs(x1);x4=abs(x2);x5=angle(x1);x6=angle(x2);figure('NumberTitle','off','Name','電信1206李嘉辛');subplot(3,1,1)stem(m,xn);title('輸入的離散序列')subplot(3,1,2)stem(m,x3);title('經(jīng)過DIF_FFT后得到的頻譜的幅度')subplot(3,1,3)stem(m,x5);title('經(jīng)過DIF_FFT后得到的頻譜的相位')figure('NumberTitle','off','Name','電信1206李嘉辛');subplot(3,1,1)stem(m,xn);title('輸入的離散序列')subplot(3,1,2)stem(m,x4);title('經(jīng)過庫函數(shù)fft后得到的頻譜的幅度')subplot(3,1,3)stem(m,x6);title('經(jīng)過庫函數(shù)fft后得到的頻譜的相位')圖4.1DIF-FFT結(jié)果圖4.2庫函數(shù)FFT結(jié)果通過觀察比較，得到的序列各點的值以及直觀的通過圖形，可以得到自己編寫的DIF_FFT函數(shù)實現(xiàn)對序列進(jìn)行FFT變換得到的結(jié)果及庫函數(shù)FFT得到的結(jié)果是一樣的。說明DIF_FFT子程序是正確的。從圖中也可以看出有限長序列通過FFT后得到的頻域為離散的。從理論講，有限長序列經(jīng)過離散傅里葉變換后，得到的頻譜為離散的，從而也說明了FFT是DFT的優(yōu)化方法，也屬于DFT。這個程序可以實現(xiàn)基2-FFT,但是如果想在運行時直接輸入要變換的點數(shù)就不行，必須在調(diào)用FFT函數(shù)前現(xiàn)將要算的序列定義好，這是這個程序的不足之處。但是該程序可以計算不是2的整數(shù)次冪的序列。所以在主程序中，輸入序列必須給出才能進(jìn)行FFT變換。當(dāng)使用編寫的程序進(jìn)行16點的DIF-FFT計算時結(jié)果如下：》xn=[12345678];N=8;DIF_FFT(xn,N)>>xn=[0:15];N=16;DIF_FFT(xn,N)ans=1.0e+02*Columns1through41.2000+0.0000i-0.0800+0.4022i-0.0800+0.1931i-0.0800+0.1197iColumns5through8-0.0800+0.0800i-0.0800+0.0535i-0.0800+0.0331i-0.0800+0.0159iColumns9through12-0.0800+0.0000i-0.0800-0.0159i-0.0800-0.0331i-0.0800-0.0535iColumns13through16-0.0800-0.0800i-0.0800-0.1197i-0.0800-0.1931i-0.0800-0.4022i當(dāng)調(diào)用matlab自帶的FFT程序進(jìn)行相同的16點的FFT計算時結(jié)果如下：>>xn=[0:15];fft(xn)ans=1.0e+02*Columns1through41.2000+0.0000i-0.0800+0.4022i-0.0800+0.1931i-0.0800+0.1197iColumns5through8-0.0800+0.0800i-0.0800+0.0535i-0.0800+0.0331i-0.0800+0.0159iColumns9through12-0.0800+0.0000i-0.0800-0.0159i-0.0800-0.0331i-0.0800-0.0535iColumns13through16-0.0800-0.0800i-0.0800-0.1