基于GPU的圓跡視頻SAR實(shí)時(shí)成像算法

22

(1.南京航空航天大學(xué)電子信息工程學(xué)院， 江蘇南京 211106；2.南京航空航天大學(xué)雷達(dá)成像與微波光子技術(shù)教育部重點(diǎn)實(shí)驗(yàn)室， 江蘇南京 211106)

0 引言

視頻合成孔徑雷達(dá)(視頻SAR)[1-2]以動(dòng)態(tài)的方式再現(xiàn)場景信息，有效擴(kuò)展了時(shí)間維的信息，能夠?qū)狳c(diǎn)區(qū)域持續(xù)監(jiān)測，并能直觀地反映出運(yùn)動(dòng)目標(biāo)的位置、速度及運(yùn)動(dòng)趨勢等，具有重要的軍事和民用價(jià)值。視頻SAR系統(tǒng)要求高分辨率和實(shí)時(shí)成像，這對成像算法和硬件平臺(tái)提出很高的要求。為了同時(shí)滿足這兩項(xiàng)要求，一方面要設(shè)計(jì)一個(gè)高效精確的成像算法，另一方面還要通過高效的硬件處理平臺(tái)對大量回波數(shù)據(jù)加速處理，實(shí)現(xiàn)實(shí)時(shí)成像。

在保證高分辨率成像的前提下達(dá)到視頻SAR成像的幀率要求是研究的關(guān)鍵。目前，視頻SAR成像主要采用后向投影(BP)算法，BP算法是一種時(shí)域算法，需要逐點(diǎn)遍歷得到精確聚焦的圖像，計(jì)算量相當(dāng)龐大。文獻(xiàn)[3-4]通過加大相鄰幀間數(shù)據(jù)重疊率，減少每幀圖像需要加入的新脈沖分量，解決BP計(jì)算量大的問題，以達(dá)到實(shí)時(shí)成像的要求。該方法依賴于相鄰幀間數(shù)據(jù)重疊率，適用性受到一定的限制。文獻(xiàn)[5]針對圓跡視頻SAR的特殊回波模型提出了一種改進(jìn)的去斜Chirp Scaling算法。該方法處理速度雖然不依賴于幀間數(shù)據(jù)重疊率，但利用中央處理器(CPU)處理成像速度較慢。隨著硬件處理器的發(fā)展，基于硬件平臺(tái)的實(shí)時(shí)成像處理算法日趨成熟。當(dāng)前主流的硬件平臺(tái)有：DSP、FPGA和GPU。文獻(xiàn)[6-7]利用DSP分別實(shí)現(xiàn)BP算法和距離多普勒(RD)算法用于SAR的實(shí)時(shí)成像，但DSP始終受到串行指令流的限制使其并不真正適用于實(shí)時(shí)成像。文獻(xiàn)[8-9]在現(xiàn)場可編程門陣列(FPGA)上實(shí)現(xiàn)PFA以滿足SAR的實(shí)時(shí)處理需求。該方法相比CPU處理速度有所提升，但是依然無法滿足視頻SAR實(shí)時(shí)成像的要求，并且FPGA時(shí)序難規(guī)劃，處理程序調(diào)試周期長，算法實(shí)現(xiàn)困難。文獻(xiàn)[10-12]利用可編程圖形處理器(GPU)實(shí)現(xiàn)BP算法加大成像處理速率，該方法雖然利用GPU得到了很高的加速比，但仍未滿足視頻SAR成像的幀率要求。GPU已經(jīng)由以前的專用圖形處理器演化成高并行度、多線程、擁有強(qiáng)大計(jì)算能力和極高存儲(chǔ)器帶寬的多核處理器，在解決計(jì)算密集型問題上具有很高的性價(jià)比[13-14]，并且基于CUDA C的編程起點(diǎn)低，非常適合復(fù)雜算法的快速實(shí)現(xiàn)。

基于上述問題，本文提出了一種基于GPU的圓跡視頻SAR實(shí)時(shí)成像算法。首先建立了圓跡視頻SAR的回波模型，并分析了視頻幀率與回波重疊率的關(guān)系。其次在GPU上實(shí)現(xiàn)了基于Chirp Scaling操作的PFA。最后對實(shí)測數(shù)據(jù)做并行處理成像，驗(yàn)證了本文所提算法的正確性和高效性，所以該算法可以適用于視頻SAR實(shí)時(shí)成像。

1 信號(hào)成像處理算法

1.1 圓跡視頻SAR模式

視頻SAR工作于圓跡模式，其成像幾何如圖1所示。機(jī)載雷達(dá)在高度為H的平面里以固定區(qū)域中心為圓心作圓周運(yùn)動(dòng)，飛行中波束始終照射場景中心O；雷達(dá)觀測角θ的值可達(dá)到360°；雷達(dá)俯仰角φ0和雷達(dá)到場景中心的距離R0均為常數(shù)，在整個(gè)運(yùn)動(dòng)過程中不變。傳統(tǒng)的圓跡模式是通過增大觀測角度來提高方位向分辨率，而視頻SAR是通過對圓跡模式下回波數(shù)據(jù)的合理分割實(shí)現(xiàn)視頻成像。

圖1 圓跡視頻SAR幾何模型

1.2 圓跡視頻SAR回波數(shù)據(jù)分割

為了保證視頻SAR圖像顯示的流暢性，視頻SAR圖像幀率應(yīng)不小于5 Hz[15]。由文獻(xiàn)[15]的分析可知，視頻SAR幀率與相鄰幀數(shù)據(jù)重疊率的關(guān)系滿足：

(1)

式中，ρa(bǔ)為方位向分辨率，V為雷達(dá)平臺(tái)移動(dòng)速度，R為視頻SAR到目標(biāo)區(qū)域的最短斜距，c為光速，fc為載波頻率，M為相鄰幀數(shù)據(jù)重疊率。因此，在給定幀率和載頻的前提下，對圓跡模式下回波數(shù)據(jù)進(jìn)行分割時(shí)，相鄰幀之間的數(shù)據(jù)就不可避免地要有重疊，相鄰兩幀之間的重疊關(guān)系如圖2所示。

圖2 相鄰幀之間重疊關(guān)系示意圖

1.3 基于Chirp Scaling操作的PFA

針對圓跡視頻SAR的成像特點(diǎn)，在回波數(shù)據(jù)劃分為若干幀后，用傳統(tǒng)的直線聚束SAR對其進(jìn)行近似，并對每一幀數(shù)據(jù)用基于Chirp Scaling操作的PFA成像?；贑hirp Scaling操作的PFA流程如圖3所示。

圖3 基于Chirp Scaling的PFA處理流程

根據(jù)直線聚束SAR的成像幾何，解調(diào)后去斜信號(hào)的回波表達(dá)式為

(2)

式中，τ為距離向時(shí)間，t為方位向時(shí)間，r(t)為天線相位中心到目標(biāo)的瞬時(shí)距離，ra(t)為天線相位中心到場景中心的瞬時(shí)距離，k為距離向調(diào)頻率，c為光速，fc為載波頻率。

PFA距離向處理包括去除去調(diào)頻引起的殘余視頻相位誤差(RVP)、距離向尺度變換和距離向FFT，其中距離向尺度變換處理實(shí)現(xiàn)了距離向重采樣。距離向處理中的尺度變換函數(shù)分別如下：

(3)

(4)

(5)

式中：t為方位向慢時(shí)間；τ為距離向快時(shí)間；δr為距離向尺度變換因子；α為一個(gè)虛引數(shù)，其在一定范圍內(nèi)具有取值任意性；T為整個(gè)采樣時(shí)間；fτ為距離向頻域；k為距離向調(diào)頻率；r0為場景中心到雷達(dá)航跡的垂線長度。

方位向處理包括方位向尺度變換和方位向FFT，其中方位向尺度變換實(shí)現(xiàn)方位向重采樣。方位向處理中的尺度變換函數(shù)分別如下：

h1(t)=exp[jπkat2]

(6)

(7)

h2(t)=exp[-jπδakat2]

(8)

(9)

式中，ka為多普勒頻率，δa為方位向尺度變換因子，ft為方位向頻域。通過上述處理，采樣點(diǎn)排列格式從極坐標(biāo)轉(zhuǎn)換到了直角坐標(biāo)。最后經(jīng)過距離向傅里葉變換將數(shù)據(jù)變換到圖像域，單幀SAR成像完成。

2 圓跡視頻SAR實(shí)時(shí)成像的GPU實(shí)現(xiàn)

本文提出的基于GPU的圓跡視頻SAR實(shí)時(shí)成像算法如圖4所示，首先根據(jù)載頻和幀率確定重疊率截取回波數(shù)據(jù)，并將其傳輸至GPU等待處理；然后在GPU上利用PFA處理數(shù)據(jù)；最后將成像結(jié)果傳輸至CPU即可視頻成像。

圖4 基于CPU+GPU的視頻SAR實(shí)時(shí)成像流程

2.1 GPU成像處理關(guān)鍵技術(shù)

根據(jù)上節(jié)提出的視頻SAR成像算法，對于PFA運(yùn)算模塊，需要執(zhí)行CPU→GPU數(shù)據(jù)復(fù)制、數(shù)據(jù)處理、GPU→CPU數(shù)據(jù)復(fù)制3個(gè)步驟；對于PFA中方位向數(shù)據(jù)處理之前進(jìn)行先排距離向到先排方位向的轉(zhuǎn)置操作，以及在數(shù)據(jù)處理之后進(jìn)行先排方位向到先排距離向的轉(zhuǎn)置操作均利用GPU執(zhí)行；對于PFA中的Chirp Scaling操作需利用GPU執(zhí)行點(diǎn)乘運(yùn)算。為了利用 CUDA 技術(shù)將這種實(shí)現(xiàn)方案在GPU上部署，并充分利用GPU設(shè)備的運(yùn)算資源，給出了3種優(yōu)化技術(shù)對程序進(jìn)行優(yōu)化處理。

2.1.1 異步并行技術(shù)

將PFA的運(yùn)算模塊部署在GPU上時(shí)，需要在CPU內(nèi)存與GPU顯存之間傳輸數(shù)據(jù)。為了讓內(nèi)存與顯存之間的數(shù)據(jù)傳輸和PFA模塊的核(kernel)函數(shù)執(zhí)行并行實(shí)現(xiàn)，本文采用異步并行技術(shù)對數(shù)據(jù)處理加速。異步并行技術(shù)即創(chuàng)建若干個(gè)流，將雷達(dá)數(shù)據(jù)平均分配到每個(gè)流上，不同流中的數(shù)據(jù)傳輸與處理不會(huì)相互干擾，這樣可以使內(nèi)存與顯存之間的數(shù)據(jù)交互與kernel函數(shù)并行執(zhí)行。異步并行方案的實(shí)現(xiàn)如圖5所示，圖中同時(shí)比較了串行執(zhí)行與4個(gè)流異步執(zhí)行的運(yùn)行時(shí)間，可以看出異步并行執(zhí)行保證了GPU運(yùn)算核心絕大部分時(shí)間處于忙碌狀態(tài)，有效掩蓋了數(shù)據(jù)在內(nèi)存與顯存中的傳輸時(shí)間。

圖5 異步并行執(zhí)行與串行執(zhí)行

2.1.2 CUDA的兩層并行技術(shù)

PFA中完成距離向處理后，雷達(dá)數(shù)據(jù)以先排距離向的方式存儲(chǔ)于CPU內(nèi)存，因此在進(jìn)行方位向操作時(shí)，需要先對數(shù)據(jù)作轉(zhuǎn)置處理。針對矩陣轉(zhuǎn)置，本文采用CUDA的兩層并行技術(shù)。兩層并行運(yùn)算：一是在同一個(gè)線程塊(block)中利用線程實(shí)現(xiàn)的需要進(jìn)行數(shù)據(jù)交換和通信的細(xì)粒度并行，二是在各個(gè)block間實(shí)現(xiàn)不需要進(jìn)行數(shù)據(jù)交換的粗粒度并行。如圖6所示，先用合并訪問方式將數(shù)據(jù)從全局存儲(chǔ)器讀入共享存儲(chǔ)器，然后每個(gè)線程與其按對角線對稱的線程交換數(shù)據(jù)完成細(xì)粒度并行(效果如圖中數(shù)字1與2位置交換)，再按合并訪問方式將結(jié)果寫到全局存儲(chǔ)器中完成粗粒度并行。這種技術(shù)利用共享存儲(chǔ)器能極大提高程序的運(yùn)行效率，同時(shí)有效避免了存儲(chǔ)體沖突。

圖6 兩層并行完成矩陣轉(zhuǎn)置示意圖

2.1.3 分塊點(diǎn)乘技術(shù)

整個(gè)算法的核心部分是PFA中的Chirp Scaling操作，該操作部署在GPU上時(shí)依賴于點(diǎn)乘實(shí)現(xiàn)，所以需要對點(diǎn)乘的實(shí)現(xiàn)進(jìn)行優(yōu)化。為了結(jié)合異步并行技術(shù)，本文對雷達(dá)數(shù)據(jù)進(jìn)行分塊點(diǎn)乘運(yùn)算。兩個(gè)二維矩陣的分塊點(diǎn)乘運(yùn)算如圖7所示，矩陣A為雷達(dá)數(shù)據(jù)，矩陣B為尺度變換函數(shù)。雷達(dá)數(shù)據(jù)分流后距離向和方位向的點(diǎn)數(shù)分別為Nr和sub_Na，每個(gè)線程塊(block)的維數(shù)分別為Ta和Tr，則總線程塊的行數(shù)和列數(shù)分別為|(Nr+Tr-1)/Tr|和|(sub_Na+Ta-1)/Ta|，|x|表示取不大于x的整數(shù)。分塊的數(shù)目與線程塊的分塊相同，各塊的數(shù)據(jù)量與調(diào)用的線程(thread)數(shù)相同。按照分塊的方法，線程塊將分塊后的點(diǎn)乘因子讀入到線程塊內(nèi)的共享存儲(chǔ)器中，然后將分流后雷達(dá)數(shù)據(jù)中的每個(gè)元素與對應(yīng)尺度變換函數(shù)的點(diǎn)乘因子進(jìn)行并行點(diǎn)乘運(yùn)算即可。這種技術(shù)利用線程塊內(nèi)的共享存儲(chǔ)器節(jié)省了數(shù)據(jù)拷貝時(shí)間。

圖7 線程分塊點(diǎn)乘示意圖

2.2 基于GPU的視頻SAR成像流程

2.2.1 回波數(shù)據(jù)截取與傳輸

在成像處理開始之前，首先需要根據(jù)載頻與幀率的關(guān)系確定數(shù)據(jù)重疊率劃分子孔徑獲取單幀回波數(shù)據(jù)，并為回波數(shù)據(jù)進(jìn)行內(nèi)存及顯存分配，在所有步驟執(zhí)行完畢之后，再將所分配的內(nèi)存及顯存釋放。為了保證該內(nèi)存始終駐留在物理內(nèi)存中，且提高內(nèi)存與顯存間復(fù)制數(shù)據(jù)的速率，使用cudaHostAlloc分配頁鎖定主機(jī)內(nèi)存。

在進(jìn)行數(shù)據(jù)傳輸時(shí)，使用異步并行技術(shù)分成4個(gè)流異步并行執(zhí)行，并為每個(gè)流在GPU設(shè)備顯存上分配存儲(chǔ)空間用于完成PFA模塊的各個(gè)計(jì)算步驟。因?yàn)镻FA模塊涉及到的距離向處理和方位向處理中都包括FFT(或IFFT)運(yùn)算，所以有必要為每個(gè)流根據(jù)數(shù)據(jù)塊二維尺寸大小調(diào)用CUDA函數(shù)中的cufftPlan1d構(gòu)造FFT句柄(plan)，并將所創(chuàng)建plan用cufftSetStream函數(shù)與所屬流關(guān)聯(lián)。模塊執(zhí)行完畢后立即調(diào)用cufftDestroy函數(shù)銷毀plan，調(diào)用cudaStreamDestroy等待指定流中所有之前的任務(wù)完成，釋放流并將控制權(quán)返回給主機(jī)線程。

2.2.2 基于GPU的PFA成像

PFA距離向處理包括FFT、IFFT以及雷達(dá)數(shù)據(jù)與尺度變換函數(shù)相乘三項(xiàng)。FFT和IFFT調(diào)用CUDA庫函數(shù)。雷達(dá)數(shù)據(jù)與尺度變換函數(shù)相乘涉及到兩個(gè)二維矩陣的點(diǎn)乘運(yùn)算，即兩個(gè)矩陣的每一點(diǎn)對應(yīng)元素相乘，這里采用線程分塊點(diǎn)乘技術(shù)來實(shí)現(xiàn)，如圖7所示。尺度變換函數(shù)的生成通過在每個(gè)并行節(jié)點(diǎn)上利用GPU并行計(jì)算出其對應(yīng)數(shù)值。

為了避免不必要的時(shí)間損失，調(diào)用cudaMemcpyDevice函數(shù)在顯存中高效地完成距離向處理后數(shù)據(jù)的存儲(chǔ)，并調(diào)用cudaDeviceSynchronize使用CUDA同步技術(shù)，強(qiáng)制運(yùn)行時(shí)等待所有流中的任務(wù)都完成再進(jìn)入下一步驟。

PFA方位向處理同樣包括FFT、IFFT以及雷達(dá)數(shù)據(jù)與尺度變換函數(shù)相乘三項(xiàng)。因?yàn)榫嚯x向處理結(jié)束后方位向雷達(dá)數(shù)據(jù)在內(nèi)存中是離散存儲(chǔ)的，為了在方位向處理中能用異步并行分流處理節(jié)省，首先需要進(jìn)行矩陣轉(zhuǎn)置，將方位向數(shù)據(jù)變?yōu)閮?nèi)存中的連續(xù)存儲(chǔ)形式，矩陣轉(zhuǎn)置利用CUDA的兩層并行技術(shù)如圖6所示構(gòu)造轉(zhuǎn)置kernel。矩陣轉(zhuǎn)置后沿距離向均分為Cr塊并分流，再進(jìn)行方位向處理，方位向處理的GPU實(shí)現(xiàn)方式與距離向處理相同。

方位向處理通過CUDA同步技術(shù)保證所有流執(zhí)行完畢后，再次調(diào)用轉(zhuǎn)置kernel將顯存中距離向數(shù)據(jù)變?yōu)檫B續(xù)存儲(chǔ)，然后對距離向數(shù)據(jù)作IFFT，PFA模塊即完成。

最后將結(jié)果由cudaMemcpyAsync函數(shù)異步寫回主機(jī)內(nèi)存的對應(yīng)位置，并將處理結(jié)果寫入圖像文件，第一幀圖像處理結(jié)束。陸續(xù)讀取之后的子孔徑數(shù)據(jù)進(jìn)行處理并顯示，最終將以視頻的形式顯示成像結(jié)果。

3 實(shí)驗(yàn)結(jié)果與分析

為了驗(yàn)證本文所提方案的正確性和有效性，下面針對實(shí)測數(shù)據(jù)進(jìn)行了處理。成像環(huán)境如下：處理器Intel(R) Xeon(R) CPU E5-1650 3.2 GHz；內(nèi)存32 GB；顯卡型號(hào)NIVIDA TeslaC2075 1.15 GHz；顯存4 GB；顯卡運(yùn)算能力2.0；操作系統(tǒng)Win 7；軟件環(huán)境Visual Studio 2010及CUDA 5.5。實(shí)驗(yàn)所用原始數(shù)據(jù)是微型SAR按圓跡模式逆時(shí)針飛行獲得，數(shù)據(jù)大小總共為2 048×56 112，每幀成像數(shù)據(jù)大小為 2 048×2 048，實(shí)驗(yàn)參數(shù)如表1所示。因?yàn)檩d頻為9.7 GHz，由式(1)計(jì)算得重復(fù)率為99%時(shí)，幀率能達(dá)到5 Hz，即滿足視頻SAR成像要求，所以本實(shí)驗(yàn)中回波數(shù)據(jù)重疊率選為99%。

表1 視頻SAR實(shí)驗(yàn)系統(tǒng)參數(shù)

從實(shí)驗(yàn)結(jié)果中截取四幀成像結(jié)果如圖8(a)、8(b)、8(c)和8(d)所示，為了清晰地看出圓跡SAR的運(yùn)行軌跡，本文每隔500幀截取一幀。從成像結(jié)果中可以明顯看出載機(jī)圍繞中心目標(biāo)按圓跡模式逆時(shí)針飛行，也清楚地看到凹凸的地形、5個(gè)建筑物和周圍的樹木，具有非常高的分辨率。

(a) 第100幀

(b) 第600幀

(c) 第1 100幀

(d) 第1 600幀圖8 視頻SAR成像結(jié)果

本文利用Chirp Scaling操作替換PFA中傳統(tǒng)的二維插值，相同性能的GPU分別實(shí)現(xiàn)這兩種方式，并處理同一幀2 048×2 048大小的數(shù)據(jù)。如表2所示，基于插值的PFA處理一幀圖像用時(shí)0.25 s，而基于Chirp Scaling操作的PFA處理一幀圖像用時(shí)0.18 s，Chirp Scaling操作比插值節(jié)省38.8%的時(shí)間?？梢?，基于Chirp Scaling操作的PFA是更高效的處理算法。

提高計(jì)算速度的關(guān)鍵是應(yīng)用 GPU達(dá)到CPU不能實(shí)現(xiàn)的處理速度，下面將GPU上的處理時(shí)間與CPU上的處理時(shí)間進(jìn)行比較，如表2所示。處理2 048×2 048大小的數(shù)據(jù)，CPU的處理時(shí)間為6.27 s，明顯不能達(dá)到視頻SAR的成像幀率要求；GPU的處理時(shí)間為0.18 s，與單線程的CPU處理相比GPU處理的加速比為34.83倍。可見，選用 GPU作為視頻SAR數(shù)據(jù)處理平臺(tái)具有很大的優(yōu)勢。

最后，基于GPU的視頻SAR成像處理算法能否滿足實(shí)時(shí)性要求，取決于成像幀率。利用本文所用算法處理一幀圖所需時(shí)間約為0.18 s。對于距離采樣點(diǎn)數(shù)為2 048的視頻SAR系統(tǒng)，一部Tesla C2075便可滿足視頻SAR 5 Hz的視頻幀率要求。因此利用GPU并行處理高分辨率頻域處理算法進(jìn)行視頻SAR實(shí)時(shí)成像處理是可行的。另外，在此基礎(chǔ)上，還可以通過增加GPU設(shè)備個(gè)數(shù)擴(kuò)充設(shè)備的處理能力，以達(dá)到增加距離采樣點(diǎn)數(shù)的要求。

表2 實(shí)測數(shù)據(jù)成像時(shí)間比較

4 結(jié)束語

本文提出了一種基于GPU的視頻SAR實(shí)時(shí)成像算法，給出了算法的詳細(xì)流程，其中主要包括回波數(shù)據(jù)劃分方法以及PFA的CUDA實(shí)現(xiàn)方法。為了讓系統(tǒng)更高效，對整個(gè)成像算法的并行處理進(jìn)行了優(yōu)化，極大地發(fā)揮了GPU并行計(jì)算的優(yōu)勢。最后在GPU并行平臺(tái)上進(jìn)行了視頻SAR成像實(shí)驗(yàn)，驗(yàn)證了并行算法性能，有效解決了視頻SAR高分辨率成像的實(shí)時(shí)性問題。