穆文爭(zhēng) 王啟智 朱子平

(中國(guó)電子科技集團(tuán)公司第三十八研究所 合肥 230088)

0 引言

在雷達(dá)信號(hào)處理中，F(xiàn)FT作為是頻域脈壓和譜分析的基礎(chǔ)，是較為常用的一種變換，一般由數(shù)字信號(hào)處理器完成。在工程應(yīng)用中，當(dāng)需要做FFT變換時(shí)，設(shè)計(jì)師可以調(diào)用開(kāi)發(fā)環(huán)境自帶的庫(kù)函數(shù)，也可以將FFT做成公用基礎(chǔ)模塊，供設(shè)計(jì)者調(diào)用。此類模塊一般采用匯編語(yǔ)言設(shè)計(jì)，針對(duì)所使用處理器的結(jié)構(gòu)特征進(jìn)行優(yōu)化，具有較高的運(yùn)算效率。但是當(dāng)需要變換的FFT點(diǎn)數(shù)較大時(shí)，該類模塊不再適用，因?yàn)檩斎霐?shù)據(jù)和中間緩存都要占據(jù)很大的內(nèi)存空間，而一般DSP的內(nèi)存空間有限且被劃分為多個(gè)不連續(xù)的block，這就破壞了該類模塊的并行運(yùn)算結(jié)構(gòu)，降低了運(yùn)算效率，如果采用與外部存儲(chǔ)器交換數(shù)據(jù)的辦法，由于與外存交換數(shù)據(jù)的速度較慢，也會(huì)降低運(yùn)算效率。因此，在工程上需要做大點(diǎn)數(shù)FFT的場(chǎng)合，就要考慮如何實(shí)現(xiàn)的問(wèn)題。

本文首先介紹了大點(diǎn)數(shù)FFT變換的數(shù)學(xué)原理，即將一維FFT拆成二維FFT實(shí)現(xiàn)的算法原理，然后以96K點(diǎn)為例，介紹其在通用處理器ADSP-TS201上的實(shí)現(xiàn)過(guò)程。

1 數(shù)學(xué)原理

設(shè)有限長(zhǎng)序列x(n)的長(zhǎng)度為N，則其DFT為

令N=N1N2，可以將x(n)分解為N2個(gè)長(zhǎng)度為N1的序列，將這些序列用陣列x'表示，則

令n和k的序號(hào)映射定義如下:

則N點(diǎn)DFT可以表示為:

令

則G(n2，k1)為序列x(N2n1+n2)的N1點(diǎn) DFT，即上式二維陣列N2行元素的DFT。計(jì)算陣列x'每一行N1點(diǎn)DFT就得到另一個(gè)矩陣:

矩陣的元素為復(fù)數(shù)G(n2，k1)。因n2行的數(shù)據(jù)在計(jì)算完x(N2n1+n2)的N1點(diǎn)DFT后不再需要，所以G(n2，k1)可以存儲(chǔ)在同一行。計(jì)算X(k)，還應(yīng)乘以旋轉(zhuǎn)因子Wk1n2N形成新的陣列

這樣DFT變換結(jié)果就可由二維陣列得出。

由此可見(jiàn)，將一維FFT拆成二維FFT的實(shí)現(xiàn)過(guò)程就是將N點(diǎn)FFT分解為N2個(gè)N1點(diǎn)行變換與N1個(gè)N2點(diǎn)列變換的級(jí)聯(lián)［2］，經(jīng)過(guò)二維拆分以后，設(shè)計(jì)師就可以調(diào)用已有的FFT函數(shù)。

2 工程實(shí)現(xiàn)面臨的問(wèn)題

除了輸入數(shù)據(jù)和中間緩存要占據(jù)很大的內(nèi)存空間外，做大點(diǎn)數(shù)FFT遇到的另外一個(gè)問(wèn)題就是存儲(chǔ)旋轉(zhuǎn)因子。當(dāng)變換的點(diǎn)數(shù)很大時(shí)，在FPGA或DSP中存儲(chǔ)這些旋轉(zhuǎn)因子需要很大的內(nèi)存空間，這在工程實(shí)現(xiàn)上不方便。下面說(shuō)明如何對(duì)旋轉(zhuǎn)因子進(jìn)行分解，以降低對(duì)內(nèi)存的需求量。

如果將N分解為N1×N2，k分解為N1×α+β，其中α為倍數(shù)，β為余數(shù)。則有

因?yàn)閗＜N，即N1·α+β＜N1·N2，又β是k對(duì) N1的余數(shù)，有0≤β≤N1-1 ，所以 N1·α ＜ N1·N2，即α＜N2。

以1M點(diǎn)FFT為例，N可以分解為N1×N2=1024×1024，這樣存儲(chǔ)的旋轉(zhuǎn)因子個(gè)數(shù)就由1M點(diǎn)減少為1024+1024=2048點(diǎn)。

3 96K點(diǎn)FFT工程實(shí)現(xiàn)

以上對(duì)大點(diǎn)數(shù)FFT的實(shí)現(xiàn)方法做了原理性介紹，下面通過(guò)一個(gè)具體的例子，說(shuō)明其在通用處理器ADSP-TS201上的實(shí)現(xiàn)過(guò)程。

由上述一維FFT拆成二維FFT的算法原理，可以歸納出大點(diǎn)數(shù)FFT實(shí)現(xiàn)的基本步驟［3］:

a.將N點(diǎn)輸入數(shù)據(jù)表示成L行C列矩陣;

b.對(duì)矩陣所有列分別做L點(diǎn)FFT，得到G(l1，c);

c.將G(l1，c)乘以中間旋轉(zhuǎn)因子得到GG(l1，c);

d.對(duì)GG(l1，c)所有行做C點(diǎn)FFT，得到X(c+l1×C);

e.將X(c+l1×C)整序?yàn)閄(k)，變換完成。

圖1為大點(diǎn)數(shù)FFT實(shí)現(xiàn)流程圖。

需要注意的是，數(shù)據(jù)輸入分解和因計(jì)算中間旋轉(zhuǎn)因子而進(jìn)行的分解是兩種不同的分解，它們可以不同。例如本例中，數(shù)據(jù)輸入可以分解為96k=3×32k，計(jì)算中間旋轉(zhuǎn)因子可以分解為96k=256×384。在對(duì)96k=N1×N2分解時(shí)，需要遵循一個(gè)原則，就是其中一個(gè)因子要為2的整冪次，因?yàn)檫@樣便于利用ADSP-TS201的匯編指令FEXT快速提取α和β。

可知，旋轉(zhuǎn)因子只需存儲(chǔ)一半即可，另一半可由共軛對(duì)稱性得到;在存儲(chǔ)旋轉(zhuǎn)因子時(shí)，由于 β 的范圍是0～(N1-1)，所以這N1個(gè)旋轉(zhuǎn)因子要全部存儲(chǔ)。這樣需要存儲(chǔ)的旋轉(zhuǎn)因子個(gè)數(shù)就可以由(N1+N2)點(diǎn)減少為點(diǎn)，從而進(jìn)一步減少了所需的內(nèi)存空間。有了這兩種旋轉(zhuǎn)因子，計(jì)算大點(diǎn)數(shù)FFT時(shí)任何一種旋轉(zhuǎn)因子都可以通過(guò)(9)式實(shí)時(shí)生成。

表1是96K點(diǎn)FFT的N1和N2不同分解方式對(duì)應(yīng)的旋轉(zhuǎn)因子存儲(chǔ)量，可見(jiàn)256×384是存儲(chǔ)量最小的一種分解方式。

表1 不同分解方式對(duì)應(yīng)的旋轉(zhuǎn)因子存儲(chǔ)量

由于所使用的處理器為ADSP-TS201，下面對(duì)該處理器的內(nèi)存加以說(shuō)明。ADSP-TS201具有24Mbit內(nèi)存，被劃分為 6個(gè) block，每個(gè) block為4Mbit(128K ×32bit)，組織方式如下［4］。

在采用開(kāi)發(fā)環(huán)境(VisualDSP++5.0)默認(rèn)的鏈接描述文件時(shí)，block0用來(lái)存放代碼，block4的部分空間用作堆棧，這樣可供用戶使用的完整的128K字內(nèi)存塊就只有4塊，即 block2、block6、block8和block10。

處理96K點(diǎn)FFT時(shí)，由于一個(gè)block最多存儲(chǔ)64K點(diǎn)復(fù)數(shù)，所以需要1.5個(gè)block作為輸入數(shù)據(jù)緩存。可以考慮將96k點(diǎn)復(fù)數(shù)分為三段、每段32K點(diǎn)進(jìn)行存儲(chǔ)，其中一種存儲(chǔ)方式如圖3所示。

圖2 ADSP-TS201的內(nèi)存組織

圖3 96K點(diǎn)輸入數(shù)據(jù)的緩存方式

下面按照前述的實(shí)現(xiàn)步驟，對(duì)96K點(diǎn)復(fù)數(shù)進(jìn)行FFT變換。

第一步:

對(duì)輸入的3段數(shù)據(jù)，先對(duì)每列做3點(diǎn)DFT，然后對(duì)做過(guò)DFT的數(shù)據(jù)乘以中間旋轉(zhuǎn)因子。需要注意的是，每列3點(diǎn)數(shù)據(jù)所乘的旋轉(zhuǎn)因子均不同，旋轉(zhuǎn)因子的上標(biāo)為行號(hào)×列號(hào)。乘積結(jié)果覆蓋原數(shù)據(jù)見(jiàn)圖4。

圖4 每列3點(diǎn)DFT

每列旋轉(zhuǎn)因子的上標(biāo)見(jiàn)表2。

表2 旋轉(zhuǎn)因子上標(biāo)

第二步:

數(shù)據(jù)經(jīng)過(guò)第一步處理后，再對(duì)每一行做FFT，結(jié)果仍然覆蓋原數(shù)據(jù)，如圖5所示。

圖5 每行做FFT

第三步:

將結(jié)果重新排序。由于DSP沒(méi)有足夠的內(nèi)存用作緩存，所以重排只能分段進(jìn)行，這樣重排就需要兩輪，重排示意圖見(jiàn)圖6。

圖6 數(shù)據(jù)重排

第一輪排序:

申請(qǐng)32K空間作為重排的緩存，對(duì)這32K緩存，只用到其中的24K。數(shù)據(jù)重排的程序由匯編指令設(shè)計(jì)，為了提高效率，存取數(shù)據(jù)同時(shí)進(jìn)行。先取數(shù)據(jù)至寄存器，存的時(shí)候以需要的順序從寄存器讀出，這樣將數(shù)據(jù)存至緩存的時(shí)候，就已經(jīng)重排好了。存完24K數(shù)據(jù)(每段取8K)后，再將數(shù)據(jù)搬回，將原數(shù)據(jù)覆蓋，這一過(guò)程如圖7。

經(jīng)過(guò)第一輪重排后，數(shù)據(jù)分布規(guī)律如圖8所示。

經(jīng)過(guò)第一輪重排后，每段數(shù)據(jù)的首地址如圖9所示。

可見(jiàn)，經(jīng)過(guò)第一輪重排后，數(shù)據(jù)雖然得到了調(diào)整，但沒(méi)有達(dá)到圖6所示的順序，還要經(jīng)過(guò)第二輪重排，才能完成整個(gè)排序。

第二輪排序:

第二輪排序就是對(duì)圖8的12個(gè)數(shù)據(jù)段進(jìn)行交換，直至達(dá)到圖6的順序。交換時(shí)以整個(gè)數(shù)據(jù)段進(jìn)行，第二輪排序的交換順序如圖10所示。

第二輪排序時(shí)，每次交換的2個(gè)數(shù)據(jù)段的首地址見(jiàn)表3。

圖7 第一輪重排示意圖

圖8 第一輪重排后的數(shù)據(jù)分布

圖9 第一輪重排后每段數(shù)據(jù)的首地址

圖10 第二輪排序數(shù)據(jù)交換過(guò)程

表3 第二輪排序時(shí)每次交換的2個(gè)數(shù)據(jù)段的首地址

至此，經(jīng)過(guò)第二輪排序以后，數(shù)據(jù)在內(nèi)存中的順序就達(dá)到了期望的順序，即完成了整個(gè)排序，從而完成了整個(gè)FFT變換。

對(duì)于大點(diǎn)數(shù)IFFT，可以在FFT的基礎(chǔ)上實(shí)現(xiàn)，步驟如圖11所示。

圖11 大點(diǎn)數(shù)IFFT實(shí)現(xiàn)步驟

運(yùn)行時(shí)間統(tǒng)計(jì):在工程實(shí)現(xiàn)時(shí)，96K點(diǎn)復(fù)數(shù)FFT模塊由四個(gè)函數(shù)組成:a.對(duì)每列數(shù)據(jù)做3點(diǎn)DFT并乘以中間旋轉(zhuǎn)因子的函數(shù)DFT_3point_multMidT-wid.asm;b.對(duì)每行做32K點(diǎn)FFT的函數(shù)bw_cvfft_c.asm，該函數(shù)為通用的參數(shù)化FFT模塊;c.第一輪排序函數(shù)ReArrange_dat_1stRound.asm;d.第二輪排序函數(shù)ReArrange_dat_2ndRound.asm。當(dāng)ADSPTS201的核時(shí)鐘運(yùn)行于480MHz時(shí)，每個(gè)函數(shù)模塊的運(yùn)行時(shí)間統(tǒng)計(jì)見(jiàn)表4。

表4 函數(shù)運(yùn)行時(shí)間統(tǒng)計(jì)

可見(jiàn)運(yùn)行完96K點(diǎn)復(fù)數(shù)FFT需要的時(shí)間為:1.66+1.6*3+0.35+0.41=7.22ms。經(jīng)測(cè)試，做完96K點(diǎn)頻域脈壓并對(duì)結(jié)果求模取對(duì)數(shù)，共計(jì)用時(shí)16.7ms?？梢?jiàn)，這種處理大點(diǎn)數(shù)FFT的方法不僅解決了內(nèi)存需求量大的問(wèn)題，而且函數(shù)模塊具有較好的實(shí)時(shí)性。

4 結(jié)束語(yǔ)

針對(duì)雷達(dá)信號(hào)處理遇到的大時(shí)寬脈壓?jiǎn)栴}，給出了大點(diǎn)數(shù)FFT的變換方法。該方法大大降低了因存儲(chǔ)旋轉(zhuǎn)因子對(duì)內(nèi)存的需求，在不與外部存儲(chǔ)器交換數(shù)據(jù)的情況下，能夠快速完成FFT變換。文章首先介紹了大點(diǎn)數(shù)FFT變換的數(shù)學(xué)原理，然后以96K點(diǎn)為例，闡述了其在通用處理器ADSP-TS201上的實(shí)現(xiàn)過(guò)程，為工程上實(shí)現(xiàn)其它點(diǎn)數(shù)的大點(diǎn)數(shù)FFT提供了有益的借鑒。

［1］丁玉美，高西全.數(shù)字信號(hào)處理［M］.西安:西安電子科技大學(xué)出版社，2000:68-77.

［2］蘇濤，莊德靖.大點(diǎn)數(shù)FFT算法的改進(jìn)及其實(shí)現(xiàn)［J］. 現(xiàn)代雷達(dá)，2005，27(7):23-26.

［3］王曉君，龍騰，周希元.二維級(jí)聯(lián)流水結(jié)構(gòu)大點(diǎn)數(shù)FFT運(yùn)算器實(shí)現(xiàn)研究［J］.無(wú)線電工程，2010，40(11):19-22.

［4］ 劉書(shū)明，羅勇江.ADSP TS20XS系列DSP原理與應(yīng)用設(shè)計(jì)［M］.北京:電子工業(yè)出版社，2007:115-118.