艾 紅，常青青，鄧大偉

（北京信息科技大學(xué) 自動(dòng)化學(xué)院，北京 100192）

0 引言

快速傅立葉變換（FFT）在雷達(dá)、通信、電子對(duì)抗和電力系統(tǒng)等領(lǐng)域有廣泛應(yīng)用，特別是在電力系統(tǒng)的諧波檢測(cè)中，F(xiàn)FT幾乎是唯一可行的檢測(cè)方法。通常提高FFT運(yùn)算速度有兩種途徑：改進(jìn)FFT算法本身和改進(jìn)運(yùn)算工具?，F(xiàn)階段提高FFT算法本身非常困難，一般方法致力于改進(jìn)運(yùn)算工具。數(shù)字信號(hào)處理器DSP 是一種可編程的高性能處理器。文中充分利用TMS320F2812 DSP強(qiáng)大的數(shù)據(jù)處理能力，實(shí)現(xiàn)了FFT運(yùn)算，并提高了運(yùn)算速度。

1 系統(tǒng)硬件結(jié)構(gòu)

系統(tǒng)設(shè)計(jì)以TMS320F2812處理器為核心，輔以外圍電路構(gòu)成。DSP負(fù)責(zé)模擬輸入信號(hào)數(shù)據(jù)采集以及FFT算法實(shí)現(xiàn)。外圍電路包括電源轉(zhuǎn)換電路，時(shí)鐘電路，復(fù)位電路以及外部RAM等。系統(tǒng)的硬件整體結(jié)構(gòu)如圖1所示[1]。對(duì)信號(hào)進(jìn)行FFT變換，首先要對(duì)模擬信號(hào)采樣將其轉(zhuǎn)換為數(shù)字信號(hào)。輸入的連續(xù)模擬信號(hào)經(jīng)信號(hào)調(diào)理電路后輸出到DSP的ADC模擬輸入通道，經(jīng)過(guò)ADC數(shù)據(jù)采集，模數(shù)轉(zhuǎn)換的結(jié)果存放于ADC結(jié)果寄存器中。信號(hào)調(diào)理電路主要是為了信號(hào)的抗混疊濾波以及電路阻抗的匹配 。信號(hào)調(diào)理電路對(duì)輸入信號(hào)進(jìn)行調(diào)理處理，包括信號(hào)的濾波、跟隨輸出以及信號(hào)的穩(wěn)定。DSP對(duì)采集的數(shù)字信號(hào)進(jìn)行FFT運(yùn)算處理，同時(shí)對(duì)運(yùn)算結(jié)果進(jìn)行相應(yīng)的數(shù)據(jù)顯示和數(shù)據(jù)存儲(chǔ)。

圖1 系統(tǒng)的硬件整體結(jié)構(gòu)圖

2 FFT算法原理

FFT是離散傅立葉變換（DFT）的快速運(yùn)算，是數(shù)字信號(hào)處理的基礎(chǔ)。因?yàn)橛行┬盘?hào)在時(shí)域很難看出特性，使用FFT將其變換到頻域，就會(huì)很容易看出其特性。DFT算法的基本公式為[3]：

FFT算法是不斷地把長(zhǎng)序列的DFT分解成幾個(gè)短序列的DFT，并利用的周期性和對(duì)稱性減少DFT的運(yùn)算次數(shù)。設(shè)序列x(n)的長(zhǎng)度為N(N=2M,M為任意正整數(shù))，按n的奇偶把x(n)分解成兩個(gè)N/2點(diǎn)的子序列：

由此可見(jiàn)，若將任何一偶數(shù)點(diǎn)序列按下標(biāo)的奇偶性分成兩個(gè)子序列，則原序列的DFT可由兩子序列的DFT線性組合得到。運(yùn)算流圖如圖2所示。

圖2 運(yùn)算流圖

圖3 蝶式運(yùn)算流圖

其中，A和B的距離稱為翅尖距。這種方法和直接進(jìn)行DFT計(jì)算相比較，運(yùn)算量減少一半。按照這種分解運(yùn)算的思想，將X1(k)和X2(k)繼續(xù)向下分解，直到最后變?yōu)橐稽c(diǎn)序列，此時(shí)的運(yùn)算量大大減小。以8點(diǎn)序列X(n)的FFT運(yùn)算為例：

蝶式運(yùn)算流圖如圖3所示。比較FFT和DFT的運(yùn)算量。假設(shè)序列點(diǎn)數(shù)為N，且有N=2L,當(dāng)使用FFT進(jìn)行運(yùn)算時(shí)，計(jì)算過(guò)程中只有蝶式運(yùn)算，它的復(fù)數(shù)乘法運(yùn)算量為，復(fù)數(shù)加法運(yùn)算量為；直接進(jìn)行DFT運(yùn)算時(shí)，復(fù)數(shù)乘法和復(fù)數(shù)加法的運(yùn)算量均為N2。由此可見(jiàn)，F(xiàn)FT運(yùn)算確實(shí)大大減小了DFT的運(yùn)算量。

由圖3可以看出，若想得到順序正確的頻域序列X(k)，必須對(duì)時(shí)域序列進(jìn)行重新排序。這里先說(shuō)明比特逆序的概念：設(shè)存儲(chǔ)地址m,其二進(jìn)制數(shù)（設(shè)m=2L）為m=(m1m2m3),若有=（m3m2m1），則稱為m的L位比特逆序列。若將圖4中的輸入時(shí)域序列下標(biāo)轉(zhuǎn)換成二進(jìn)制，依次是：000,100,010,110,001,101,011,111；而原時(shí)域序列下標(biāo)轉(zhuǎn)換成二進(jìn)制依次是:000,001,010,011,100,101,110,111，將兩者對(duì)比可發(fā)現(xiàn)，輸入時(shí)域序列下標(biāo)就是原時(shí)域序列下標(biāo)的比特逆序。所以FFT蝶式運(yùn)算的第一步就是對(duì)時(shí)域序列進(jìn)行比特排序。

3 系統(tǒng)軟件實(shí)現(xiàn)以及結(jié)果分析

3.1 FFT程序?qū)崿F(xiàn)

FFT算法主要包括比特逆序、蝶式權(quán)值的計(jì)算、各級(jí)的蝶式運(yùn)算以及FFT序列的輸出。程序中涉及到復(fù)數(shù)的運(yùn)算，由于計(jì)算機(jī)無(wú)法處理復(fù)數(shù)，故而將復(fù)數(shù)拆為實(shí)部和虛部，只需計(jì)算出實(shí)部和虛部的數(shù)值即可。程序設(shè)計(jì)將蝶式權(quán)w[i]拆分為實(shí)部pr[i]和虛部pi[i]，將參與各級(jí)蝶式運(yùn)算的序列拆分為實(shí)部fr和虛部fi。復(fù)數(shù)進(jìn)行加法運(yùn)算時(shí)需將兩復(fù)數(shù)的實(shí)部與實(shí)部、虛部與虛部對(duì)應(yīng)相加即可，例如w[i]+w[i+1],所得的和的實(shí)部為pr[i]+pr[i+1]，虛部為pi[i]+pi[i+1]；而復(fù)數(shù)的乘法，相對(duì)來(lái)說(shuō)有些復(fù)雜，不再是簡(jiǎn)單的虛部和實(shí)部各自相乘,例如w[i]×w,所得的乘積的實(shí)部為pr[i]×pr+pi[i]×pi，虛部為:(pr[i]+pi[i])×(pr+pi)-( pr[i]×pr+pi[i]×pi)。該程序執(zhí)行的是128個(gè)點(diǎn)的FFT運(yùn)算，共有7級(jí)蝶式運(yùn)算，每個(gè)蝶式權(quán)都可看作是的n次方，所以在一個(gè)循環(huán)內(nèi)即可計(jì)算出所有的蝶式權(quán)。每級(jí)蝶式運(yùn)算計(jì)算出新的序列值代替該級(jí)的輸入序列值，這樣避免了開(kāi)辟新的存儲(chǔ)空間，有效地節(jié)省了存儲(chǔ)空間。

FFT算法程序設(shè)計(jì)如下所示：

Uint16 ConversionCount;

Uint16 px[128];

Uint16 pz[128];

void kfft(pr,pi,n,k,fr, fi ,l,il)

{ Uint16 n,k,l,il;

double pr[],pi[],fr[], fi [];

int it,m,is,i,j,nv,l0;

double p,q,s,vr,vi,poddr,poddi;

for (it=0; it＜=n-1; it++)

{ m=it;

is=0;

for (i=0; i＜=k-1; i++)

{ j=m/2;

is=2*is+(m-2*j);

m=j;}

fr[it]=pr[is];

fi[it]=pi[is]; //此循環(huán)為比特逆序的實(shí)現(xiàn)，pr、pi是原序列的實(shí)部與虛部

} //fr、 fi 是排完序后的序列

pr[0]=1.0; pi[0]=0.0;

p=6.283185306/(1.0*n);

pr[1]=cos(p);

if (l!=0) pi[1]=-pi[1];

for (i=2; i＜=n-1; i++)

{ p=pr[i-1]*pr[1];

q=pi[i-1]*pi[1];

s=(pr[i-1]+pi[i-1])*(pr[1]+pi[1]);

pr[i]=p-q;

pi[i]=s-p-q;} //此循環(huán)為權(quán)值運(yùn)算實(shí)現(xiàn)

for (it=0; it＜=n-2; it=it+2) //it是序列下標(biāo)，

在第一級(jí)蝶式運(yùn)算中

//相鄰的兩個(gè)運(yùn)算蝶下標(biāo)相差2

{ vr=fr[it]; vi= fi [it];

fr[it]=vr+fr[it+1]; //第一級(jí)蝶式運(yùn)算中相鄰兩個(gè)數(shù)值進(jìn)行加減運(yùn)算

fi [it]=vi+ fi [it+1];

fr[it+1]=vr-fr[it+1];

fi [it+1]=vi- fi [it+1];} //此循環(huán)為第一級(jí)蝶式運(yùn)算，該級(jí)中的蝶式權(quán)均為1

m=n/2; nv=2;

for (l0=k-2; l0＞=0; l0--) //除去第一級(jí)，還有六級(jí)蝶式運(yùn)算

{ m=m/2; nv=2*nv; // 同一級(jí)相鄰的兩個(gè)蝶式運(yùn)算中數(shù)值的下標(biāo)相差為nv

for (it=0; it＜=(m-1)*nv; it=it+nv) //it仍為數(shù)值的下標(biāo)

for(j=0;j＜=(nv/2)-1;j++)

{ p=pr[m*j]*fr[it+j+nv/2]; // nv/2為蝶式運(yùn)算的翅間距

q=pi[m*j]*f i [it+j+nv/2];

s=pr[m*j]+pi[m*j];

s=s*(fr[it+j+nv/2]+f i [it+j+nv/2]);

poddr=p-q; poddi=s-p-q;

fr[it+j+nv/2]=fr[it+j]-poddr;

fi [it+j+nv/2]=f i [it+j]-poddi;

r[it+j]=fr[it+j]+poddr;

fi [it+j]=f i [it+j]+poddi; //計(jì)算出的新序列值代替輸入序列值,節(jié)省了存儲(chǔ)空間

}}} // fr、f i 分別為計(jì)算出的FFT序列的實(shí)部和虛部

3.2 運(yùn)行結(jié)果

用FFT的計(jì)算公式對(duì)單一頻率的正余弦信號(hào)進(jìn)行變換，可以得出單一頻率的信號(hào)在頻域表現(xiàn)為在其正負(fù)頻率點(diǎn)上的兩個(gè)脈沖。根據(jù)單頻信號(hào)的這一特性可驗(yàn)證上述程序是否正確。將單一頻率的正弦波信號(hào)接入硬件系統(tǒng)，打開(kāi)圖形觀察窗口，將觀察點(diǎn)數(shù)設(shè)置為128點(diǎn)，即可觀察到如圖4所示的波形。

圖4 時(shí)域正弦波的FFT運(yùn)行結(jié)果

圖5 改善后的FFT運(yùn)行結(jié)果

圖4中上半部分波形為輸入的時(shí)域正弦波，下半部分波形為經(jīng)過(guò)變換后的頻域波形。在圖4中，由于產(chǎn)生的模擬信號(hào)伴有噪音，導(dǎo)致信號(hào)頻率不單一，故而在信號(hào)的頻域中，除了信號(hào)對(duì)應(yīng)的頻率點(diǎn)有脈沖外，其它頻率點(diǎn)處也有較弱的脈沖。

對(duì)輸入信號(hào)的頻率進(jìn)行調(diào)整，使信號(hào)的組成頻率盡可能單一，可以得到更好的波形，如圖5所示?？梢钥吹皆谛盘?hào)對(duì)應(yīng)的頻率點(diǎn)有脈沖，其它頻率點(diǎn)的振幅為0，與理論推導(dǎo)結(jié)果一致，由此可驗(yàn)證程序設(shè)計(jì)的正確性。

4 結(jié)論

FFT是聲學(xué)、圖像和信號(hào)處理等領(lǐng)域中一種重要的分析工具，文中闡述了硬件結(jié)構(gòu)圖和信號(hào)調(diào)理電路。詳細(xì)介紹了FFT算法原理，采用C語(yǔ)言編寫(xiě)程序?qū)崿F(xiàn)了FFT算法。程序運(yùn)行結(jié)果表明TMS320F2812 DSP實(shí)現(xiàn)FFT運(yùn)算速度快，精度高[4]。

[1] 賈瑋,楊錄,張艷花.基于TMS320VC5416的FFT算法的實(shí)現(xiàn)[J].山西電子技術(shù),2009,2:11-13.

[2] Rright Hert. Rapid algorithms of digital signal processing[M].Beijing,Electronic Industrial Press,2002.

[3] 胡廣書(shū).數(shù)字信號(hào)處理[M].北京:清華大學(xué)出版社,2003.

[4] 蘇奎峰,呂強(qiáng),耿慶鋒,陳圣儉.TMS320F2812原理與開(kāi)發(fā)[M].北京:電子工業(yè)出版社,2005.