基于FPGA的二維譜峰搜索算法硬件架構設計

房海松，司偉建

（哈爾濱工程大學，信息與通信工程學院，黑龍江 哈爾濱 150001）

0 引 言

經典的MUSIC算法是陣列測向技術中所用到的DOA估計方法之一，是空間譜估計測向理論的重要基石［1－2］。當前 大 多采 用 數字 信 號 處 理 （digital signal processing，DSP）芯片來實現二維或更高維的MUSIC算法，而同DSP相比，FPGA具有內部并行運行、可根據具體要求進行設計、運行速度更快等優點［3－4］。當前FPGA主要用于實現一維MUSIC算法，相比而言，二維MUSIC算法更適合應用于實際環境的DOA估計，研究二維 MUSIC算法更具有實際意義［5－6］。

譜峰搜索的FPGA實現的主要難點是計算量比較大，消耗的硬件資源比較多，因此在減少資源消耗的同時滿足估計精度和實時性要求成為相關研究的熱點［7］。在信號處理技術的發展過程中，譜峰搜索的FPGA實現的最初應用是針對一維信號進行處理［8－9］。針對不同的成本和速度要求，文獻［10］采用一維脈動結構實現空間譜的計算，加速了譜峰搜索，在200 MHz工作頻率下耗時57μs實現了精度為0．5°的一維搜索。文獻［11］提出了基于Chirp－z變換的搜索方法，適合對一定區域進行精細的搜索，能有效提高搜索效率。文獻［12］給出了適用于均勻線陣的實數化預處理方法和實用的空間譜定義，提高了算法的實時性，但是譜峰搜索的范圍小，精度有限。

二維譜峰搜索包含更多的有用信息，更具有實用性。文獻［13］在進行二維譜峰搜索的粗搜索時使用了二維數據插值，能有效地減少運算量。文獻［14］在二維譜峰搜索中使用雞群算法進行優化，尋優能力強，能快速搜索到譜峰。文獻［15］設計了一種新的方案，能夠減少執行時間和節省FPGA內部資源，完成二維譜峰搜索需要的時鐘周期為31 500，方位角和仰角偏差均小于0．1°，估計精度高，滿足實時性要求。文獻［16］采用分步搜索法實現譜峰搜索，設計專用硬件電路實時計算方向向量以節省存儲資源，具有精度高、速度快、資源消耗少的優勢，在100 MHz工作頻率的FPGA芯片上，實現估計精度為0．1°的二維譜峰搜索耗時2．989 ms。文獻［17］通過空域劃分實現了二維譜峰的并行搜索，有效地提高了搜索效率，并且能估計空域交界處的目標，但在計算譜值時使用等效值代替平方和，降低了估計精度。

本文針對實時性要求更高的二維MUSIC算法的譜峰搜索部分，合理設計出各個步驟的實現結構來解決實現時間與硬件資源消耗的平衡問題。二維譜峰搜索涉及的計算量較大，硬件資源消耗也較大，因此在保證資源消耗和估計精度滿足要求的前提下，需使譜峰搜索的實現時間更少。本文所采用的硬件架構設計方法，方位角的搜索范圍為0°～360°，仰角的搜索范圍為30°～90°，當工作頻率為100 MHz、估計精度為1°時，二維譜峰搜索實現時間僅需幾微秒，具有很高的實時性。

1 譜峰搜索原理與實現步驟

本文采用8陣元均勻分布的平面圓陣，陣元分布如圖1所示。

圖1 陣元分布Fig．1 The distribution of array elements

采用分步搜索法來進行搜索，具體操作為：首先采用較大的步長（3°）在整個搜索范圍進行粗搜索，確定譜峰的大致位置；然后在粗搜索的搜索結果附近采用較小的步長（1°）進行精搜索，得到準確的信號入射角（方位角θ和仰角φ）。這樣可以在保證精度的前提下，大幅度降低計算量［18］。

譜峰搜索包括粗搜索的陣列流型矢量模塊、譜值模塊和譜峰粗搜索模塊以及精搜索的陣列流型矢量模塊、譜值模塊和譜峰精搜索模塊。譜峰搜索的流程框圖如圖2所示。

圖2 譜峰搜索的流程框圖Fig．2 The flowchart of spectrum peak search

其實現步驟如下：

1）第一步，選定一組方位角θ和仰角φ，根據8個陣元的坐標（x，y），可得到陣列流型矢量a（θ，φ）（1×8矩陣），如式（1）所示

其中，λ為入射信號的波長。

2）第二步，將陣列流型矢量a（θ，φ）和噪聲子空間UN（8×8矩陣，整體輸入）按照式（2）進行內積累加運算，形成譜值。

3）第三步，對所有的譜值進行粗搜索。

4）第四步，對粗搜索選取的譜峰及其周圍的譜值進行精搜索。精搜索同樣需要經過第一步和第二步，然后再對譜值進行精搜索。

2 譜峰搜索模塊分析

2．1 粗搜索

2．1．1 陣列流型矢量模塊

此模塊的目的是計算式（1），選定一組方位角和仰角，得到8×1矩陣，分別對應8個陣元。記新矩陣X，使其滿足式（3）；并令新矩陣Y滿足式（4），可得

為了減小計算量，需要進行預處理，完成復數域到實數域的轉換，且不改變最終結果［20－21］。因此，對于第1、2、3、4陣元，有

根據仰角取值30°、33°、…、87°、89°，可將陣列流型矢量模塊分成21個小模塊，對每個小模塊制定查找表數據X。在每個小模塊中，根據方位角取值0°、3°、…、357°將每個小模塊分成120部分，每一部分均實施如圖3所示的5個步驟。按照流水線結構實現這120部分，5個步驟根據時鐘依次完成，從而有效提高了資源利用率，降低了資源消耗。

求取X時共制定21個并行的查找表，分別對應仰角30°、33°、…、87°、89°；每個查找表都有120個地址，分別對應方位角0°、3°、…、357°。值得注意的是，21個查找表并行消耗的是硬件資源，依次讀取120個地址中的數據消耗的是時鐘，這是經過權衡時鐘消耗和硬件資源消耗得出的結果。求取正弦值和余弦值時為了減少查找表的資源消耗，需將第1、2、3、4陣元對應的余弦值轉化為正弦值，并且將其他象限的正弦值均轉化為第一象限的正弦值。由此可使該查找表的資源消耗降低到原來的八分之一，且不會增加誤差。因此，制定84個并行的查找表即可。

圖3 陣列流型矢量模塊實現步驟Fig．3 The realization steps of array flow pattern vector module

2．1．2 譜值模塊

譜值模塊本質是矩陣計算，在計算式（2）時，可先計算式（8），再計算式（9），得到譜值。其中新矩陣b為1×8矩陣。

同陣列流型矢量模塊類似，譜值模塊可分成22個小模塊（仰角90°占用一個小模塊），每個小模塊同樣采用流水線結構來實現譜值的計算。

現以仰角30°為示例說明對應譜值的計算過程。對于方位角0°、3°、…、357°來說，矩陣相乘中的運算均可分成4個步驟，如圖4所示。步驟1、2與步驟3、4分別按照兩條獨自進行的流水線結構實現，如表1所示，采用流水線結構極大地減少了乘法器資源的消耗。仰角30°、33°、…、89°、90°（共22個）對應的譜值計算需并行進行。

圖4 矩陣相乘的運算步驟Fig．4 The operational steps of matrix multiplication

表1 矩陣相乘的運算步驟的流水線結構Tab．1 The pipelined structure of the operational steps of matrix multiplication

2．1．3 譜峰粗搜索模塊

由于仰角對應的譜值計算是并行進行的，所以每個時鐘周期可得到22個并行的譜值。隨著時鐘周期的增加，依次可得到方位角0°、3°、…、357°（共120個）對應的譜值，共計2 640（22×120）個譜值。這些譜值按照圖5的方式進行排列。分別以這些譜值為中心，與其周圍的8個譜值相比較，判斷這些中心譜值是否為極小值并記錄下來。由于同一時鐘兩個相鄰的仰角對應的譜值不可能同時為極小值，因此可將這些譜值分成11組，如表2所示。這樣一來，譜峰粗搜索模塊就可分為11個小模塊。

將這11組得到的極小值（記為偽譜峰）分別按照從小到大的順序排列，并分別選取前4個偽譜峰，共計得到44個偽譜峰。另外，由于仰角90°對應的譜值均相同，只需選取1個即可，因此譜峰粗搜索模塊得到45個偽譜峰，然后根據這些偽譜峰進行精搜索。

圖5 譜值的排列方式Fig．5 The arrangement of spectral values

表2 譜峰粗搜索中譜值的分組情況Tab．2 The grouping of spectral values in spectrum peak coarse search

而在邊緣搜索中，對于仰角30°，可與中心譜值相比較的譜值只有5個；對于仰角90°，在同一方位角的條件下，對應的譜值與仰角89°對應的譜值相比較。

2．2 精搜索

精搜索以1°為搜索步長，以粗搜索得到的44個偽譜峰（記為中心譜值，仰角90°對應的偽譜峰暫不考慮）為中心，重新計算其周圍3×3區域的譜值（記為區域譜值）。每個區域有9個譜值，找到其中的最小值，共計44個，再將這44個最小值和仰角90°對應的偽譜峰按照從小到大的順序排列，并選取前6個，即為最終譜峰。

為了保持時鐘消耗和資源消耗的平衡，需合理分配這44個偽譜峰的先后次序。將粗搜索得到的第1、2組偽譜峰隨時鐘以串行的方式先后進行精搜索，同理，其他偽譜峰的分組情況如表3所示。因此需要對這6組偽譜峰進行并行操作，精搜索的執行過程如圖6所示。

表3 偽譜峰的分組情況Tab．3 The grouping of pseudo spectrum peaks

圖6 精搜索的執行過程Fig．6 The execution of fine search

2．2．1 陣列流型矢量模塊

根據粗搜索得到的6組偽譜峰，可將陣列流型矢量模塊分成6個小模塊，每個小模塊和粗搜索的陣列流型矢量模塊類似，也分為5個步驟，不同的是前者需根據粗搜索得到的偽譜峰對應的方位角和仰角確定譜值區域（3×3區域）。確定譜值區域需結合求取X的查找表實現，根據中心譜值的仰角與方位角可確定區域譜值的仰角與方位角，從而確定區域譜值對應的查找表及其地址。精搜索和粗搜索在陣列流型矢量模塊求取sin值時復用同樣的查找表，節省了大量的查找表資源。

2．2．2 譜值模塊

本模塊可分為6個小模塊，和陣列流型矢量的6個小模塊相對應，分別對應A、B、C、D、E、F組偽譜峰，且并行計算譜值。和粗搜索的譜值模塊相同，精搜索的譜值模塊采用相同的流水線結構來實現譜值的計算。

2．2．3 譜峰精搜索模塊

譜峰精搜索模塊只有一個模塊。以A組的第一個偽譜峰為例，其對應的區域譜值有9個，在這9個譜值中找到最小值。同理，分別在B、C、D、E、F組的第一個偽譜峰對應的區域譜值中找到最小值。然后將這6個最小值按照從小到大的順序排列。隨著時鐘的進行，分別在A、B、C、D、E、F組的第二個偽譜峰對應的區域譜值中找到最小值，并將這些最小值和之前得到的6個最小值按照從小到大的順序排列。如此，直到把所有偽譜峰對應的區域譜值的最小值進行排序。最后，將粗搜索得到的仰角90°對應的偽譜值與這些最小值進行排序，得到排序后的6個譜值，即為最終結果。

3 仿真驗證與性能分析

3．1 FPGA實現

為了對本文提出的二維譜峰搜索設計方案進行性能測試，首先采用Verilog語言，在Vivado工具中進行綜合仿真，然后將設計下載到Xilinx公司Virtex－7系列FPGA的XC7V690T開發板中進行實際驗證。完成所有模塊需要377個時鐘周期，估計精度為1°，當工作頻率為100 MHz時，實現時間為3．77μs。硬件資源消耗如表4所示。

表4 二維譜峰搜索硬件資源消耗Tab．4 The hardware resource consumption of two－dimensional spectrum peak search

3．2 仿真驗證

采用8陣元均勻圓陣的陣元模型，入射信號的頻率為3 GHz，采樣快拍數為128，要求估計精度為1°且實現時間為5μs以內。信噪比分別為4～15 dB，在100 MHz工作頻率的FPGA芯片上完成了12組數據仿真實驗，每個信噪比下分別進行了500次實驗，得到方位角和仰角，在單信號源和雙信號源條件下其均方根誤差分別如圖7和圖8所示?？梢钥闯?，單信號源和雙信號源對應的方位角和仰角均滿足1°的精度要求。隨著信噪比的增加，誤差越來越小。

圖7 單信號源的仿真結果Fig．7 The simulation results of single signal source

圖8 雙信號源的仿真結果Fig．8 The simulation results of dual signal source

3．3 性能分析

對于二維譜峰搜索，實現時間是衡量其性能的重要依據。本文與其他二維譜峰搜索的實現方案性能對比如表5所示。從表5中可以看出，本文譜峰搜索的速度明顯高于其他文獻給出方案的搜索速度，實時性更高。本文方案與文獻［15］、文獻［16］的方案都采用分步搜索法進行搜索，都采用查找表計算陣列流型矢量。但是同本文方案相比，文獻［15］、文獻［16］的方案流水化程度較低，導致硬件資源的消耗和算法的實現時間大幅度增加。而且，本文方案優于其他文獻方案之處還在于增加信號源不會增加譜峰搜索的實現時間。

表5 二維譜峰搜索方案性能對比Tab．5 The performance comparison of two－dimensional spectrum peak search

4 結束語

本文針對二維譜峰搜索的實時性進行了深入研究，提出了一種基于FPGA的二維譜峰搜索實現方案。該方案采用分步搜索法實現譜峰搜索，合理設計了譜峰搜索各個步驟的實現結構，大幅度減少了計算量和計算復雜度，從而在滿足估計精度的同時解決了實現時間與硬件資源消耗的平衡問題。陣列流型矢量計算、譜值計算和譜峰粗（精）搜索3個模塊完全流水化，模塊交叉重疊進行，用硬件資源交換實現時間，并且通過模塊復用降低資源消耗，從而大大提高了計算效率，極大地縮短了整體時間。與現有其他方案相比，該二維譜峰搜索實現結構的主要優點在于速度快，實時性高，具有較高的工程應用價值。