FPGA實時多普勒中心頻率估計

甘翼，李向陽

(中國電子科技集團公司第十研究所，四川 成都 610036)

0 引言

合成孔徑雷達(synthetic aperture radar,SAR)受氣象條件的影響較小，具有全天候工作的能力，在測繪、軍事等領域具有廣范的應用[1-2]。彈載SAR能夠實現對目標的快速高精度成像和持續跟蹤，是導彈實現遠距離精確打擊的重要傳感器之一[3]。彈載SAR相對其他平臺體積小、質量輕、功耗低，并且對處理實時性和大場景處理能力提出了更高的要求。因此，在彈載SAR實現方面，基于純FPGA的SAR成像實現方法，可以提高大場景處理能力、減小成像時間、提高系統集成度、降低系統功耗[4]，在并行處理能力和成像時間方面相對DSP處理器具備極大的優勢[5]。

SAR成像算法包括距離向壓縮、距離走動矯正、方位向壓縮以及多普勒中心頻率估計、多普勒調頻率估計等模塊。其中距離向壓縮、距離走動矯正和方位向壓縮主要完成FFT變換、點乘和IFFT變換，具有單一的數據處理結構，無迭代操作，因而適合采用FPGA(field-programmable gate array)進行流水化處理。多普勒中心頻率(frequency of Doppler centroid,FDC)估計、多普勒調頻率(frequency of Doppler rate,FDR)估計[6]，因為算法相對較為復雜，使用FPGA實現具有較高的難度[7]，該部分通常采用DSP處理器完成[8]。這種折中方案會引入的數據交換延遲和FPGA處理等待延遲，大幅增加系統成像處理時間。以常用的DSP實現多普勒參數估計流程為例：首先，為了在DSP上實現參數估計，需要FPGA向DSP通過SRIO發送脈壓數據，而距離向4096、方位向2048的IQ單精度浮點數據高達64 MB，通過2Lane@5Gbps的SRIO單次傳輸需要64 ms；其次，在DSP進行FDC估計時，因為后續距離走動校正需要多普勒中心頻率估計結果，所以FPGA需要暫停當前流程，等待DSP的FDC處理結果才能進行后續操作，通常會產生高達80 ms的處理等待延遲；因此，基于FPGA+DSP的SAR處理流程在彈載條件硬件資源受限的情況下，無法滿足大場景實時成像的要求，只能通過減小成像區域，來滿足系統對實時性的要求。而通過將FDC估計采用FPGA實現即可大幅降低系統處理時間，提高系統的實時性。文獻[9]提出了一種基于Zynq系列FPGA和RD算法的SAR成像實現方法，其中運動參數估計、多普勒參數估計、多普勒調頻率估計等放在ARM處理器上進行處理，而距離向壓縮、方位向壓縮以及距離走動校正則是放在FPGA上進行。文獻[10-12]提出了3種基于FPGA的SAR實現方法，其中文獻[10]沒有對參數估計部分進行闡述，文獻[11-12]分別提出了基于FPGA和CS算法的SAR成像實現方法，并且采用FPGA實現多普勒參數估計。但是在多普勒參數估計部分均采用一次曲線進行多普勒中心頻率擬合。而在彈載大斜視系統中，為了提高多普勒中心頻率估計的準確性，通常采用2次或者3次曲線擬合進行逼近[13]。

本文提出一種基于FPGA的FDC估計方法，通過調整FDC數據處理流程、實現了SAR信號處理機的“零等待”。在FDC估計算法實現過程中，提出了基于FPGA的流水線式相位解纏繞方法和基于最小二乘的參數曲線擬合算法，實現相位解纏繞的連續流水輸出以及中心頻率估計參數精確插值。

1 多普勒中心頻率估計算法

多普勒中心頻率的估計有2種途徑：①通過慣導獲得[14]；②通過算法從回波數據中提取。由于慣導本身存在測量誤差，為了獲得更好的聚焦性能，通常從回波數據中通過算法提取多普勒中心頻率參數15]。利用回波數據對多普勒中心頻率進行估計有多種方法[16]。其中，能量均衡法、匹配相關法和最大似然法需要把距離向功率譜與某一參考函數做相關，把相關函數的零點作為多普勒中心[17]；相關函數法和數據序列符號自相關法則利用回波數據中的相位信息通過自相關提取多普勒中心，不需要額外的參考函數[18]。基于FPGA的可實現性，本文采用相位函數法實現多普勒中心估計[19]。由維納-辛欽定理可知，相關函數r(τ)可以表示為其功率譜函數E{S(f)}的IDFT變換。

r(τ)=F-1{E[S(f)]}.

(1)

對于雷達發射系統，假設距離向功率譜函數E{S(f)}為天線功率方向圖W(f)，則

r(τ)=F-1{W(f-fdc)}=r0(τ)exp{j2πfdcτ}.

(2)

從式(2)中可以看出自相關函數r(τ)是多普勒中心fdc的函數，當采用與k個方位向上的數據進行自相關時，也即τ=kT，其中T=1/PRF為重頻周期，fdc可以表示為

(3)

不同方位向上回波S(nT)的自相關函數可以表示為

r(kT)=E{s*(nT)s(nT+kT)}∝

(4)

從式(3)和式(4)通過不同方位向上相同距離門的回波函數相關，即可得到fdc。因為非明顯場景中無法提取正確的多普勒信息，需根據回波能量對距離門進行篩選，截取有效的場景數據進行多普勒相關估計。同時為了降低計算量，在進行方位向相關之前，會對距離門上的數據進行抽取，最終只能得到部分多普勒頻率估計結果。為了得到完整的多普勒中心頻率，提高多普勒中心頻率估計的正確性，需要對多普勒頻率估計結果進行二次曲線擬合和插值處理，通過求均值最終得到多普勒中心頻率，用于距離走動校正。

由于FPGA內部存儲資源的限制，無法保存所有的距離向數據，因而本文采用流水線的方位向相關累積方法，無需保存距離向數據，大幅降低對存儲資源的消耗。在基于最小二乘的二次曲線擬合中，需要進行矩陣求逆操作，本文采用矩陣分解的形式，在FPGA上直接實現3×3矩陣的求逆操作。下面將詳細介紹基于FPGA的多普勒中心估計處理流程和方法。

2 “零等待”的多普勒中心估計處理流程

多普勒中心頻率估計包含距離向積累抽取、方位向相關積累、相位解纏繞和插值擬合等。FPGA的處理流程和DSP的處理流程[20]對比如圖1所示。在DSP中距離向脈壓塊數據的讀寫和多普勒中心頻率估計需順序進行，相比FPGA實現，需等待一個孔徑或子孔徑的距離向脈壓數據全部存儲后，再完整讀取進行多普勒參數估計，隨后還須增加一次DDR數據讀取時間，才能進行RCMC處理。而RCMC依賴多普勒中心頻率的估計結果。而此時多普勒中心頻率估計尚未完成，因此，一個成像孔徑或子孔徑的距離向脈壓數據需要緩存。以孔徑大小為4 096×2 048個處理單元為例，距離向脈壓后的單精度IQ數據為64 MB，由于數據量較大，FPGA內部BRAM資源無法滿足要求，需要轉儲到外部大容量DDR中，與前述DSP實現方法不同，FPGA可以發揮其并行處理的優勢，實現邊存儲邊處理。如圖1所示，在對脈壓結果進行存儲的過程中，即可同步進行多普勒中心頻率的估計。距離向壓縮結果存儲完成時，即可同步完成多普勒中心頻率估計，隨即開展RCMC處理，實現系統數據處理的“零等待”，降低了系統的處理時間，同時也簡化了系統設計。

圖1 FPGA與DSP對比處理流程Fig.1 Comparison of DSP and FPGA processes

如圖2所示，基于FPGA的多普勒中心頻率估計可以采用流水線的方式進行數據處理。傳統非流水線處理方式需要的時間為t=N(tddr+tfdc)，而流水線的處理方式其時間消耗為t=Ntddr+tfdc，其中tddr為DDR讀取一個距離向數據所需要的時間，tfdc為處理一個距離向所需要的時間，N為CPI周期內方位向積累點數。在相同tddr和tfdc的情況下，通過流水線的處理方式，可以“隱藏”N-1次FDC處理時間，對于4 096個距離向脈壓數據，在數據接收完成后，只需額外的一個處理周期即可外出多普勒中心頻率的估計，亦即“隱藏”了4 095次FDC處理時間。

圖2 基于FPGA流水與非流水對比處理流程Fig.2 Comparison between FPGA-based pipelined and non-pipelined processes

以上是基于FPGA的多普勒中心估計采用的基本流程和優化方法，下面將從距離向積累抽取、相位解纏繞、二次曲線擬合和插值3個方面詳細介紹實現過程。

2.1 距離向數據的積累抽取

為了提高多普勒中心頻率估計的魯棒性，降低系統噪聲和場景突變對參數估計的影響，采用鄰域積累實現距離向數據的平滑濾波。雖然會降低多普勒中心頻率估計的精度，但是在彈載SAR面臨的復雜電磁環境中可以提高算法的魯棒性。同時為了降低方位向相關的數據處理量，對積累后的數據進行1/4抽取。如圖3所示，在FPGA的實現過程中，采用定點累加的方式完成四鄰域數據的積累。通過2 bit計數器對數據進行循環處理，當循環計數值為‘0’時，使用當前距離門數據對累加器進行初始化，其他情況下進行距離門數據累加，在累加完成4組數據后取走處理結果。雖然理論上應該采用均值來表示累積的結果，但是每個距離門上的鄰域累積具有相同的長度，所以直接采用累加結果表示平滑濾波結果，從而可以減少邏輯資源的消耗。

圖3 距離向平滑抽取Fig.3 Smooth extraction of range dimension

根據FPGA數據處理的特點，以及距離向脈沖壓縮的實現過程，脈壓后的數據與回波數據一樣是按照固定重復周期采用流水線的方式依次接收和處理。在多普勒中心頻率估計中，距離向累積抽取可以直接在數據流上直接完成。而方位向相關累積數據則需要緩存。

數據緩存在FPGA中可以通過LE(logic elements)資源實現，也可以通過內部BRAM實現。LE寄存器具有操作簡單、讀取延時低等特點。但是對于4 096點的64 bit復數據，其資源消耗是巨大的，會占用xc7v690t芯片20%的LE資源。通過BRAM進行緩存可以使用FPGA內部的專用存儲資源，存儲容量相對較多，但是BRAM會產生1～2個時鐘周期的讀寫延遲。為了實現數據的流水線處理，需要進行預讀取。對于距離向脈壓數據，在4鄰域平滑累加的4個時鐘周期內，需要完成一次讀取和一次寫入操作。對于采用的2時鐘延時的BRAM，數據讀取操作需要在4鄰域求和計算輸出完成2個時鐘周期之前進行，寫入操作則在求和計算結果之后，為了避免讀寫操作的沖突，通常安排在求和結束后的第2個時鐘周期寫入。而對于方位向相關積累來說，則在5個時鐘周期內需完成讀取、累加和寫入操作。其處理時間為3個時鐘周期，預讀取需要2個時鐘周期。

如圖4所示，方位向相關積累需要在CPI內的所有PRT的每個距離門上進行，抽取后的每個距離門需要獨立的寄存器用于緩存方位向相關積累結果，對于4 096的距離向長度，需要2組1 024點深度的寄存器組用于分別存放IQ數據的相關累加結果。

圖4 方位向相關積累Fig.4 Correlation accumulation of azimuth range

方位向相關累積完成后，對存放在BRAM中的各個距離門累加結果依次讀取，作為cordic模塊的輸入，進行arctan操作。cordic生成的只是相關累積結果的瞬時相位信息，范圍在[-π,π)之間，因而瞬時相位信息和多普勒頻率相位之間需要進行解纏繞操作。

2.2 解纏繞

多普勒中心頻率反映了脈沖發射單元與目標之間的徑向速度差。不同PRT之間的相關性反映了不同距離門上的多普勒變化率，這是一種相位累積效應，這種累積結果通過反正切變換為角度信息時，映射在[-π,π)的范圍內。解纏繞操作可以恢復相位對應于角度上的累加信息。當角度信息在-π和-π附近變化時，通過判斷鄰域數據的變化率，增加或者減去2π，來達到相位連續累加的目的。在得到多普勒相位值后，根據其對應圓周角度的比例，通過PRT換算即可得到多普勒頻率。由于PRT的離散數據采集造成的多普勒模糊[21]，則通過慣導信息來提取多普勒模糊數進行修正，從而可以得到實際的多普勒頻率值。解纏繞的FPGA實現偽代碼如表1所示。

表1 解纏繞FPGA實現偽代碼Table 1 Unwrapping pseudo code

2.3 基于最小二乘的曲線擬合

如前所述，為了降低計算量，多普勒中心頻率估計是1/4抽取后進行處理的，但是在方位向壓縮中，需要計算每個距離門上的多普勒頻率信息，為了得到完整的多普勒頻率估計值，需要對多普勒頻率的估計值進行插值。本文采用了二次曲線擬合的方法來實現缺失部分多普勒中心頻率的插值。基于最小二乘法的曲線擬合方法提供基于均方誤差最小的最優解，可以有效地降低多普勒中心頻率估計誤差，提高估計精度。也能在一定程度上抑制由于平臺抖動等帶來的多普勒中心頻率劇烈變化，對多普勒中心頻率估計結果具有濾波的作用。多普勒中心頻率在標準掃描SAR場景中可以用二次曲線對多普勒中心頻率進行近似描述。

下面詳細介紹基于最小二乘的二次曲線擬合方法中矩陣求逆的FPGA簡化實現方法[22]。對于線性函數y=Ax，若A,y為已知量，x為待解量，則x可以表示為：x=inv(ATA)ATy。其中，由于矩陣A可能不是方陣，通過左乘AT得到方陣ATA。

對于二次曲線：

(5)

其矩陣A：

(6)

x=(a,b,c).

對于3×3矩陣，矩陣求逆操作可以通過代數余子式來實現[23]：

(7)

定義矩陣M=ATA，亦即

(8)

(9)

|A|=a11|A11|+a12|A12|+a13|A13|，

(10)

(11)

式中：

(12)

在FPGA實現過程中，矩陣M的行列式可以通過其代數余子式進行簡化實現：

|M|=a11·detA11+a12·detA12+a13·detA13=

fx4·detA11+fx3·detA12+fx2·detA13.

(13)

在得到矩陣M的代數余子式和行列式后，即可以通過cordic模塊實現除法操作，實現矩陣的求逆，從而可以求解方陣ATA的逆矩陣，從而可以實現二次曲線系數(a，b，c)。再得到二次曲線的系數后，可以通過動態生成x及其相關x2，從而完成二次曲線的插值。

綜上所述，在FPGA實現基于最小二乘的二次曲線擬合中起主要關鍵點在于：①動態生成采用流水線方式生成x通過乘法器的流水線操作從x中得到x2，隨后通過流水線方式直接輸出fx,fx2,fx3,fx4。②通過代數余子式的形式求解方陣ATA的逆矩陣inv(ATA)，從而為使用FPGA實現最小二乘算法提供有力保證。

3 實驗結果

基于FPGA的多普勒中心頻率估計為彈載SAR純FPGA實現的一部分，已經在某項目上實現，并完成了各項測試，所有指標均達到或優于相關要求。SAR成像模塊采用一片Xilinx公司的Virtex 7系列xc7vx690t芯片外擴2片DDR3芯片，完成多普勒中心頻率估計、RCMC、多普勒調頻率估計以及方位向壓縮，生成4 096×2 048大小的圖像。系統信號處理模塊主頻為200 MHz，系統LE資源占用約占芯片總資源的58%。信號處理板數量由原理樣機的4塊減少為2塊，處理功耗下降85 W以上。FDC模塊資源使用情況如表2所示。其中arcTan函數為采用cordic完成相關數據的反正切計算；CorrComplex函數計算復數據相關；corrMemRe函數為使用LUT資源例化相關累積。lsqFDC函數使用最小二乘法實現二次曲線擬合系數估計。genFDC函數使用lsqFDC計算的系數，實時生成2 048點有效場景的FDC值；thetaWarp函數完成相位角解纏繞。FDC估計模塊在RC數據接收過程中直接流水完成，不消耗額外的運行時間，接收完數據可以直接得到相關累積結果。為了使參數擬合和插值模塊不產生等待，FDC相關累積模塊拋棄最后一個方位向累積，直接轉入多普勒參數擬合狀態，因基于最小二乘的二次曲線擬合和插值模塊以及后續的多普勒模糊數疊加模塊總處理時間為4 124個時鐘周期，在200 MHz處理時鐘下，為2.1 μs，遠小于重頻周期。因而可以在RC數據接收完最后一包距離向脈壓數據前直接輸出FDC估計結果，直接轉入RCMC處理，實現系統“0等待”，因此FDC估計模塊系統消耗時間為0。

表2 FDC模塊資源使用統計Table 2 FDC module resource usage statistics

4 結束語

本文根據彈載SAR對信號處理軟硬件實現的需求，提出了一種基于FPGA的快速多普勒中心頻率估計方法，采用流水線并行處理的模式，實現FDC估計時間的零消耗。為了實現基于FPGA的FDC估計，創新性地提出了基于FPGA的流水線式相位解纏繞方法和基于最小二乘的參數曲線擬合算法，實現相位解纏繞的連續流水輸出以及中心頻率估計的二次取下擬合和插值，并在實際工程中進行了有效驗證。測試結果說明，該方法充分發揮了FPGA并行處理的優勢，大大降低了SAR信號處理對硬件資源的需求，在彈載SAR和無人機載SAR等對結構、質量和功耗限制較嚴格的條件下，有著非常廣闊的應用前景。