基于ADSP-TS201的大點數FFT實現

穆文爭 王啟智 朱子平

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

0 引言

在雷達信號處理中，FFT作為是頻域脈壓和譜分析的基礎，是較為常用的一種變換，一般由數字信號處理器完成。在工程應用中，當需要做FFT變換時，設計師可以調用開發環境自帶的庫函數，也可以將FFT做成公用基礎模塊，供設計者調用。此類模塊一般采用匯編語言設計，針對所使用處理器的結構特征進行優化，具有較高的運算效率。但是當需要變換的FFT點數較大時，該類模塊不再適用，因為輸入數據和中間緩存都要占據很大的內存空間，而一般DSP的內存空間有限且被劃分為多個不連續的block，這就破壞了該類模塊的并行運算結構，降低了運算效率，如果采用與外部存儲器交換數據的辦法，由于與外存交換數據的速度較慢，也會降低運算效率。因此，在工程上需要做大點數FFT的場合，就要考慮如何實現的問題。

本文首先介紹了大點數FFT變換的數學原理，即將一維FFT拆成二維FFT實現的算法原理，然后以96K點為例，介紹其在通用處理器ADSP-TS201上的實現過程。

1 數學原理

設有限長序列x(n)的長度為N，則其DFT為

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

令n和k的序號映射定義如下:

則N點DFT可以表示為:

令

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

矩陣的元素為復數G(n2，k1)。因n2行的數據在計算完x(N2n1+n2)的N1點DFT后不再需要，所以G(n2，k1)可以存儲在同一行。計算X(k)，還應乘以旋轉因子Wk1n2N形成新的陣列

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

由此可見，將一維FFT拆成二維FFT的實現過程就是將N點FFT分解為N2個N1點行變換與N1個N2點列變換的級聯［2］，經過二維拆分以后，設計師就可以調用已有的FFT函數。

2 工程實現面臨的問題

除了輸入數據和中間緩存要占據很大的內存空間外，做大點數FFT遇到的另外一個問題就是存儲旋轉因子。當變換的點數很大時，在FPGA或DSP中存儲這些旋轉因子需要很大的內存空間，這在工程實現上不方便。下面說明如何對旋轉因子進行分解，以降低對內存的需求量。

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

因為k＜N，即N1·α+β＜N1·N2，又β是k對 N1的余數，有0≤β≤N1-1 ，所以 N1·α ＜ N1·N2，即α＜N2。

以1M點FFT為例，N可以分解為N1×N2=1024×1024，這樣存儲的旋轉因子個數就由1M點減少為1024+1024=2048點。

3 96K點FFT工程實現

以上對大點數FFT的實現方法做了原理性介紹，下面通過一個具體的例子，說明其在通用處理器ADSP-TS201上的實現過程。

由上述一維FFT拆成二維FFT的算法原理，可以歸納出大點數FFT實現的基本步驟［3］:

a.將N點輸入數據表示成L行C列矩陣;

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

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

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

e.將X(c+l1×C)整序為X(k)，變換完成。

圖1為大點數FFT實現流程圖。

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

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

表1是96K點FFT的N1和N2不同分解方式對應的旋轉因子存儲量，可見256×384是存儲量最小的一種分解方式。

表1 不同分解方式對應的旋轉因子存儲量

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

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

處理96K點FFT時，由于一個block最多存儲64K點復數，所以需要1.5個block作為輸入數據緩存。可以考慮將96k點復數分為三段、每段32K點進行存儲，其中一種存儲方式如圖3所示。

圖2 ADSP-TS201的內存組織

圖3 96K點輸入數據的緩存方式

下面按照前述的實現步驟，對96K點復數進行FFT變換。

第一步:

對輸入的3段數據，先對每列做3點DFT，然后對做過DFT的數據乘以中間旋轉因子。需要注意的是，每列3點數據所乘的旋轉因子均不同，旋轉因子的上標為行號×列號。乘積結果覆蓋原數據見圖4。

圖4 每列3點DFT

每列旋轉因子的上標見表2。

表2 旋轉因子上標

第二步:

數據經過第一步處理后，再對每一行做FFT，結果仍然覆蓋原數據，如圖5所示。

圖5 每行做FFT

第三步:

將結果重新排序。由于DSP沒有足夠的內存用作緩存，所以重排只能分段進行，這樣重排就需要兩輪，重排示意圖見圖6。

圖6 數據重排

第一輪排序:

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

經過第一輪重排后，數據分布規律如圖8所示。

經過第一輪重排后，每段數據的首地址如圖9所示。

可見，經過第一輪重排后，數據雖然得到了調整，但沒有達到圖6所示的順序，還要經過第二輪重排，才能完成整個排序。

第二輪排序:

第二輪排序就是對圖8的12個數據段進行交換，直至達到圖6的順序。交換時以整個數據段進行，第二輪排序的交換順序如圖10所示。

第二輪排序時，每次交換的2個數據段的首地址見表3。

圖7 第一輪重排示意圖

圖8 第一輪重排后的數據分布

圖9 第一輪重排后每段數據的首地址

圖10 第二輪排序數據交換過程

表3 第二輪排序時每次交換的2個數據段的首地址

至此，經過第二輪排序以后，數據在內存中的順序就達到了期望的順序，即完成了整個排序，從而完成了整個FFT變換。

對于大點數IFFT，可以在FFT的基礎上實現，步驟如圖11所示。

圖11 大點數IFFT實現步驟

運行時間統計:在工程實現時，96K點復數FFT模塊由四個函數組成:a.對每列數據做3點DFT并乘以中間旋轉因子的函數DFT_3point_multMidT-wid.asm;b.對每行做32K點FFT的函數bw_cvfft_c.asm，該函數為通用的參數化FFT模塊;c.第一輪排序函數ReArrange_dat_1stRound.asm;d.第二輪排序函數ReArrange_dat_2ndRound.asm。當ADSPTS201的核時鐘運行于480MHz時，每個函數模塊的運行時間統計見表4。

表4 函數運行時間統計

可見運行完96K點復數FFT需要的時間為:1.66+1.6*3+0.35+0.41=7.22ms。經測試，做完96K點頻域脈壓并對結果求模取對數，共計用時16.7ms。可見，這種處理大點數FFT的方法不僅解決了內存需求量大的問題，而且函數模塊具有較好的實時性。

4 結束語

針對雷達信號處理遇到的大時寬脈壓問題，給出了大點數FFT的變換方法。該方法大大降低了因存儲旋轉因子對內存的需求，在不與外部存儲器交換數據的情況下，能夠快速完成FFT變換。文章首先介紹了大點數FFT變換的數學原理，然后以96K點為例，闡述了其在通用處理器ADSP-TS201上的實現過程，為工程上實現其它點數的大點數FFT提供了有益的借鑒。

［1］丁玉美，高西全.數字信號處理［M］.西安:西安電子科技大學出版社，2000:68-77.

［2］蘇濤，莊德靖.大點數FFT算法的改進及其實現［J］. 現代雷達，2005，27(7):23-26.

［3］王曉君，龍騰，周希元.二維級聯流水結構大點數FFT運算器實現研究［J］.無線電工程，2010，40(11):19-22.

［4］ 劉書明，羅勇江.ADSP TS20XS系列DSP原理與應用設計［M］.北京:電子工業出版社，2007:115-118.