基于GPU的SAS成像算法并行實(shí)現(xiàn)研究＊

（海軍工程大學(xué)電子工程學(xué)院 武漢 430033）

1 引言

合成孔徑聲納（Synthetic Aperture Sonar，SAS）是一種高分辨率水下成像聲納，基本原理是［1～3］：利用小尺寸基陣沿空間勻速直線運(yùn)動來合成大孔徑基陣，獲得不依賴于距離和波長大小的高方位向分辨率圖像。SAS被廣泛應(yīng)用在海底地形測繪、海洋地質(zhì)考察、沉底物打撈以及水雷的探測和識別。隨著SAS成像技術(shù)的發(fā)展和實(shí)際應(yīng)用需求，對回波數(shù)據(jù)實(shí)時(shí)處理的要求也越來越迫切。由于SAS成像算法本身的復(fù)雜性和數(shù)據(jù)量大，傳統(tǒng)的處理方式難以滿足實(shí)時(shí)成像的要求，采用并行技術(shù)是實(shí)現(xiàn)高性能SAS成像的必經(jīng)之路。文獻(xiàn)［4］中指出，RD［5］算法能夠兼顧成熟、簡單、高效和精確等因素，至今仍是使用最廣泛的成像算法。本文首先給出了多子陣［9］合成孔徑聲納的RD算法，然后研究了基于GPU 的RD成像算法的并行處理方法。

2 多子陣RD 算法簡介

假定發(fā)射信號為線性調(diào)頻信號，易知多子陣合成孔徑聲吶第i個(gè)接收陣在距離為r處的點(diǎn)目標(biāo)的回波信號基帶形式［9］可寫為

式中A0為幅度函數(shù)，可視為常數(shù)，ωa（·）為方位向包絡(luò)函數(shù)，k為發(fā)射信號調(diào)頻斜率，c為聲速，f0為信號中心頻率，R0（t，r）為延時(shí)距離誤差。

首先對式（1）距離向快變時(shí)間和方位向慢變時(shí)間作傅立葉變換，進(jìn)入fr-fa域：

其中A2為常數(shù)，Wr（·）為距離頻譜包絡(luò)，Wa（·）為方位頻譜包絡(luò)，fr為距離向瞬時(shí)頻率，fa為方位向瞬時(shí)頻率，Δhi表示為第i接收子陣中心與發(fā)射陣中心的間距。

將上述公式中的根式對fr做泰勒展開并舍棄二次以上項(xiàng)［8］并作距離向進(jìn)行反傅立葉變換，進(jìn)入RD 域：

這里ks為距離多普勒域新的距離調(diào)頻率，D＝

從式（5）中可以發(fā)現(xiàn)，第一項(xiàng)為方位調(diào)制項(xiàng)，它包含兩部分內(nèi)容：方位脈壓函數(shù)φa1和方位走動函數(shù)φa2。

第二項(xiàng)為距離彎曲項(xiàng)，時(shí)間徙動量為Δτ：

其中Δτi為停走停近似及各個(gè)接收陣接收路徑不同造成的延遲誤差，與方位向無關(guān)，因而可以在距離脈壓后的回波數(shù)據(jù)中分別予以補(bǔ)償。

式（5）第三項(xiàng)為停走停近似及不同接收陣引入的固定相位φi，也與方位向無關(guān)，也可以在各接收陣回波數(shù)據(jù)中分別予以補(bǔ)償。

式（5）最后一項(xiàng)為只與各接收陣的方位位置有關(guān)的相位項(xiàng)，其物理意義可表征為各接收陣數(shù)據(jù)方位向的間隔為Δhi／2順序排列。

至此完成了多接收陣的距離多普勒算法，其流程圖如圖1所示。

圖1 多接收陣RD 成像算法流程圖

3 基于CUDA的GPU 并行技術(shù)簡介

CUDA 采用了統(tǒng)一處理架構(gòu)，降低了編程的難度，并且NVIDIA GPU 引入了片內(nèi)共享存儲器，提高了效率。這兩項(xiàng)改進(jìn)使CUDA 架構(gòu)更加適合進(jìn)行GPU 通用計(jì)算。CUDA 中線程按照兩個(gè)層次進(jìn)行組織、共享存儲器和柵欄同步［6］。這些關(guān)鍵特性使得CUDA 擁有了兩個(gè)層次的并行：線程級并行實(shí)現(xiàn)的細(xì)粒度數(shù)據(jù)并行，任務(wù)級并行實(shí)現(xiàn)的粗粒度并行。

CUDA編程模型是將CPU 作為主機(jī)端（Host），將GPU 作為設(shè)備端（Device），同時(shí)在一個(gè)系統(tǒng)中可以存在一個(gè)主機(jī)和多個(gè)設(shè)備或者多個(gè)主機(jī)和多個(gè)設(shè)備［7］，即GPU 集群。CUDA 程序的核心是內(nèi)核（Kernel）函數(shù)，其在設(shè)備端執(zhí)行并通過主機(jī)端調(diào)用，主要任務(wù)是完成數(shù)據(jù)的并行處理。一個(gè)完整的CUDA程序是由若干個(gè)運(yùn)行在Device端的Kernel函數(shù)并行處理步驟和Host端的串行處理步驟共同組成。在主機(jī)端調(diào)用時(shí)需要設(shè)置內(nèi)核函數(shù)的執(zhí)行參數(shù)。執(zhí)行參數(shù)設(shè)置主要是指線程塊、線程網(wǎng)格維度的設(shè)計(jì)、共享存儲空間的分配以及流的設(shè)計(jì)等。

4 RD 算法的GPU 并行實(shí)現(xiàn)

本文選擇RD 算法進(jìn)行SAS的成像［9］，首先對距離進(jìn)行脈沖壓縮，再進(jìn)行距離徙動校正，然后進(jìn)行方位向脈沖壓縮。距離脈壓和方位脈壓在這里使用頻域匹配濾波實(shí)現(xiàn)，主要涉及的是傅里葉變換、向量和矩陣的點(diǎn)乘和逆傅里葉變換。在CUDA中包含了CUFFT 庫，可以利用CUFFT 庫實(shí)現(xiàn)傅里葉變換［10］和逆傅里葉變換在GPU 上的并行計(jì)算。然后，通過編寫的核函數(shù)實(shí)現(xiàn)向量和矩陣的點(diǎn)乘和轉(zhuǎn)置運(yùn)算。具體算法流程如圖2所示。

圖2 基于CUDA 的RD 成像并行算法流程圖

利用CUDA開發(fā)的RD成像算法實(shí)現(xiàn)流程如下：

1）將回波數(shù)據(jù)和參數(shù)傳入到顯存。

2）調(diào)用CUFFT 庫實(shí)現(xiàn)距離向的傅里葉變換并行運(yùn)算，同時(shí)計(jì)算距離匹配函數(shù)、方位匹配函數(shù)和距離徙動參數(shù)。

3）編寫的點(diǎn)乘核函數(shù)實(shí)現(xiàn)距離向數(shù)據(jù)和距離匹配函數(shù)的點(diǎn)乘，然后調(diào)用CUFFT 庫并行實(shí)現(xiàn)距離向的逆傅里葉變換，完成距離向脈沖壓縮。

4）利用距離徙動參數(shù)在GPU 上實(shí)現(xiàn)距離徙動校正，然后用核函數(shù)實(shí)現(xiàn)數(shù)據(jù)的并行轉(zhuǎn)置。

5）調(diào)用CUFFT 庫實(shí)現(xiàn)方位向數(shù)據(jù)的傅里葉變換，然后用點(diǎn)乘核函數(shù)實(shí)現(xiàn)與方位匹配函數(shù)的點(diǎn)乘，最后調(diào)用CUFFT 庫并行實(shí)現(xiàn)方位向的逆傅里葉變換，完成方位向脈沖壓縮。

6）將圖像數(shù)據(jù)傳回主機(jī)端。

5 實(shí)驗(yàn)與結(jié)果

為了驗(yàn)證本文所提出的SAS成像方法的有效性，本文進(jìn)行了五個(gè)點(diǎn)的SAS 成像仿真實(shí)驗(yàn)。本實(shí)驗(yàn)采用的計(jì)算機(jī)CPU 為Intel（R）G840CPU 雙核2．8G，內(nèi)存4G，顯卡為GeForce 210，顯存為1G，運(yùn)算能力為1．2，運(yùn)算精度為單精度。軟件使用的是Windows XP下的CUDA 4．0．和VS2008。數(shù)據(jù)來源于合成孔徑聲納接收回波的計(jì)算機(jī)仿真。仿真參數(shù)如表1所示。

表1 SAS信號仿真參數(shù)

針對不同數(shù)據(jù)大小的回波，在CPU 和GPU 上分別進(jìn)行RD成像仿真效率的比較，其結(jié)果如表2所示。

表2 CPU 和GPU 上RD 算法效率比較

從表2的計(jì)算效率比較結(jié)果可以看出，隨著原始回波數(shù)據(jù)量的增加，并行RD 算法的加速比也在逐步增加，可見對于高維數(shù)的SAS回波，更能夠體現(xiàn)并行算法的優(yōu)勢。

圖3和圖4分別是RD 算法在GPU 和CPU上局部放大后的1024×1024個(gè)復(fù)數(shù)成像的結(jié)果。

圖3 GPU 成像結(jié)果

圖4 CPU 成像結(jié)果

從圖3和圖4可以看出，基于GPU 的仿真結(jié)果和基于CPU 的仿真結(jié)果在成像結(jié)果上的一致性，但是從表2可以看出兩者的運(yùn)算速度卻有明顯的差距。

由測試的結(jié)果可以知道，基于GPU 的RD成像算法，具有較高的效率，能夠滿足實(shí)時(shí)成像的要求。

6 結(jié)語

本文設(shè)計(jì)了基于GPU 的SAS成像并行算法，并對該并行成像算法進(jìn)行了計(jì)算機(jī)仿真，驗(yàn)證了算法的可行性和高效性。實(shí)驗(yàn)結(jié)果表明，基于CUDA開發(fā)出來的SAR 成像并行算法，實(shí)現(xiàn)了實(shí)時(shí)成像，可滿足下一步處理的需要。

［1］張學(xué)波，唐勁松，張森，等．多接收陣合成孔徑聲納線頻調(diào)變標(biāo)成像算法［J］．系統(tǒng)工程與電子技術(shù)，2013，35（7）：1415-1420．

［2］張學(xué)波，唐勁松，鐘何平．合成孔徑聲納多接收陣數(shù)據(jù)融合CS成像算法［J］．哈爾濱工程大學(xué)學(xué)報(bào)，2013，34（2）：240-244．

［3］楊海亮，張森，唐勁松．寬測繪帶多接收陣合成孔徑聲納成像處理方法［J］．系統(tǒng)仿真學(xué)報(bào)，2011，23（7）：1424-1428．

［4］禹衛(wèi)東．合成孔徑雷達(dá)信號處理研究［D］．南京：南京航空航天大學(xué)博士論文，1997．

［5］保錚，邢孟道，王彤．雷達(dá)成像技術(shù)［M］．北京：電子工業(yè)出版社，2006．

［6］張舒，諸艷利．GPU 高性能運(yùn)算之CUDA［M］．北京：水利水電出版社，2009．

［7］桑德斯（美）．GPU 高性能編程CUDA 實(shí)戰(zhàn)［M］．聶雪軍，等譯．北京：機(jī)械工業(yè)出版社，2011．

［8］Cumming I G，Wong F H．合成孔徑聲納雷達(dá)成像-算法與實(shí)現(xiàn)［M］．洪文，等譯．北京：電子工業(yè)出版社，2007，10．

［9］楊海亮．多接收陣合成孔徑聲納成像算法研究［J］．博士論文．海軍工程大學(xué)，2009．

［10］趙麗麗，張勝兵，張萌，等．基于CUDA 的高速FFT 計(jì)算［J］．計(jì)算機(jī)應(yīng)用研究，2011，28（4）：1556-1559．