基于自適應(yīng)高程約束的TomoSAR三維成像

任子帥 張 照 高雨欣 郭 睿*

①(西北工業(yè)大學(xué)自動化學(xué)院 西安 710129)

②(西北工業(yè)大學(xué)電子信息學(xué)院 西安 710129)

1 引言

層析合成孔徑雷達成像(Synthetic Aperture Radar Tomography,TomoSAR)技術(shù)是干涉SAR技術(shù)的延伸，通過對同一目標(biāo)的多角度觀測，可以反演其在高程向(垂直于方位向和斜距向平面)上不同高度的散射值，從而還原真實的SAR成像三維場景[1–5]。近年來，隨著機載、星載SAR系統(tǒng)的不斷成熟，高質(zhì)量多基線SAR影像被成功獲取，使得TomoSAR技術(shù)對地物高度向結(jié)構(gòu)的重建成為可能。TomoSAR成像可解決由疊掩效應(yīng)引起的SAR圖像解譯模糊問題[6,7]，并在許多方面都有成功應(yīng)用，例如城市監(jiān)測、生物量估計、目標(biāo)檢測、冰川變化等[8,9]。TomoSAR成像技術(shù)因其獨特的高程向成像能力，在還原建筑三維結(jié)構(gòu)方面具有顯著優(yōu)勢。TomoSAR能夠沿高程向重建單個像元內(nèi)的散射剖面，并提取各個散射體的高程位置和強度信息，進而還原建筑整體的三維結(jié)構(gòu)并獲得場景范圍內(nèi)的高精度建筑三維模型[10]。2010年以來，國內(nèi)外學(xué)者對TomoSAR技術(shù)在還原建筑三維結(jié)構(gòu)方面的應(yīng)用進行了廣泛的研究并取得了進展[11–13]。

TomoSAR高程反演的傳統(tǒng)譜估計算法包括自適應(yīng)波束形成算法(Capon)[14]、多重信號分類算法(Multiple Signal Classification,MUSIC)[15]等。由于受到了高程瑞利分辨率的限制，傳統(tǒng)譜估計算法的改進空間有限。隨著稀疏重構(gòu)技術(shù)與算法概念的提出和發(fā)展，壓縮感知(Compressed Sensing,CS)[16]已經(jīng)成為TomoSAR高程反演的主流算法。壓縮感知能突破奈奎斯特采樣定理的限制，使用較少的非均勻觀測樣本恢復(fù)出原始稀疏信號。在數(shù)據(jù)集的獲取方面，主要包括陣列SAR和多航過SAR。中國科學(xué)院空天信息創(chuàng)新研究院發(fā)布了機載陣列干涉SAR系統(tǒng)在山西運城市和四川峨眉山市飛行試驗獲取的部分數(shù)據(jù)，兩個數(shù)據(jù)中都包含疊掩和陰影現(xiàn)象較嚴重的建筑區(qū)域[17]。德國宇航中心則利用TerraSAR-X多航過數(shù)據(jù)開展了大量的關(guān)于三維建筑重構(gòu)的研究[18,19]。

然而TomoSAR往往存在嚴重的因相位纏繞導(dǎo)致的高程模糊問題。在高程反演過程中，若高程搜索范圍超過最大不模糊高程，高程反演會出現(xiàn)周期性的峰值，產(chǎn)生偽影，影響成像結(jié)果。機載Tomo-SAR成像最大不模糊高程較小，高程模糊明顯[20]；多航過星載TomoSAR成像最大不模糊高程較大，但噪聲水平更高，會存在其他偽影。目前專門針對TomoSAR高程解模糊的相關(guān)研究較少，但高程模糊是TomoSAR處理中不可忽略的問題之一。文獻[21]利用圓跡SAR的全方位觀測、解疊掩和解高程模糊優(yōu)勢，進行無先驗?zāi)Ｐ偷腟AR三維成像；文獻[22]采用互質(zhì)采樣子陣列提高TomoSAR處理中的觀測并增大基線孔徑解決采樣不足引起的模糊問題，提高TomoSAR三維成像信噪比；文獻[23]針對山地區(qū)域的TomoSAR成像，利用三維點云分割重組的方法進行高程解模糊；文獻[24]中提到建筑高度超過最大不模糊高度時，傳統(tǒng)TomoSAR三維點云會存在嚴重模糊，通過對三維重建的解空間進行約束提高三維重建效果，同時結(jié)合圖像語義信息進行三維點云聚類，實現(xiàn)解模糊。

針對TomoSAR中的高程模糊本文提出一種自適應(yīng)調(diào)整高程搜索范圍的方法，該方法通過自適應(yīng)的高程采樣中心和搜索半徑設(shè)置，使得高程向的模糊峰值被隔離，從而改善高程模糊和偽影問題。在實驗中，針對機載陣列SAR數(shù)據(jù)和星載TerraSAR-X多航過數(shù)據(jù)，該方法能改善TomoSAR三維點云重建結(jié)果，并提出平均鄰域高度差和散射體完整度評價指標(biāo)驗證該方法的有效性和廣泛適用性。

2 TomoSAR成像原理

SAR影像每個像元的復(fù)數(shù)值是同一距離的各個高度目標(biāo)后向散射信號在距離-方位成像平面的投影，換句話說，散射信號沿高程向的積分組成了像元值[25]。TomoSAR成像幾何原理如圖1，不同航過或陣列天線在垂直于距離向的方向上的距離稱為基線，選取一個基線作為參考基線，第n景SAR影像對應(yīng)的基線記為bn，?b為基線跨度。目標(biāo)在距離向與雷達天線的距離為斜距，記為r。雷達天線運動方向為方位向，記為x。與距離向和方位向正交的方向為高程向，高程向距離記為s。設(shè)參考基線對應(yīng)的SAR天線S1獲取的影像為主影像，天線Sn對應(yīng)SAR影像的任一像素gn為[26]

圖1 TomoSAR三維成像幾何模型Fig.1 TomoSAR 3D geometric model

其中，ξn=2bn/(λr)為高程向頻率，λ代表波長。γ(s)為沿著高程向分布的后向散射強度分量，也是需要計算恢復(fù)的目標(biāo)變量。從不同航過獲取的N景SAR影像中，可以構(gòu)造一個高程向的合成孔徑。沿高程向的采樣范圍 ?s進行L次均勻離散采樣，則式(1)可以轉(zhuǎn)化為

考慮N景SAR影像，將每幅SAR影像上同名像素點的復(fù)數(shù)值作為向量，可以得到：

其中，ε為噪聲分量，Gnl代表觀測矩陣元素：

此時，式(3)可簡寫為[6,27]

TomoSAR處理是根據(jù)預(yù)處理后的N景SAR影像中復(fù)數(shù)值向量g和構(gòu)造的測量矩陣G，逆向求出沿高程向分布的目標(biāo)散射強度γ，進而可以根據(jù)散射強度峰值、峰值個數(shù)和峰值位置計算散射體的強度、數(shù)目和高程向位置。

完成圖像配準、通道誤差補償、大氣相位誤差補償?shù)鹊腟AR數(shù)據(jù)預(yù)處理后，可以依據(jù)成像模型(5)通過Capon算法[14,15]、MUSIC算法[16]等譜估計算法或CS算法解算出高程向分布的散射體，實現(xiàn)Tomo-SAR三維成像。CS算法由Zhu和Bamler[11]首次引入TomoSAR，并利用TerraSAR-X聚束模式下的高分辨率數(shù)據(jù)驗證該方法在欠采樣下的超分辨率成像能力，CS算法得到的稀疏層析譜中可以分離高度相近的疊掩散射體。在實際應(yīng)用過程中普遍使用L1范數(shù)的凸優(yōu)化模型[28,29]：

其中，λk為嶺估計參數(shù)，與噪聲等級有關(guān)。相比于傳統(tǒng)譜估計算法，CS算法重建的后向散射剖面不是連續(xù)的而是稀疏的，而且無需進行多視處理即可實現(xiàn)超分辨率成像。

3 高程模糊

3.1 高程模糊分析

由式(2)可知，SAR影像任一像素的相位部分包含散射剖面的高程信息：

式(7)相位項被纏繞在2π里并以2kπ為周期(k ∈Z)，導(dǎo)致在利用相位解算高程散射剖面時出現(xiàn)偽影現(xiàn)象之一的高程模糊。且以下特點：高程向散射強度呈現(xiàn)周期性重復(fù)譜峰，其中一個譜峰為真實譜峰，其余為模糊譜峰；單個像素存在多個散射體的情況下(即疊掩散射體)，不同散射體的真實譜峰和模糊譜峰會相互影響；三維點云成像結(jié)果中出現(xiàn)分層的偽目標(biāo)，如圖2。

圖2 高程向的偽目標(biāo)示意圖Fig.2 Pseudo target in the elevation

圖2中偽目標(biāo)與實像之間的高程距離在本文中作為高程模糊間隔或者最大不模糊高程。根據(jù)式(7)，可以得到高程向模糊的間隔周期[30]：

其中，samb與載波波長λ、斜距r以及垂直基線bn有關(guān)。并且式(7)的相位項存在以下關(guān)系：

根據(jù)式(9)，由于高程模糊間隔和譜峰周期性的存在，在TomoSAR提取散射體高程信息時會遇到3種問題：

(1) 多散射體高程距離差?s′小于高程模糊間隔，并且小于瑞利分辨率ρs的情況下(即?s′<ρs)，一般的譜估計方法難以分離疊掩散射體。對于該問題可以利用CS算法來突破瑞利分辨率的限制。

(2) 多散射體高程距離差?s′小于高程模糊間隔而且大于等于瑞利分辨率的情況下，但高程采樣的搜索范圍大于高程模糊間隔時(即ρs≤?s′samb)，高程反演算法會重構(gòu)出模糊峰值。在通過模型選擇判斷散射體數(shù)目之后，該問題可能會導(dǎo)致將模糊峰值位置錯誤地估計為散射體位置的情況，從而忽視了真實散射體的位置。針對該問題，需要對高程采樣范圍進行優(yōu)化調(diào)整。高程向采樣范圍 ?s通常需要滿足以下限制[11]：

其中，?bmin為最小基線間隔。但這是一種浮動的約束，只提供了采樣范圍的最大跨度，沒有提供一個采樣中心。在面對實際數(shù)據(jù)時，特別是SAR天線相位中心高度較低的情況(如機載數(shù)據(jù))，建筑高度非常接近最大不模糊高度，該約束條件并不能完全滿足一個準確的高程采樣范圍。本文基于實際應(yīng)用過程中獲得的經(jīng)驗，提出一種自適應(yīng)高程采樣范圍定義來改善層析三維點云中的高程模糊現(xiàn)象。

(3) 當(dāng)多散射體的高程距離差大于等于一個高程模糊間隔時(即?s′ ≥samb)，不同散射體的高程散射峰值會相互影響，給估計散射體正確的高程位置帶來困難。針對該問題，一般的解決方式是優(yōu)化基線分布，使得高程模糊間隔遠大于目標(biāo)場景的建筑高度；或者采用互質(zhì)基線設(shè)計，在不影響瑞利分辨率和基線數(shù)量的情況下改善高程模糊。從式(10)中可以發(fā)現(xiàn)高程模糊間隔可以直接通過調(diào)整垂直基線而發(fā)生改變。優(yōu)化基線分布同時對于SAR傳感器的控制和分布有一定的要求，而且對于機載和星載等常用的SAR成像平臺，散射體高程差不小于高程模糊間隔的問題一般不會存在，在本文不考慮。

3.2 高程模糊及其他偽影仿真

本節(jié)對TomoSAR單個像素的層析結(jié)果進行仿真分析，表1為仿真參數(shù)，針對基線短、基線均勻、斜距小的情況。根據(jù)表1中的仿真參數(shù)和均勻等距的基線分布進行TomoSAR仿真實驗。實驗假設(shè)高程向有兩個散射體，散射體的歸一化散射強度為1.0和0.6，分別設(shè)置5 m,60 m和120m的散射體間距，并加入了信噪比(Signal to Noise Ratio,SNR)為15 dB的高斯噪聲。采用CS算法對單個像素的高程向散射剖面進行比較和分析。根據(jù)表1中的參數(shù)和式(8)進行計算，可以得到仿真數(shù)據(jù)的高程模糊間隔為100 m；根據(jù)瑞利分辨率的計算公式ρs=λr/(2?bmax)，可得仿真數(shù)據(jù)的瑞利分辨率為10 m。3種散射體間距的高程反演仿真結(jié)果如圖3所示。

表1 仿真實驗1參數(shù)Tab.1 Experimental parameters for the first simulation data

圖3 仿真數(shù)據(jù)1的高程散射剖面反演結(jié)果Fig.3 Reconstruction results of elevation reflectivity profile of the first simulation data

在圖3(a)中，散射體5 m的間隔已經(jīng)小于仿真數(shù)據(jù)的瑞利分辨率，兩個散射體得到了識別，模糊位置和真實散射體位置相距較大，分離比較方便。在圖3(b)中，在散射體間隔60 m的情況下，兩個散射體位置處均出現(xiàn)了峰值。但是在高程向采樣搜索范圍大于高程模糊間隔的情況下，會產(chǎn)生嚴重的高程模糊現(xiàn)象，每一個散射體都在間隔100 m的高程位置產(chǎn)生模糊峰值，模糊峰值與真實峰值的高度相當(dāng)，而且兩者的高程距離較近。在圖3(c)中，散射體間隔120 m的情況下，散射體間距大于高程模糊間隔，真實散射體峰值和模糊峰值相互纏繞在一起，難以分辨。因此，如果高程采樣的搜索范圍設(shè)置不合適，會導(dǎo)致反演的高程散射剖面出現(xiàn)多余的峰值，導(dǎo)致錯誤地將局部模糊峰值當(dāng)作目標(biāo)散射體峰值，從而造成散射體高程位置估計錯誤。

此外，本文也對基線長、基線不均勻、斜距大的TomoSAR單像素層析進行仿真，表2為其仿真實驗參數(shù)。該仿真實驗假設(shè)中假設(shè)高程間距為120 m的兩個散射體，歸一化散射強度分別為1.0和0.6。該仿真數(shù)據(jù)的理論最大不模糊高程為428.6 m，高程瑞利分辨率為20 m。星載多航過SAR仿真數(shù)據(jù)的基線采用不均勻排列方式，最小基線間隔為21 m，最大基線跨度為450 m。由于最大不模糊高程較大，不易產(chǎn)生非常明顯的高程模糊。但是由于不均勻的基線分布和因長距離探測導(dǎo)致的較低SNR，TomoSAR處理可能也會產(chǎn)生一些不規(guī)則的偽影。在仿真信號中加入了SNR為15 dB,10 dB和5 dB的高斯噪聲，以模擬噪聲影響。不同SNR下高程散射剖面重構(gòu)結(jié)果如圖4所示。

表2 仿真實驗2參數(shù)Tab.2 Experimental parameters for the second simulation data

圖4 仿真數(shù)據(jù)2的高程散射剖面反演結(jié)果Fig.4 Reconstruction results of elevation reflectivity profile of the second simulation data

在圖4(a)中，SNR=15 dB時可以明顯看出高程向有兩個明顯峰值，與真實散射體位置幾乎重合。在圖4(b)中，SNR=10 dB高程散射剖面出現(xiàn)多個較大的峰值，其中最高的兩個峰值在真實散射體的位置。因為不模糊高程較大，仿真數(shù)據(jù)2中的偽影現(xiàn)象有一定程度緩解，在較高信噪比情況下沒有出現(xiàn)大于真實散射體峰值的情況。在圖4(c)，隨著信噪比降至SNR=5 dB，偽影現(xiàn)象已經(jīng)非常明顯，最高的兩個峰值與真實散射體位置出現(xiàn)了偏差，表示噪聲較大會加劇TomoSAR的偽影現(xiàn)象，因此也需要在TomoSAR處理過程中優(yōu)化高程采樣的搜索范圍，以減少偽影的影響。

4 自適應(yīng)高程搜索范圍設(shè)置方法

根據(jù)上文的分析，在TomoSAR三維成像過程中的合適的高程搜索范圍約束是有必要的。高程向搜索范圍需要滿足式(10)才能保證在層析模型的相位項被約束在一個周期內(nèi)。但是式(10)的約束只提供了一個最寬松的條件，而且沒有確定約束在高程向上的絕對位置，在實現(xiàn)中通常根據(jù)經(jīng)驗確定，得到的三維成像結(jié)果可能會存在高程模糊。因此，本文提出一種自適應(yīng)高程搜索范圍方法，以改善TomoSAR結(jié)果中的高程模糊現(xiàn)象。

自適應(yīng)高程搜索范圍方法結(jié)合場景高度等信息，確定高程搜索的采樣中心和搜索半徑，并根據(jù)每個像素的斜距方位索引，自適應(yīng)設(shè)置高程搜索范圍。方法具體實現(xiàn)如下：

(1) 場景最大高度預(yù)估計。一般TomoSAR中依據(jù)式(10)設(shè)置場景的高程范圍，本文結(jié)合SAR影像中的疊掩信息[31]，獲取場景中的最高目標(biāo)高度估計值，并將此作為場景的最大高度預(yù)估計hmax。

(2) 地面最低高度預(yù)估計。一般選擇SAR影像中的一個地面像素點作為最低高度參考像素，即認為所選像素中散射體高度不高于其他散射體高度?？梢栽谛本嘞?高程向的坐標(biāo)系中過參考點構(gòu)建一條斜率為cotθ的直線，作為地面最低高度的約束。但是在實際TomoSAR處理過程中，地形地貌起伏等的影響會使得參考像素的高度與場景真實的最低高度不完全一致，導(dǎo)致高程搜索范圍的約束下邊界并不準確。因此，在實際應(yīng)用中地面最低高度的估計還可以通過以下方式獲得：先利用層析算法估計所有散射體的高度信息，重建未解模糊的三維點云；再將點云投影到斜距向-高程向的坐標(biāo)系中；計算投影的點云密度，將密度較大的點作為地面點；根據(jù)與零高度線的距離確定非模糊地面點；提取地面點的最低點作為地面最低高度的估計hmin。

(3) 構(gòu)建高程采樣中心線性函數(shù)。設(shè)SAR影像斜距方向的坐標(biāo)索引為i(i=1,2,...,I)，I為斜距向像素數(shù)；方位向上的坐標(biāo)索引為j(j=1,2,...,J)，J為方位向像素數(shù)；近端斜距為rnear，距離向采樣間隔為ρr。選擇SAR影像中某一像素點作為零高程參考像素，參考像元的坐標(biāo)為(iref,jref)，參考點的斜距為rref，高程位置s與高度h的對應(yīng)關(guān)系為

場景中最大的建筑的高度值為hmax，再利用式(11)轉(zhuǎn)換成高程值，將其設(shè)置為高程向搜索范圍上邊界。場景中地面最低高度為hmin，也可以轉(zhuǎn)換成高程值并設(shè)置為高程向搜索范圍下邊界。根據(jù)高程向搜索范圍上下邊界構(gòu)建線性函數(shù)，要求與上下邊界等間距，斜率為cotθ。將線性函數(shù)上的點作為每個像素的高程采樣中心，于是每個像素的高程采樣中心為

(4) 計算搜索半徑并設(shè)置各像素的高程搜索范圍。高程向搜索范圍的搜索半徑為

其中，?h=hmax-hmin為高度差。因此可以根據(jù)每個像素的坐標(biāo)索引自適應(yīng)地確定高程向搜索范圍：

(5) 高程反演。最后，利用CS算法逐像素在高程搜索范圍內(nèi)反演散射剖面，在高程反演結(jié)果提取最高的M個峰值(M為像素內(nèi)疊掩散射體個數(shù))。最后還原散射體的三維信息，得到三維點云。

該方法的流程圖如圖5所示。

圖5 自適應(yīng)高程搜索范圍方法流程圖Fig.5 The flowchart of adaptive elevation search range method (blue: The boundary before constraints;red: Boundaries after constraints;yellow: Elevation sampling center)

5 實測數(shù)據(jù)驗證與分析

5.1 機載陣列SAR數(shù)據(jù)

首先采用中國科學(xué)院空天信息創(chuàng)新研究院在《雷達學(xué)報》發(fā)布的峨眉機載陣列SAR數(shù)據(jù)集[17]進行驗證，該數(shù)據(jù)集的相關(guān)參數(shù)如表3所示，Google Earth光學(xué)圖和SAR幅度圖如圖6所示。數(shù)據(jù)集由11景SAR影像組成，空間基線均勻分布。按照參數(shù)計算得到理論高程瑞利分辨率為18.7 m，高程模糊間隔為189.4 m，則高度上的模糊間隔為99.2 m。場景高度分布在[–22 m,65 m]，場景高度差約為87 m，場景高度差接近高度向上的模糊間隔，由此引起的高程模糊可能會影響三維成像結(jié)果。

表3 峨眉數(shù)據(jù)參數(shù)Tab.3 Parameters of Emei data

圖6 峨眉區(qū)域的Google Earth光學(xué)圖和SAR幅度圖Fig.6 Google Earth optical map and SAR amplitude map of Emei area

圖7(a)為未進行高程搜索范圍自適應(yīng)約束的三維成像結(jié)果。因為高程搜索范圍超過了高程模糊間隔，該三維點云存在高程模糊，部分散射體的高程向位置被錯誤估計，不同高度層出現(xiàn)相似的建筑點云。圖7(b)為采用自適應(yīng)調(diào)整高程搜索范圍方法的三維成像結(jié)果，高程搜索范圍的大致邊界如圖7(a)紅色虛線所示。在每個像元的散射體個數(shù)M已知的情況下，首先選擇高程散射剖面中最高的M個峰值。然后在采樣中心與搜索半徑的約束下，每個像素的高程搜索范圍被約束到一個更合適的區(qū)間。最后通過高程反演和三維信息還原得到了解模糊的三維點云。本方法既保證了真實點云的正常顯示，也完成了高程模糊的大部分濾除。圖7(c)為將斜距-方位-高程坐標(biāo)系轉(zhuǎn)換為地距-方位-高度坐標(biāo)系的結(jié)果，圖7(d)為峨眉數(shù)據(jù)全部區(qū)域的三維點云結(jié)果。峨眉數(shù)據(jù)的三維成像結(jié)果顯示最高建筑的高度約為65 m，與前期研究[32]和高度預(yù)估計相符，建筑結(jié)構(gòu)清晰，高程模糊現(xiàn)象被大幅度改善。

為了進一步分析高程模糊的影響以及本文方法在高程解模糊上的效果，實驗中取SAR圖像兩條沿斜距向的剖線進行分析，并標(biāo)號為剖線1和剖線2，如圖8(a)所示。圖8(b)為所選剖線在光學(xué)圖中的大致位置。剖線1包含有一個建筑物的散射體，而剖線2包含有多個建筑物的散射體。圖8(c)和圖8(d)分別為所選剖線1和所選剖線2的歸一化高程向散射強度譜，可以看到圖中存在結(jié)構(gòu)相似的高亮區(qū)域分布在不同的分層。其中紅色橢圓內(nèi)的區(qū)域是真實散射體的散射強度譜，而其他高亮區(qū)域代表高程模糊。在某些斜距向像素中，高程模糊的散射強度甚至高于真實散射體的散射強度，易造成散射體高程向位置的錯誤估計。例如在圖8(d)中，上下兩個綠色圓框分別位于模糊散射體位置和真實散射體位置，上方的綠色圓框內(nèi)的散射強度明顯高于下方的綠色圓框內(nèi)的散射強度。取圖8(c)中的紅色箭頭位置處對應(yīng)的單個像素，得到該像素處高程剖面圖如圖8(e)和圖8(f)所示。其中藍色曲線為該像素的高程向散射強度譜，散射強度譜的峰值位置可以大概反映散射體的實際位置和高程模糊的位置。紅色曲線是基于CS算法的高程稀疏反演。圖8(e)是未采用自適應(yīng)高程范圍約束的高程反演結(jié)果，根據(jù)圖8(c)可以確定該像素的真實散射體的大概位置和高程模糊的位置，并在圖8(e)中標(biāo)注出。圖8(e)中左側(cè)高程模糊的散射強度甚至高于真實散射體的散射強度，且該像素的CS高程反演也錯誤估計了真實散射體位置。圖8(f)為采用本文方法進行高程向搜索范圍約束后的高程反演結(jié)果，其中黑色虛線為該像素的高程向搜索范圍的邊界?？梢钥吹?，CS反演的估計結(jié)果正確地反映了真實散射體的位置。

圖8 峨眉數(shù)據(jù)的散射剖面反演結(jié)果Fig.8 Inversion results of scattering profile from Emei data

為了評估本文方法對高程模糊的改善，參考鄰域均值這一概念，提出“平均鄰域高度差”這一定義來衡量三維點云在高度向上的集中度和連續(xù)性。平均鄰域高度差被定義如下：

其中，MS為估計的散射體總數(shù)，h(i,j)為像素(i,j)中所有疊掩散射體的高度之和，h(i,j,m)為對應(yīng)像素(i,j)中第m個散射體的高度，M(i,j)為對應(yīng)像素(i,j)中疊掩散射體的數(shù)目，MN(i,j)為對應(yīng)像素(i,j)相鄰像素的疊掩散射體數(shù)目。平均鄰域高度差?hE可以概括為相鄰像素包含的散射體高度差的均值，數(shù)值越小表示三維點云在高度向上的集中度和連續(xù)性越好。表4列出未采用和采用自適應(yīng)高程搜索范圍的三維點云平均鄰域高度差，可以看出本文方法較大地提高點云分布的集中度和連續(xù)性。點云在高度向上分布更集中，則高程模糊產(chǎn)生的相鄰像素高度突變現(xiàn)象更少。

表4 峨眉數(shù)據(jù)的平均鄰域高度差(m)Tab.4 ?hE of Emei data (m)

為了進一步分析高度預(yù)估計對本文方法和三維點云的影響，本文還定義了“散射體完整度”來反映三維點云的完整性。散射體完整度定義如下：

其中，Mref為參考散射體數(shù)目。散射體完整度C為一定條件下得到的三維點云中散射體的數(shù)目占參考散射體數(shù)目的比值。散射體完整度數(shù)值越高表示三維點云的散射體數(shù)量越多，建筑結(jié)構(gòu)展現(xiàn)得越完整。以1 m為高度網(wǎng)格，將高度預(yù)估計結(jié)果hmax和hmin附近的高度值分別作為本文方法中高程采樣中心和搜索半徑的參數(shù)，得到三維點云，并計算平均鄰域高度差?hE和散射體完整度C，結(jié)果如圖9所示。圖9(a)是場景最高建筑高度預(yù)估計結(jié)果hmax浮動時的變化，將hmax=70 m時的點云散射體數(shù)量作為參考散射體數(shù)量。當(dāng)hmax小于65 m時，散射體完整度急劇下降；當(dāng)hmax大于65 m時，散射體完整度幾乎不變，但是平均鄰域高度差上升，點云集中度和連續(xù)性變差。這表明了最高建筑高度預(yù)估計結(jié)果hmax=65 m的準確性。圖9(b)是場景地面點最低高度預(yù)估計結(jié)果hmin浮動時的變化，將hmin=-27 m時的點云散射體數(shù)量作為參考散射體數(shù)量。當(dāng)hmin大于–22 m時，平均鄰域高度差緩慢下降，但是散射體完整度急劇下降；當(dāng)hmax小于–22 m時，散射體完整度幾乎不變，但是平均鄰域高度差急劇上升。這表明了場景地面點最低高度預(yù)估計結(jié)果hmin=-22 m的準確性。

圖9 不同高度估計下的平均鄰域高度差和散射體完整度Fig.9 ?hE and C at different height estimation

5.2 TerraSAR-X數(shù)據(jù)

本文同時采用高分辨率TerraSAR-X數(shù)據(jù)集進行驗證。衛(wèi)星軌道為升軌模式，工作模式為條帶模式，數(shù)據(jù)以多航過方式獲得。此數(shù)據(jù)集包含巴塞羅那城區(qū)23景高分辨率影像數(shù)據(jù)，該數(shù)據(jù)集的入射角、基線、斜距等具體參數(shù)詳見表5，Google Earth光學(xué)圖和SAR幅度圖如圖10所示。巴塞羅那數(shù)據(jù)的理論高程瑞利分辨率為39.1 m，最大不模糊高程為1207.4 m。在第3節(jié)中已說明星載多航過SAR雖然具有較大的高程模糊間隔，但高程散射剖面仍有許多旁瓣影響成像效果，低信噪比造成的偽影現(xiàn)象仍然不能忽視，如3.2節(jié)中仿真所示，在對巴塞羅那數(shù)據(jù)進行三維成像時也有必要采用自適應(yīng)高程搜索范圍方法。

表5 巴塞羅那數(shù)據(jù)參數(shù)Tab.5 Parameters of Barcelona data

圖10 巴塞羅那區(qū)域的Google Earth光學(xué)圖和SAR幅度圖Fig.10 Google Earth optical map and SAR amplitude map of Barcelona area

圖11(a)為未進行自適應(yīng)高程搜索范圍約束的建筑三維重建結(jié)果，圖11(b)為采用了自適應(yīng)調(diào)整高程搜索范圍方法的建筑三維重建結(jié)果?？梢钥吹綀D11(a)中的較高區(qū)域和較低區(qū)域存在偽影(紅圈內(nèi))，而圖11(b)中偽影點云被重建到了正確位置，驗證了本文算法在星載多航過SAR數(shù)據(jù)層析成像中的適用性。

圖11 巴塞羅那數(shù)據(jù)的三維重建結(jié)果Fig.11 3D reconstruction results of Barcelona data

選擇圖11(a)紅圈內(nèi)某一散射體所在的SAR影像像元，得到該像素高程向剖面如圖12(a)所示，并將高程映射為高度。在圖12(b)中，藍色曲線為未采用自適應(yīng)高程搜索范圍的高度反演結(jié)果，所選像素的散射體位置被估計在高度215 m左右，但根據(jù)實際經(jīng)驗和光學(xué)圖可以判斷該高度被錯誤估計。紅色曲線為采用本文自適應(yīng)高程搜索范圍的高度反演結(jié)果，黑色虛線為所選像素的高程搜索范圍邊界?？煽吹?，所選像素的散射體位置被估計在高度50 m附近，Google Earth光學(xué)圖顯示所選像素的所在實際地理位置為中等高度建筑，根據(jù)實際經(jīng)驗可以判斷該高度估計得到了改進。

進一步采用平均鄰域高度差來定量判斷本文方法在改善點云偽影方面的性能。表6列出未采用自適應(yīng)高程搜索范圍方法的三維點云結(jié)果和采用本文自適應(yīng)高程搜索范圍的三維點云結(jié)果的平均鄰域高度差?hE。本文方法得到的三維點云的?hE更小，表明點云的偽影現(xiàn)象得到了改善。

6 結(jié)語

本文研究了TomoSAR中高程模糊現(xiàn)象產(chǎn)生的原因和特性，針對這一問題提出一種自適應(yīng)調(diào)整高程搜索范圍的方法。該方法結(jié)合場景的高度預(yù)估計，構(gòu)建以零高程像元為參考的高程采樣中心線性函數(shù)，從而實現(xiàn)每個像素的高程搜索中心自適應(yīng)變化，并估計高程向搜索半徑，以此得到根據(jù)每個像素索引自適應(yīng)變化的高程向搜索范圍約束，最后通過高程反演和三維信息還原得到場景三維點云?；跈C載及星載實測數(shù)據(jù)的三維重建實驗中，本文所提方法能夠?qū)C載平臺飛行高度下的高程模糊問題進行改善，除此，該方法能改善星載平臺下的偽影現(xiàn)象。通過定義的平均鄰域高度差和散射體完整度指標(biāo)，本文所提方法在TomoSAR三維點云重建中的性能得以驗證，但對于更低高度的成像平臺，本文方法的性能需進一步的討論。

致謝感謝中國科學(xué)院空天信息創(chuàng)新研究院將機載陣列干涉SAR數(shù)據(jù)公開于雷達學(xué)報，促使了本文的研究。