張冰塵王萬影*③畢 輝③趙 曜洪 文

①(微波成像技術(shù)重點實驗室 北京 100190)

②(中國科學(xué)院電子學(xué)研究所 北京 100190)

③(中國科學(xué)院大學(xué) 北京 100190)

基于壓縮多信號分類算法的森林區(qū)域極化SAR層析成像

張冰塵①②王萬影*①②③畢 輝①②③趙 曜①②洪 文①②

①(微波成像技術(shù)重點實驗室 北京 100190)

②(中國科學(xué)院電子學(xué)研究所 北京 100190)

③(中國科學(xué)院大學(xué) 北京 100190)

該文研究了一種基于壓縮多信號分類算法的森林區(qū)域極化SAR層析成像方法。其具體步驟包括：全極化的SAR接收成像區(qū)域的反射回波，利用各極化通道的信號建立多觀測向量模型；應(yīng)用小波基對高程向結(jié)構(gòu)進行稀疏表示，采用壓縮多信號分類算法對觀測區(qū)域的高程向后向散射系數(shù)進行重建，實現(xiàn)對森林區(qū)域?qū)游龀上?。最后，通過仿真實驗、PolSARpro仿真數(shù)據(jù)和德宇航E-SAR的P-波段數(shù)據(jù)驗證了該方法在同等測量精度的要求下可以有效減少SAR層析成像所需的航過數(shù)，同時降低了虛假目標(biāo)的出現(xiàn)概率。

極化SAR；層析；壓縮多信號分類；小波基

1 引言

3D-SAR成像技術(shù)既繼承了傳統(tǒng)SAR系統(tǒng)所具備的全天時、全天候、高分辨率成像等優(yōu)點，又避免了2維成像中處于同一散射單元內(nèi)的目標(biāo)散射點與雷達(dá)間斜距相等時存在的疊掩效應(yīng)，能夠?qū)⒛繕?biāo)的高度向與距離向完全分離，實現(xiàn)對目標(biāo)的3維分辨能力[1,2]。

SAR層析成像具有3維成像能力。SAR層析成像是沿垂直于視線的法線方向排列多個天線或利用同一天線在法線方向的不同軌跡高度對同一目標(biāo)區(qū)域成像，來獲得高度向上的分辨力。SAR層析成像技術(shù)不僅能夠獲得目標(biāo)散射體的高程信息，同時還可以獲得散射體在高度向上的分布，能完全恢復(fù)真實3維場景。SAR層析成像能夠有效、快捷地實現(xiàn)3維成像，極大地擴大了SAR的應(yīng)用范圍，因此，已成為SAR技術(shù)發(fā)展的重要方向。SAR層析成像的方法包括譜估計和壓縮感知，其中，壓縮感知被廣泛應(yīng)用于雷達(dá)成像[3]。

1998年，德國宇航局(DLR)利用機載E-SAR系統(tǒng)進行了層析成像飛行試驗，并采用傅里葉變換聚焦的算法成功實現(xiàn)3維成像。2000年，Reigber等人[4]提出了機載SAR層析模型，介紹了多基線SAR層析成像的原理，并采用頻譜估計的方法實現(xiàn)了多基線L波段的層析成像。在SAR層析成像中，常用的譜估計算法有快速傅里葉變換(FFT)、多信號分類(MUSIC)以及CAPON等[5]。

傳統(tǒng)的頻譜估計方法信號序列長度受限，所能達(dá)到的高程向分辨率受限于高程向合成孔徑的大小，所以高程向分辨率相對較低。為了解決這一問題，壓縮感知的方法被應(yīng)用于SAR層析成像[6,7]。壓縮感知可以應(yīng)用于SAR層析成像的前提是高程向的散射體分布是稀疏的，或者經(jīng)過稀疏基處理之后滿足稀疏性[8]。這種方法實現(xiàn)了高程向的超分辨率成像，然而仍然需要大量的航過來保證重建精度，這在實際實驗中需要較高的成本。為了減少SAR層析成像所需的航過數(shù)提出了多信號的壓縮感知處理方法[9,10]。

文獻[11]提出了壓縮多信號分類(CS-MUSIC)算法，該算法采用壓縮感知的方法來識別部分支撐集，剩余的部分則利用一種廣義的MUSIC算法進行恢復(fù)。該方法結(jié)合多觀測向量實現(xiàn)對場景的重建，通過仿真實驗證明該算法在保證恢復(fù)精度的同時，又可以大大地減少所需的航過數(shù)。本文結(jié)構(gòu)如下，第2節(jié)給出森林區(qū)域SAR層析的多觀測問題模型，通過小波基對森林區(qū)域場景進行稀疏化，并在后面的實驗驗證部分對實測P-波段數(shù)據(jù)結(jié)果進行多視處理；第3節(jié)介紹了CS-MUSIC算法的具體實施步驟；第4節(jié)利用仿真實驗說明該算法應(yīng)用于SAR層析成像中與其他算法相比的優(yōu)越性；最后，第5節(jié)對全文進行總結(jié)。

2 信號模型

SAR層析成像的目的在于對同一目標(biāo)場景在不同位置進行多次重復(fù)觀測，以獲得一個目標(biāo)的多個投影，利用這些投影來計算物體的3維圖像。SAR層析成像利用垂直于斜距方向的不同基線形成高程向的一個合成孔徑，來獲得在高程方向上合適的分辨率。

SAR層析成像沿著每一個軌道獲得的回波信號都是一幅由方位-距離像素點組成的2維SAR圖像，每一個像素點代表該方位-距離分辨單元內(nèi)高程向所有目標(biāo)點散射回波的疊加，表示如下[4]：

其中，a, r和v分別代表方位向、距離向和高程向，ym(a,r)為第m條基線處獲得的2維SAR圖像中位于(a,r)處的像素值，x(a,r,v)為高程向后向散射系數(shù)，λ為載波波長，r為觀測點到場景中心的斜距，bm為第m條基線高程向高度。

對其進行高程向的離散化處理，得到SAR層析問題的數(shù)學(xué)模型，可表示為

此時，對于有M次航過的SAR層析成像情況，我們可以用下面的矩陣相乘的形式來表示HH通道方位-距離分辨單元的回波采樣向量，且。此時SAR層析成像數(shù)學(xué)模型可以表示為

其中N為高程向網(wǎng)格數(shù)。我們稱A為觀測矩陣，單觀測的SAR層析成像問題可以表示為

我們所要求解的是后向散射系數(shù)向量XHH，且XHH∈CN×1，即可獲得高程向散射體分布情況。對于文中多極化信號的情況，回波信號由3個極化通道信號組成，此時回波采樣矩陣，且滿足X∈CN×3。此時，多信號SAR層析問題數(shù)學(xué)模型表示為

由于不同極化通道獲得的回波信號具有一定的相關(guān)性，并且都是針對于同一目標(biāo)場景，所以我們可以假設(shè)不同通道高程向后向散射系數(shù)均具有相同的支撐集[9]。

此時，我們所要解決的多信號壓縮感知問題表示為

實際上，森林區(qū)域高程向一般并不具備稀疏性。所以在應(yīng)用壓縮感知的方法進行森林區(qū)域?qū)游龀上駮r，需要先將其在小波域進行稀疏表示[12]。這里令W∈CN×L代表小波基，α∈CL×3代表小波系數(shù)，則=XWα，此時，問題模型表示為

此時，觀測矩陣表示為H=AW，通過多信號的正則化求解方法直接獲得小波系數(shù)α，根據(jù)公式X=Wα，可以求得森林區(qū)域的后向散射系數(shù)。

這種多觀測問題主要是求解具有相同支撐集的輸入向量。目前，其求解方法分為兩類，分別為陣列信號處理和壓縮感知。然而，兩種方法各有其局限性，所以本文采用結(jié)合了這兩種方法的CS-MUSIC算法進行重建。

3 CS-MUSIC算法

傳統(tǒng)的SAR層析成像是對同一場景進行重復(fù)觀測，獲得多幅2維復(fù)圖像，每一幅復(fù)圖像上的方位-距離像素點可視為本點高程向后向散射系數(shù)x(a,r,v)的傅里葉變換在某點的取值。若有多幅單視復(fù)圖像，即有了x(a,r,v)的傅里葉變換在多點的取值，則可利用頻譜估計的方法進行重建。

但是，傳統(tǒng)的頻譜估計方法信號序列長度受限，所能達(dá)到的高程向分辨率受限于高程向合成孔徑的大小，所以高程向分辨率相對較低。如果要提高高程向分辨率就必須增加觀測次數(shù)，來增加信號序列的長度，因此就需要多次的航過。

為了減少航過數(shù)，利用多極化通道信號之間的相關(guān)性，并引入了處理多觀測向量(MMV)問題的CS-MUSIC算法。CS-MUSIC算法將所要求得的支撐集分為兩部分，其中，k-3個支撐集元素用壓縮感知的方法進行恢復(fù)，其余的3個元素通過一種廣義的MUSIC準(zhǔn)則求得。其具體步驟如下[11]：

(1)利用HH, HV以及VV通道的信號，構(gòu)建多信號的像素點矩陣；

(2)確定像素點矩陣的列秩r =3，預(yù)估場景的稀疏度k，設(shè)置像素點矩陣和支撐集初值分別為

(3)利用同步的正交基追蹤(SOMP)[13]算法求出k-3個場景支撐集的索引值，更新支撐集；其中，該步驟SOMP算法又可以分為以下子步驟：

(a)J次迭代之后，更新支撐集和殘差信號

其中，I為單位矩陣，lj為第j次迭代求得的索引值，PSJ為回波信號在alj上的正交投影，為觀測矩陣A的第lj列。

(b)選擇1Jl+，滿足

更新索引集，使

(c) 確定迭代終止條件，當(dāng)J+1=k-3 時，終止迭代，并令I(lǐng)k-3=SJ+1，否則，繼續(xù)步驟(a)和(b) 。

(4)利用廣義MUSIC算法估計出其余的3個支撐集元素，并將其與步驟(3)獲得的支撐集進行合并，得到我們所需要的完整支撐集。其中，該步驟又可以分為以下子步驟：

(a)求得步驟(3)之后殘差信號的自相關(guān)矩陣，根據(jù)目標(biāo)點個數(shù)劃分信號子空間U和噪聲子空間Q；

(b)構(gòu)造譜搜索函數(shù)

(c)將上一步搜索函數(shù)的結(jié)果進行升序排列，取前3個值對應(yīng)的索引，并將其與Ik的元素一起放入S中，則S為我們所要求得的支撐集。

CS-MUSIC算法結(jié)合了壓縮感知和譜估計兩種方法，將支撐集分為兩部分進行恢復(fù)，3k-個元素利用壓縮感知的方法進行恢復(fù)，其余的3個元素則采用譜估計的方法獲得。這種方法利用了壓縮感知在欠采樣的情況下可以進行精確恢復(fù)，同時也保證了MUSIC算法在已知目標(biāo)個數(shù)情況下的精確重建，大大減少了虛假峰值的出現(xiàn)概率。

4 實驗

以下通過仿真和真實數(shù)據(jù)實驗來驗證CS-MUSIC算法在森林區(qū)域極化SAR層析成像中的有效性和優(yōu)越性。主要從兩方面對算法性能進行評估：一是散射中心位置估計的準(zhǔn)確性；二是減少航過數(shù)時高程向后向散射系數(shù)恢復(fù)情況以及虛假目標(biāo)的出現(xiàn)概率。

4.1 仿真數(shù)據(jù)實驗

根據(jù)文獻[14,15]中介紹的森林場景體散射模型，給出森林區(qū)域在3個極化通道情況下的高程向能量分布情況，如圖1所示

仿真數(shù)據(jù)實驗1設(shè)計如下：選擇不同極化通道情況下的散射模型，取高程向80 m范圍內(nèi)散射體分布，其散射中心分別在10 m和50 m處，分別代表地面和樹冠散射中心，取10次航過，高程向均勻基線，基線間隔10 m，中心斜距5494 m，下視角45o，波長0.86 m，信噪比30 dB，選擇Daubechies小波對高程向進行稀疏化處理。分別在航過數(shù)為10和6的情況下，將IST(迭代軟閾值)[16]，MUSIC, CAPON[17]和CS-MUSIC算法的重建結(jié)果進行對比。

圖1 3個極化通道高程向歸一化能量分布

圖2 不同航過數(shù)時各種算法獲得的高程向后向散射能量分布與給定的仿真結(jié)果對比

仿真數(shù)據(jù)實驗的結(jié)果表明，CS-MUSIC算法可以實現(xiàn)森林區(qū)域的SAR層析成像，在航過數(shù)為10時，IST, CAPON和CS-MUSIC算法都能準(zhǔn)確的估計出10 m和50 m處散射中心的位置，并且樹冠區(qū)域后向散射系數(shù)與仿真結(jié)果基本相符，而MUSIC算法只能估計出10 m處散射中心的位置，不能估計出另一處散射中心；當(dāng)航過數(shù)減少為6時，IST和MUSIC算法均不能準(zhǔn)確估計出兩處散射中心位置，CAPON算法估計的散射中心位置發(fā)生偏移，并且估計的結(jié)果中出現(xiàn)很多虛假目標(biāo)，干擾了對散射中心的判斷，而CS-MUSIC算法在減少航過的情況下仍能準(zhǔn)確地估計出兩處散射中心的位置，并且后向散射系數(shù)變化趨勢與仿真結(jié)果基本相符。

4.2 實測P-波段SAR數(shù)據(jù)實驗

為了驗證CS-MUSIC算法的有效性，我們使用了10組已聚焦和校正的P波段SAR圖像數(shù)據(jù)，該數(shù)據(jù)是德國宇航局(DLR) 與瑞典國防研究局(Swedish Defence ResearchAgency, FOI)應(yīng)用機載E-SAR 系統(tǒng)于 2007年 3 月到 5 月在瑞典Remingstorp 獲取的，為歐空局(ESA)BioSAR 計劃的一部分，稱為 BioSAR 2007，主要樹種為挪威云杉(Picea abies)，蘇格蘭紅松(Pinus sylvestris)和樺樹(Betula spp)，圖3所示為該場景P-波段主影像的Pauli基。地形相對比較平坦，海拔從 120 m 到 145 m 之間，樹高從 10 m 到 30 m。P-波段數(shù)據(jù)水平基線間隔大約為 10 m，水平總基線大約為 70 m，飛機多次觀測幾乎處于相同的高度。其他參數(shù)如下：飛行高度大約4005 m，空間分辨率斜距約為3 m，方位向為 1 m，入射角由近距的 25°到遠(yuǎn)距的 55°變化，中心頻率為0.35 GHz。

圖3 P-波段主影像 Pauli 基

我們選取了距離向800 m處的距離門包含的1800個方位向連續(xù)點，進行高程向的后向散射系數(shù)重建。從Pauli圖中可以看到我們選擇的實驗區(qū)域(藍(lán)色矩形框標(biāo)注的距離門)以樹木為主，對于驗證CS-MUSIC算法在森林區(qū)域SAR層析成像中的應(yīng)用具有代表性。實驗中，IST算法和CS-MUSIC算法使用的小波仍然是消失矩為4的Daubechies小波，兩組實驗參數(shù)4種算法獲得的實驗結(jié)果，均進行窗口大小為64的方位向多視處理。如圖4所示，在航過數(shù)為10時，CS-MUSIC算法與迭代軟閾值(IST)算法、CAPON算法以及MUSIC算法高程向重建結(jié)果進行的對比。在航過數(shù)為10時，4種算法均可以將地面和樹冠區(qū)域區(qū)分開，相比于MUSIC和CAPON算法，CS-MUSIC算法可以獲得完整的地面信息，樹冠信息更加完整，同時出現(xiàn)的虛假目標(biāo)數(shù)目較少；與IST算法相比，CS-MUSIC算法很大程度上減少了虛假目標(biāo)的出現(xiàn)概率。圖5為航過數(shù)等于6時，CS-MUSIC算法與迭代軟閾值(IST)算法，CAPON算法以及MUSIC算法高程向重建結(jié)果進行的對比。當(dāng)航過數(shù)減少時，IST, MUSIC和CAPON算法已經(jīng)不能找到地面和冠層的散射中心位置，并且出現(xiàn)很多虛假目標(biāo)，這與前面的仿真結(jié)果一致。而CS-MUSIC算法仍然可以清晰地將地面與冠層區(qū)分開，并出現(xiàn)相對較少的虛假目標(biāo)。

5 結(jié)束語

本文將CS-MUSIC算法應(yīng)用于森林區(qū)域極化SAR層析成像，針對傳統(tǒng)的SAR層析成像需要較多航過數(shù)的問題，提出了應(yīng)用全極化SAR各極化通道信號之間的相關(guān)性構(gòu)建回波觀測向量矩陣，應(yīng)用小波基對高程向結(jié)構(gòu)進行稀疏表示，并采用多信號的CS-MUSIC算法進行高程向的后向散射系數(shù)重建，相比于常見的單信號SAR層析成像方法，該方法可以準(zhǔn)確地估計出地表和樹冠散射中心的位置，同時降低虛假目標(biāo)出現(xiàn)概率，在同等重建精度的情況下，降低森林區(qū)域極化SAR層析成像時所需的航過數(shù)。

圖4 航過數(shù)等于10時各種算法高程向重建結(jié)果(圖示為50×1800的高程-方位向2維平面，對方位向進行平滑濾波的結(jié)果)

圖5 航過數(shù)等于6時各種算法高程向重建結(jié)果(圖示為50×1800的高程-方位向2維平面，對方位向進行平滑濾波的結(jié)果)

致謝 感謝中國林業(yè)科學(xué)研究院資源信息研究所陳爾學(xué)教授提供的Biomass數(shù)據(jù)集(龍計劃項目第3期支持，ID10609)，幫助我們驗證了算法的實用性。

[1] 張福博, 劉梅. 基于頻域最小二乘APES 的非均勻多基線SAR 層析成像算法[J]. 電子與信息學(xué)報, 2012, 34(7): 1568-1573.

Zhang Fu-bo and Liu Mei. Uneven multi-baseline SAR tomography base on frequency domain least squares amplitude and phase estimation (APES)[J]. Journal of Electronics & Information Technology, 2012, 34(7): 1568-1573.

[2] 彭學(xué)明, 王彥平, 譚維賢, 等. 基于跨航向稀疏陣列的機載下視MIMO 3D-SAR三維成像算法[J]. 電子與信息學(xué)報, 2012, 34(4): 943-949.

Peng Xue-ming, Wang Yan-ping, Tan Wei-xian et al.. Airborne downward-looking MIMO 3D-SAR imaging algorithm based on cross-track thinned array[J]. Journal of Electronics & Information Technology, 2012, 34(4): 943-949.

[3] Zhang Bing-chen, Hong Wen, and Wu Yi-rong. Sparse microwave imaging: principles and applications[J]. SCIENCECHINA Information Sciences, 2012, 55(8): 1-33.

[4] Reigber A and Moreira A. First demonstration of airborne SAR tomography using multibaseline L-band data[J]. IEEE Transactions on Geoscience and Remote Sensing, 2000, 38(5): 2142-2152.

[5] Mariotti d’Alessandro M and Tebaldini S. Phenomenology of P-band scattering from a tropical forest through three-dimensional SAR tomography[J]. IEEE Geoscience and Remote Sensing Letters, 2012, 9(3): 442-446.

[6] Zhu Xiao-xiang and Bamler R. Super-resolution power and robustness of compressive sensing for spectral estimation with application to spaceborne tomographic SAR[J]. IEEE Transactions on Geoscience and Remote Sensing, 2012, 50(1): 247-258.

[7] Budillon A, Evangelista A, and Schirinzi G. Three-dimensional SAR focusing from multipass signals using compressive sampling[J]. IEEE Transactions on Geoscience and Remote Sensing, 2011, 49(1): 488-499.

[8] 王金峰, 皮亦鳴, 曹宗杰. 一種機載SAR 層析三維成像算法[J]. 電子與信息學(xué)報, 2010, 32(5): 1029-1033.

Wang Jin-feng, Pi Yi-ming, and Cao Zong-jie. An algorithm for airborne SAR tomography 3D imaging[J]. Journal of Electronics & Information Technology, 2010, 32(5): 1029-1033.

[9] Aguilera E, Nannini M, and Reigber A. Multisignal compressed sensing for polarimetric SAR tomography[J]. IEEE Geoscience and Remote Sensing Letters, 2012, 9(5): 871-875.

[10] Aguilera E, Nannini M, and Reigber A. A data-adaptive compressed sensing approach to polarimetric SAR tomography of forested areas[J]. IEEE Geoscience and Remote Sensing Letters, 2013, 10(3): 543-547.

[11] Kim J M, Lee O K, and Ye J C. Compressive MUSIC: revisiting the link between compressive sensing and array signal processing[J]. IEEE Transactions on Information Theory, 2012, 58(1): 278-301.

[12] Aguilera E, Nannini M, and Reigber A. Wavelet-based compressed sensing for SAR tomography of forested areas[J]. IEEE Transactions on Geoscience and Remote Sensing, 2013, 51(12): 5283-5295.

[13] Tropp J A, Gilbert A C, and Strauss M J. Algorithms for simultaneous sparse approximation. Part I: Greedy pursuit[J]. Signal Processing, 2006, 86(3): 572-588.

[14] Treuhaft R N and Siqueria P. Vertical structure of vegetated land surfaces from interferometric and polarimetric radar[J]. Radio Science, 2000, 35(1): 141-177.

[15] Cloude S R and Papathanassiou K P. A 3-stage inversion process for polarimetric SAR interferometry[J]. IEE Proceedings-Radar, Sonar and Navigation, 2003, 150(3): 125-134.

[16] Daubechies I, Defriese M, and Mol C D. An iterative thresholding algorithm for linear inverse problems with a sparsity constraint[J]. Communications on Pure and Applied Mathematics, 2004, 57(11): 1413-1457.

[17] Tebaldini S. Algebraic synthesis of forest scenarios from multibaseline PolInSAR data[J]. IEEE Transactions on Geoscience and Remote Sensing, 2009, 47(12): 4132-4142.

張冰塵： 男，1973年生，研究員，主要研究方向為微波遙感與雷達(dá)技術(shù)、稀疏信號處理.

王萬影： 女，1989年生，碩士生，研究方向為信號與信息處理、稀疏微波成像算法.

畢 輝： 男，1991年生，碩士生，研究方向為信號與信息處理、稀疏微波成像算法.

趙 曜： 男，1984年生，助理研究員，主要研究方向為信號與信息處理、稀疏微波成像算法.

洪 文： 女，1968年生，研究員，主要研究方向為合成孔徑雷達(dá)成像與系統(tǒng)及其應(yīng)用、極化/極化干涉合成孔徑雷達(dá)數(shù)據(jù)處理及應(yīng)用、3維微波成像新概念新體制新方法等.

Polarimetric SAR Tomography for Forested Areas Based on Compressive Multiple Signal Classification

Zhang Bing-chen①②Wang Wan-ying①②③Bi Hui①②③Zhao Yao①②Hong Wen①②

①(Science and Technology on Microwave Imaging Laboratory, Beijing 100190, China)
②(Institute of Electronics, Chinese Academy of Sciences, Beijing 100190, China)
③(University of Chinese Academy of Sciences, Beijing 100190, China)

This paper focuses on the polarimetric SAR tomography for forested areas based on compressive Multiple Signal Classification (MSC). First, full polarimetric SAR receives the reflected echo of the imaging area. Then, the signals from polarimetric channels are used to build multiple measurement vector model, and a wavelet basis is used in order to sparsely represent vertical structure. For achieving the measurement of forested area, the backscattering coefficients are reconstructed by Compressive Multiple Signal Classification (CMSC) algorithm. Simulated data from PolSARpro software and P-band data acquired by the E-SAR sensor of the German Aerospace Center validate that the method can effectively reduce the passes for SAR tomography and the probability of occurrence of spurious spikes under the same measurement accuracy.

Polarimetric SAR; Tomography; Compressive Multiple Signal Classification (CMSC); Wavelet basis

TN958

A

1009-5896(2015)03-0625-06

10.11999/JEIT140584

2014-05-06收到，2014-10-09改回

國家973計劃項目(2010CB731905)和中國科學(xué)院創(chuàng)新團隊國際合作伙伴計劃“先進微波探測與信息處理”資助課題

*通信作者： 王萬影 wangwanying0913@126.com