一種深度欠采樣的快速SAR成像算法

王 健 宗竹林

(電子科技大學(xué) 成都 611731)

1 引言

合成孔徑雷達(dá)能夠全天時(shí)、全天候?qū)h(yuǎn)距離目標(biāo)進(jìn)行二維高分辨成像，具有廣泛的應(yīng)用。隨著應(yīng)用的不斷深入，人們對(duì)成像的質(zhì)量，如圖像的分辨率，可判讀性及視覺(jué)效果，都提出了更高的要求。傳統(tǒng)的成像算法是基于奈奎斯特采樣定理，由于SAR回波原始數(shù)據(jù)量非常大，給數(shù)據(jù)的傳輸和存儲(chǔ)帶來(lái)很大困難[1，2]。

CS理論是由Donoho等人從信號(hào)稀疏分解和逼近理論進(jìn)一步發(fā)展的一種新的信號(hào)處理理論，它徹底改變了人們對(duì)信息獲取的傳統(tǒng)觀(guān)念，由于突破了奈奎斯特采樣定理的限制，該理論一經(jīng)提出，就在信號(hào)/圖像處理、醫(yī)療成像、模式識(shí)別、地質(zhì)勘探、光學(xué)/雷達(dá)成像、無(wú)線(xiàn)通信等領(lǐng)域受到高度關(guān)注，并被美國(guó)科技評(píng)論評(píng)為2007年度十大科技進(jìn)展[3]。CS在雷達(dá)領(lǐng)域的應(yīng)用也得以迅速展開(kāi)。文獻(xiàn)[4]首次提出了隨機(jī)濾波的信號(hào)處理方法，實(shí)現(xiàn)了基于CS的SAR成像;文獻(xiàn)[5]采用Alltop sequence作為雷達(dá)發(fā)射信號(hào)，構(gòu)造了一個(gè)冗余字典，并在模糊域上實(shí)現(xiàn)了高分辨率雷達(dá)成像，突破了最大不模糊距離的限制。此外許多國(guó)內(nèi)研究者[6～9]也在CS雷達(dá)成像上取得了一些研究成果。

目前CS在雷達(dá)成像上的研究主要集中在距離向或方位向的一維壓縮感知成像，這些算法需要利用一個(gè)完整的合成孔徑時(shí)間內(nèi)的回波數(shù)據(jù)才能實(shí)現(xiàn)二維成像。本文以進(jìn)一步減少SAR成像所需數(shù)據(jù)量，縮短成像時(shí)間為目的，提出了一種深度欠采樣的壓縮感知成像算法，并在此基礎(chǔ)上利用OMP和BCS算法對(duì)壓縮后的數(shù)據(jù)恢復(fù)重構(gòu)，實(shí)現(xiàn)了雷達(dá)二維成像。由于該算法所需的回波數(shù)據(jù)小于一個(gè)合成孔徑時(shí)間內(nèi)的回波數(shù)據(jù)，從而等效的減小數(shù)據(jù)獲取時(shí)間，加快了成像速度。

2 壓縮感知基本理論

CS是建立在信號(hào)稀疏表示基礎(chǔ)上的信號(hào)采樣理論。該理論指出，任意N維K稀疏信號(hào)的稀疏表示X可以通過(guò)在RN空間的M=O(K·log(N/K))，(M?N)維隨機(jī)測(cè)量值y以很高的概率精確重構(gòu)。

假設(shè)一個(gè)離散信號(hào) X=(x1，x2，…，xN)，長(zhǎng)度為N，如果信號(hào)X是稀疏的或可壓縮的，則存在一組稀疏基 Ψ =[Ψ1，Ψ2，…，ΨN]使得X可以表示為:

其中Ψ是一個(gè)N×N的稀疏矩陣，σ是N×1的稀疏系數(shù)，包含K個(gè)非零值。

隨機(jī)測(cè)量過(guò)程是將信號(hào)X投影到一組測(cè)量矩陣Φ上，使測(cè)量對(duì)象從N維降為M維，從而得到測(cè)量值。隨機(jī)測(cè)量過(guò)程的數(shù)學(xué)表示為:

式中y表示M維隨機(jī)測(cè)量值，Φ為M×N維測(cè)量矩陣，ne為零均值高斯白噪聲。

當(dāng)矩陣θ滿(mǎn)足RIP準(zhǔn)則[10]時(shí)，就可以運(yùn)用壓縮感知重構(gòu)算法恢復(fù)信號(hào)X。

2.1 BCS[12] 算法

y的似然函數(shù)為:

使用最大似然法求上式的X會(huì)導(dǎo)致過(guò)學(xué)習(xí)。為避免這種情況，對(duì)權(quán)值X賦予先驗(yàn)的條件概率分布:

其中α =[α1，…，αN]T是一個(gè)控制權(quán)值先驗(yàn)方差的N維向量。如果能求得α，就很容易得到X。α可通過(guò)最小化下面的代價(jià)函數(shù)而得到。

使用EM算法求解上式，得到的第(k+1)步的迭代公式如下:

2.2 OMP[11]算法流程

A.初始化測(cè)量次數(shù)，初始化殘差:r=y;

B.利用最小二乘法找到與y相似最高的θ的列即內(nèi)積最大的列，最小二乘法步驟如下:

a.將θ的每一列代入＜rH，θi＞，找到使內(nèi)積＜rH，θi＞ 最大的一列記為 θq1;

b.A=θq1，并將對(duì)應(yīng)的 θ的第 i列置零，θ(:，i)=0;

C.返回(2)進(jìn)行下一次迭代，此時(shí) A={θq2，θq1};

D.迭代終止。其中 A=diag(α0，…，αN)。對(duì)于稀疏信號(hào)表示的問(wèn)題，以上迭代算法可寫(xiě)成:

3 基于壓縮感知的SAR成像

3.1 深度欠采樣回波信號(hào)模型

圖1中給出了正側(cè)視模式下SAR平臺(tái)的幾何模型，X軸表示方位向，Y軸表示平臺(tái)速度方向，Z軸為高度向。

圖1 正側(cè)視SAR幾何關(guān)系

雷達(dá)以速度V沿y軸方向勻速運(yùn)動(dòng)，在一個(gè)合成孔徑時(shí)間內(nèi)，連續(xù)發(fā)射線(xiàn)性調(diào)頻脈沖信號(hào):

在慢時(shí)間τ=m×PRT時(shí)刻，雷達(dá)平臺(tái)至點(diǎn)目標(biāo)k的距離可表示為rk(τ)，雷達(dá)系統(tǒng)接收到的回波信號(hào)可以描述為:

其中σk為散射點(diǎn)k的后向散射系數(shù);f0為載頻;c為光速;Tp為線(xiàn)性調(diào)頻脈沖信號(hào)的時(shí)間寬度。

根據(jù)CS可以完成不完全數(shù)據(jù)重構(gòu)的思想，本文利用前1/2個(gè)合成孔徑時(shí)間內(nèi)的數(shù)據(jù)來(lái)獲取測(cè)量信號(hào)，這樣縮短了數(shù)據(jù)獲取的時(shí)間，加快了成像的速度。

圖2 回波采樣示意圖

圖2表示RD算法和CS算法對(duì)回波數(shù)據(jù)采樣的示意圖，“■”表示采樣點(diǎn)，“□”表示非采樣點(diǎn)。圖2(b)表示對(duì)一個(gè)合成孔徑時(shí)間內(nèi)的接收到的前50%的數(shù)據(jù)采樣，并且選擇隨機(jī)高斯矩陣[13]作為測(cè)量矩陣。測(cè)量過(guò)程如下:

其中Φ是M×N維隨機(jī)高斯矩陣，即隨機(jī)選擇N維信號(hào)X中的M維，從而得到欠采樣信號(hào)y。隨機(jī)測(cè)量矩陣可以表示為:

Φ中“1”表示隨機(jī)測(cè)量數(shù)據(jù)，“0”表示丟棄的數(shù)據(jù)，CS算法即由對(duì)應(yīng)于Φ中“1”的測(cè)量數(shù)據(jù)重構(gòu)出原始信號(hào)的過(guò)程。

3.2 基于壓縮感知的SAR成像

假設(shè)雷達(dá)平臺(tái)在一個(gè)合成孔徑時(shí)間內(nèi)掃描的二維場(chǎng)景為ξ，將二維場(chǎng)景區(qū)域均勻離散成Na×Nr個(gè)點(diǎn)散色目標(biāo)，將Na×Nr個(gè)點(diǎn)目標(biāo)回波延遲按行展開(kāi)，構(gòu)成一個(gè)一維向量，長(zhǎng)度為N=Na×Nr，N個(gè)點(diǎn)目標(biāo)的坐標(biāo)可以表示為 β =(β1，β2，…，βN)，其中 βi=(xi，yi，zi)是點(diǎn)目標(biāo)i的坐標(biāo)向量，對(duì)應(yīng)的散色系數(shù)可以表示為 σ =(σ1，σ2，…，σN)，為保證信號(hào)的稀疏性，假設(shè)只有少數(shù)的K個(gè)散色系數(shù)是非零的，即:

根據(jù)(11)式和(14)式，可以得知對(duì)于N個(gè)點(diǎn)目標(biāo)β =(β1，β2，…，βN)，在Ψ變換域下可以實(shí)現(xiàn)對(duì)雷達(dá)回波的稀疏表示。對(duì)于每個(gè)點(diǎn)目標(biāo)，變換域Ψ可以表示為 Ψ =(φ1，φ2，…，φN)，其中點(diǎn)目標(biāo) βi對(duì)應(yīng)于變換基φi，按行展開(kāi)得到N×1的向量，表示為:

式中，Ri是點(diǎn)目標(biāo)i與平臺(tái)間的斜距;yi是點(diǎn)目標(biāo)i在y軸方向與點(diǎn)目標(biāo)之間的距離;Lsar為合成孔徑長(zhǎng)度。

將上述變換基和(13)式代入(12)式，便可得到測(cè)量信號(hào)，由于隨機(jī)采樣在實(shí)現(xiàn)過(guò)程較為復(fù)雜，在實(shí)際中可以采用等間隔采樣來(lái)取代隨機(jī)高斯測(cè)量，在硬件實(shí)現(xiàn)上可以用低速率AD對(duì)回波數(shù)據(jù)采樣獲得測(cè)量信號(hào)。綜上所述，CS雷達(dá)成像處理流程可以用圖3描述。

圖3 CS雷達(dá)成像框架

圖3中的A/D是以較低的采樣率對(duì)雷達(dá)回波信號(hào)采樣，得到M×1欠采樣測(cè)量信號(hào)y=(y1，y2，…，yM)，這樣在得到了傳感矩陣和測(cè)量信號(hào)后，就可以用壓縮感知重構(gòu)算法對(duì)回波信號(hào)做脈沖壓縮。脈沖壓縮過(guò)程即(12)式中獲取測(cè)量信號(hào)的逆過(guò)程，通過(guò)測(cè)量信號(hào)y和傳感矩陣θ=ΦΨ得到散色系數(shù)σ。

本文在 CS脈沖壓縮過(guò)程采用 OMP[11]和BCS[12]兩種算法實(shí)現(xiàn)。OMP算法利用上述傳感矩陣θ=ΦΨ和測(cè)量信號(hào)y可以實(shí)現(xiàn)對(duì)欠采樣回波的脈沖壓縮。BCS算法需要對(duì)θ=ΦΨ和y做如下變換來(lái)得到傳感矩陣和測(cè)量信號(hào)。

上式將測(cè)量信號(hào)的實(shí)部虛部分開(kāi)來(lái)獲取實(shí)測(cè)量信號(hào) yr。

對(duì)應(yīng)的變換域是對(duì)Ψ中的φi做如下變換。

從 而 得 到 BCS 的 變 換 域 Ψr= (φr1，φr2，…，φrN)。

將(13)式中測(cè)量矩陣擴(kuò)展為2M×2N維，得到BCS算法隨機(jī)測(cè)量矩陣Φr，利用BCS算法做脈沖壓縮的過(guò)程表示如下:

圖3中后處理過(guò)程即將得到的一維散色系數(shù)轉(zhuǎn)化為二維圖像的過(guò)程，OMP和BCS算法的后處理過(guò)程可以表示如下:

4 仿真結(jié)果

假設(shè)雷達(dá)工作在正側(cè)視模式下，場(chǎng)景中心到載機(jī)航線(xiàn)的垂直距離Rb=2000m，距離向欠采樣率和方位向欠采樣率分別為Dr和Da，距離分辨率和方位分辨率為1m，仿真中用到的雷達(dá)參數(shù)有:

信號(hào)帶寬:150MHz 采樣率:200MHz

脈沖寬度:0.5μs 脈沖重復(fù)頻率:500Hz

載機(jī)速度:150m/s 載頻:10GHz

圖4 重構(gòu)誤差與方位向降采樣率

本文首先分析了在距離向降采樣率Dr=1/5，方位向降采樣率與重構(gòu)誤差的關(guān)系，在15個(gè)點(diǎn)目標(biāo)的場(chǎng)景區(qū)域，仿真結(jié)果如圖4所示。

由圖4可以看出，在圖中降采樣率情況下重構(gòu)誤差可以滿(mǎn)足精確重構(gòu)的要求，所以在后面的仿真中選用方位向?qū)⒉蓸勇蔇a=1/2，即半個(gè)合成孔徑時(shí)間的數(shù)據(jù)。

圖5中顯示的是在無(wú)噪聲情況下傳統(tǒng)的RD算法成像結(jié)果。圖6為無(wú)噪聲情況下壓縮感知算法重構(gòu)的目標(biāo)場(chǎng)景區(qū)域，其中圖6(a)是由OMP算法得到的成像結(jié)果，圖6(b)是BCS算法的重構(gòu)結(jié)果?？梢钥闯鯟S算法在無(wú)噪聲情況下可以完美的重建圖像信息。

由于信噪比大小對(duì)RD算法成像的影響是我們所熟悉的，下面是在SNR不同情況下的仿真給出CS方法成像的結(jié)果，對(duì)應(yīng)的RD算法成像效果與圖5并無(wú)明顯區(qū)別。

圖7 SNR=5dB情況的成像結(jié)果

圖7是在信噪比SNR=5dB情況下，對(duì)目標(biāo)空間的重構(gòu)結(jié)果。圖7(a)和圖7(b)分別是OMP算法和BCS算法的重構(gòu)結(jié)果，從圖中可以看出，出現(xiàn)了一些能量較弱的“假目標(biāo)”。

圖8 SNR=25dB情況的成像結(jié)果

圖8是在SNR=25dB情況下，OMP算法和BCS算法重構(gòu)的二維場(chǎng)景，可以看出在25dB時(shí)，噪聲對(duì)重構(gòu)結(jié)果的影響很小，場(chǎng)景區(qū)域得到了正確的重建。

圖9和圖10是RD算法與CS算法在單個(gè)點(diǎn)目標(biāo)下的三維圖像，與圖5和圖6二維圖相對(duì)應(yīng)。

通常方位向采樣數(shù)據(jù)的減少會(huì)直接降低方位向分辨率，由圖6和圖10與圖5圖9對(duì)比可知，本文提出的CS算法在降低方位向采樣數(shù)據(jù)的情況下，不僅沒(méi)有降低方位向的分辨率，并且該算法消除了目標(biāo)旁瓣，提高了成像分辨率。

對(duì)上述信噪比不同情況下，OMP算法和BCS算法的重構(gòu)結(jié)果進(jìn)行定量的分析，通過(guò)重構(gòu)誤差、重構(gòu)時(shí)間、重構(gòu)信號(hào)的LSNR幾個(gè)參量對(duì)重構(gòu)圖像進(jìn)行評(píng)估。其中重構(gòu)誤差定義為:

定義重構(gòu)信號(hào)的局部信噪比為L(zhǎng)SNR[7]，用來(lái)反映圖像的重構(gòu)效果，計(jì)算公式如下:

通過(guò)對(duì)不同SNR條件下重構(gòu)誤差的計(jì)算，得到圖11所示的曲線(xiàn)，可以看出，BCS和OMP算法的重構(gòu)誤差隨著輸入信噪比的增大，逐漸減少，并且逐漸趨近于零。

圖11 重構(gòu)誤差與SNR的關(guān)系

在重構(gòu)的過(guò)程中統(tǒng)計(jì)重構(gòu)所需時(shí)間，得到圖12中重構(gòu)時(shí)間與SNR的關(guān)系，重構(gòu)時(shí)間隨著信噪比的增大逐漸減少，在30dB之后逐漸趨于穩(wěn)定。由圖12可以看出信噪比的大小對(duì)重構(gòu)時(shí)間有較大的影響。

圖12 重構(gòu)時(shí)間與SNR的關(guān)系

圖13表示了重構(gòu)信號(hào)局部SNR和輸入信號(hào)SNR的關(guān)系，LSNR與輸入信號(hào)SNR成正比，并隨著輸入SNR增大而增大。圖13中的曲線(xiàn)表明，隨著輸入SNR的增大，重構(gòu)圖像的效果就越好。由此可知，圖11與圖13的定量分析與上述成像的視覺(jué)效果相吻合。

圖13 LSNR與SNR的關(guān)系

由以上仿真結(jié)果可以看出，在一個(gè)合成孔徑時(shí)間內(nèi)，傳統(tǒng)的RD成像方法掃描的數(shù)據(jù)量為Na×Nr，CS方法僅使用了Na/2×Nr/5的數(shù)據(jù)量，大幅度的降低了數(shù)據(jù)量，同時(shí)仿真僅用了一個(gè)合成孔徑時(shí)間內(nèi)的前50%的數(shù)據(jù)，這等效于在1/2個(gè)合成孔徑時(shí)間內(nèi)獲得的數(shù)據(jù)，因此縮短了數(shù)據(jù)獲取時(shí)間，實(shí)現(xiàn)了快速成像。與傳統(tǒng)SAR成像相比，CS成像算法消除了目標(biāo)旁瓣，提高了成像的分辨率，并且有效的降低了成像系統(tǒng)的數(shù)據(jù)率。

5 結(jié)論

本文提出了一種深度欠采樣快速SAR成像算法，在小于一個(gè)合成孔徑時(shí)間里獲得的數(shù)據(jù)進(jìn)行測(cè)量，獲取欠采樣回波數(shù)據(jù)，然后利用CS脈沖壓縮算法重構(gòu)，實(shí)現(xiàn)了快速高分辨雷達(dá)二維成像。該算法能夠有效地降低高分辨率雷達(dá)成像系統(tǒng)的數(shù)據(jù)率，減少硬件系統(tǒng)(特別是A/D轉(zhuǎn)換器)的實(shí)現(xiàn)難度，消除了目標(biāo)旁瓣，提高了成像的分辨率，同時(shí)縮短了成像時(shí)間。

[1]Sujit Bhattacharya，Thomas Blumensath，Bernard Mulgrew，et al.Fast encoding of synthetic aperture radar raw data using compressed sensing[J].IEEE Workshop on Statistical Signal Processing，2007:448－452.

[2]Gabriel Rilling，Mike davies and bernard mulgrew.Compressed sensing based compression of SAR raw data[C].Singnal Processing with A-daptive Sparse Structured Representations Workshop，Saint-Malo，F(xiàn)rance，2009.

[3]喻玲娟，謝曉春.壓縮感知理論簡(jiǎn)介[J].數(shù)字視頻，2008，32(12):16 －15.

[4]Baraniuk R and Steeghs P.Compressive radar imaging[C].IEEE Radar Conference，Boston，MA，USA，Apr.17－20，2007:128－133.

[5]Herman M and Strohmer T.Compressed sensing radar[C].IEEE International Conference on A-coustics，Speech and Signal Processing，Las Vegas，NV，USA，Mar.30 － Apr.4，2008:1509－1512.

[6]謝曉春，張?jiān)迫A.基于壓縮感知的二位雷達(dá)成像算法[J].電子與信息學(xué)報(bào)，2010，32(5):1234－1238.

[7]謝曉春.壓縮感知理論在雷達(dá)成像中的應(yīng)用研究[D].中國(guó)科學(xué)院空間科學(xué)與應(yīng)用研究中心，2010.

[8]Lei Zhang，Meng-dao Xing，Cheng-Wei Qiu，Jun Li，Jia-lian Sheng，Ya-chao Li，Zheng Bao.Resolution enhancement for inversed synthetic aperture radar imaging under low SNR via improved compressive sensing[J]，IEEE transactions on geosciences and remote sensing，2010，48(10):3824 －3838.

[9]Hongxian Wang，Yinghui Quan，Mengdao Xing，and Shouhong Zhang.ISAR imaging via sparse probing frequencies[J]，IEEE Transactions on Geosciences and Remote Sensing，2010，8(3):451－454.

[10]E Candès.The restricted isometry property and its implications for compressed sensing[J].Acadèmie des sciences，2006，346(I):598 －592.

[11]Tropp J A，Gilbert A C.Signal recovery from random measurements via orthogonal matching pursuit[J]..IEEE Transactions on Information Theory，2007，53(12):4655－4666.

[12]Ji S H，Xue Y，and Carin L.Bayesian compressive sensing[J].IEEE Transactions on Signal Processing，2008，56(6):2346 －2356.

[13]R Baraniuk.A lecture on compressive sensing[J].IEEE Signal Processing Magazine，2007，24(4):118－121.

[14]Sadegh Samadi，Müjdat Cetin，Mohammad Ali Masnadi-Shirazi.Sparse representation-based SAR imaging[J].IET Radar Sonar Nacig，2011，5(2):182 －193.