劉 鵬，李真芳，樓良盛，楊偉明，王 震

1. 西安電子科技大學(xué)雷達(dá)信號(hào)處理國(guó)家重點(diǎn)實(shí)驗(yàn)室，陜西 西安 710071； 2. 地理空間信息國(guó)家重點(diǎn)實(shí)驗(yàn)室，陜西 西安 710054； 3. 西安測(cè)繪研究所，陜西 西安 710054

干涉合成孔徑雷達(dá)(interferometric synthetic aperture radar，InSAR)技術(shù)得益于其獲取目標(biāo)高程信息的能力，被廣泛應(yīng)用于地形測(cè)繪、海洋監(jiān)測(cè)、災(zāi)害監(jiān)測(cè)與預(yù)警等領(lǐng)域[1-4]。由于合成孔徑雷達(dá)屬于主動(dòng)輻射測(cè)距成像模式的一種技術(shù)，使得疊掩和陰影成為SAR圖像中普遍存在的現(xiàn)象。在InSAR信號(hào)處理中，疊掩和陰影區(qū)域的相位奇異性導(dǎo)致相位無(wú)法正確濾波和相位解纏，錯(cuò)誤的相位信息導(dǎo)致相鄰區(qū)域的相位展開(kāi)錯(cuò)誤，并造成高程反演出現(xiàn)錯(cuò)誤。因此需要準(zhǔn)確判斷并標(biāo)記SAR圖像的陰影區(qū)和疊掩區(qū)，避免該區(qū)域影響整個(gè)相位解纏過(guò)程。

近年來(lái)，許多學(xué)者對(duì)InSAR中陰影和疊掩的檢測(cè)方法進(jìn)行了相關(guān)研究，并取得了一定成果，如基于幅度的閾值分割方法[5]、基于相干系數(shù)的閾值分割方法[6]、基于邊緣特征和模糊理論的方法[7]、基于局部頻率估計(jì)的方法[8]、基于特征值分解的方法[9]和基于恒虛警(CFAR)檢測(cè)的方法[10]等。其中，局部頻率估計(jì)和特征值分解更為有效，但其都有各自的局限性。一方面，上述方法都是逐像素進(jìn)行判別，SAR圖像實(shí)際上反映的是場(chǎng)景對(duì)電磁波的雙站散射能力，目標(biāo)的幾何和物理特性通過(guò)影響場(chǎng)景電磁雙站散射而被SAR系統(tǒng)捕獲并成像，SAR圖像上每個(gè)分辨單元是由完全隨機(jī)分布的許多散射點(diǎn)矢量疊加組成，在SAR圖像上形成相干斑噪聲，進(jìn)而使得簡(jiǎn)單閾值分割方法失效，逐像素的判別結(jié)果較為離散，即使做形態(tài)學(xué)處理，也很難反映SAR數(shù)據(jù)的紋理信息；另一方面，局部頻率估計(jì)需要對(duì)每個(gè)像素周?chē)付ù翱谶M(jìn)行二維FFT操作，對(duì)SAR數(shù)據(jù)逐像素進(jìn)行判別計(jì)算量較大，特征值分解方法需要獲取同一區(qū)域的多幅SAR圖像，通過(guò)構(gòu)造圖像陣列的方式判別疊掩和陰影，無(wú)法應(yīng)用于大部分工程需求，恒虛警檢測(cè)(constant false alarm，CFAR)基于目標(biāo)與背景對(duì)比強(qiáng)烈的準(zhǔn)則，需要估計(jì)目標(biāo)的分布函數(shù)，對(duì)于較大的疊掩和陰影區(qū)域檢測(cè)性能不佳。

針對(duì)逐像素檢測(cè)方法在準(zhǔn)確度和速度上的不足，考慮到在SAR圖像中屬于同一地物同等坡度的相鄰像素具有相似信息，疊掩或陰影區(qū)域具有散射信息近似相等且距離上相近的特點(diǎn)，如果能夠?qū)⒕哂邢噜徢蚁嗨菩畔⒌南袼亟M合成一個(gè)基元塊，則可以大大降低后續(xù)圖像處理任務(wù)的復(fù)雜性?；诖?，本文結(jié)合超像素分割方法提出一種基于多信息融合的超像素檢測(cè)算法，實(shí)現(xiàn)了快速且準(zhǔn)確的疊掩和陰影區(qū)域檢測(cè)。近年來(lái)，常用的超像素分割方法是基于圖論和梯度下降聚類(lèi)的方法。基于圖論的方法具體包括基于隨機(jī)游走的熵率的方法[11]、偽布爾優(yōu)化的方法[12]和基于圖的貪心聚類(lèi)算法[13]。基于梯度下降聚類(lèi)的超像素分割包括簡(jiǎn)單線性迭代聚類(lèi)(SLIC)方法[14]、迭代邊緣細(xì)化方法[15]、多尺度分割方法[16]和分水嶺算法[17]。這些方法由于其處理速度快、邊界擬合度較好，被廣泛應(yīng)用于光學(xué)圖像分割中。其中，SLIC方法由于在遵守圖像邊界的能力、速度、內(nèi)存效率以及分割性能上的優(yōu)越性而被廣泛使用[18]。然而，由于SAR圖像固有的相干斑噪聲效應(yīng)[19]，以及相對(duì)于光學(xué)圖像而言SAR圖像具有更多像素?cái)?shù)量，通常達(dá)到106～108數(shù)量級(jí)，使得傳統(tǒng)SLIC在超像素分割方面表現(xiàn)不佳。近年來(lái)多尺度分割和SLIC方法在SAR圖像分割領(lǐng)域得到了廣泛的應(yīng)用。盡管過(guò)去已提出了許多的多尺度分割方法，但嵌入在eCognition軟件中的多尺度分割(multiscale segmentation，MS)方法仍然是目前最有效的方法[20]。針對(duì)SAR圖像存在的相干斑噪聲問(wèn)題以及SAR數(shù)據(jù)固有的數(shù)據(jù)特征，文獻(xiàn)[21]基于概率密度函數(shù)的SLIC(PDF-based SLIC)方法，實(shí)現(xiàn)了L波段的SAR圖像超像素分割。文獻(xiàn)[22—23]通過(guò)改進(jìn)SLIC超像素算法中的距離度量方法實(shí)現(xiàn)了PolSAR圖像中目標(biāo)的分類(lèi)。然而SLIC算法應(yīng)用于SAR圖像分割中依舊有些問(wèn)題沒(méi)有很好地解決。首先，大部分研究都是集中在改進(jìn)建立有效的距離度量以區(qū)分地物類(lèi)別，并未考慮到初始化步驟對(duì)于SAR圖像分割的重要性；其次，迭代過(guò)程中復(fù)雜的距離度量使得算法在大尺寸SAR圖像中效率低，不利于實(shí)際工程應(yīng)用。因此，需要考慮SAR圖像中疊掩和陰影檢測(cè)的實(shí)際情況，設(shè)計(jì)一種快速且準(zhǔn)確的SLIC處理流程，以及將SLIC嵌入到后續(xù)處理中實(shí)現(xiàn)疊掩和陰影的快速檢測(cè)。

本文基于改進(jìn)后的SLIC算法，結(jié)合相干系數(shù)、閾值幅度分割方法和局部頻率估計(jì)方法，提出一種基于多信息融合的超像素檢測(cè)算法，識(shí)別SAR圖像中的疊掩和陰影區(qū)域，來(lái)保證InSAR后續(xù)處理的可靠性和有效性。由于缺少真實(shí)場(chǎng)景的數(shù)字高程模型，導(dǎo)致無(wú)法驗(yàn)證本文算法在實(shí)際應(yīng)用中的準(zhǔn)確性，因此本文提出一種基于電磁散射模型的快速回波仿真法，將先進(jìn)積分方程近似模型(advanced integral equation model，AIEM)[24-25]和頻域回波仿真[26]相結(jié)合，通過(guò)利用衛(wèi)星工具包(satellite tool kit，STK)軟件設(shè)計(jì)仿真所需的軌道參數(shù)，實(shí)現(xiàn)星載全極化雙站SAR的原始回波數(shù)據(jù)仿真，回波能夠反映多站、多角度、全極化、復(fù)雜地物的回波特性以及反映SAR圖像的疊掩、陰影、相干斑噪聲等物理特征，可用于檢驗(yàn)本文所提出的基于多信息融合的超像素檢測(cè)算法有效性。結(jié)合天繪二號(hào)衛(wèi)星[27-28]的2019年9月獲取的河北赤城精度檢測(cè)場(chǎng)區(qū)域和2020年6月獲取的哥倫比亞區(qū)域數(shù)據(jù)，證明本文算法的有效性。

1 SAR圖像疊掩和陰影檢測(cè)理論研究

1.1 疊掩和陰影的幾何特性

工程上利用兩個(gè)在俯仰向觀測(cè)角度略有差別的接收天線，對(duì)同一目標(biāo)進(jìn)行各自成像處理，通過(guò)干涉的方法獲取目標(biāo)高程信息，稱(chēng)其為干涉合成孔徑雷達(dá)(interferometric synthetic aperture radar，InSAR)技術(shù)。如圖1所示，S1和S2代表裝載雷達(dá)天線的兩部衛(wèi)星，當(dāng)雷達(dá)觀測(cè)具有陡峭地形的山區(qū)或城市時(shí)，具有相同斜距的目標(biāo)點(diǎn)P1和P2在錄取回波時(shí)會(huì)疊加到同一個(gè)像素內(nèi)，將這一現(xiàn)象稱(chēng)為疊掩，將電磁波無(wú)法照射的區(qū)域成為陰影區(qū)域。由圖1中的幾何關(guān)系可以得到局部入射角θ的具體表達(dá)式為

θ=β-α

(1)

式中，β為雷達(dá)下視角；α為目標(biāo)點(diǎn)的坡度角。當(dāng)θ=0°時(shí)，斜坡的坡面與雷達(dá)視線垂直，在該目標(biāo)點(diǎn)兩側(cè)的點(diǎn)將具有相同的斜距，使得不同層次高程的目標(biāo)被疊加到同一距離單元內(nèi)形成疊掩；當(dāng)θ<0°時(shí)，此時(shí)α>β，底部的目標(biāo)點(diǎn)到雷達(dá)的斜距大于頂部的目標(biāo)點(diǎn)到雷達(dá)的斜距，使得SAR圖像上出現(xiàn)頂?shù)椎怪玫默F(xiàn)象，進(jìn)而形成疊掩現(xiàn)象。陰影區(qū)域一部分是由于電磁波被其他目標(biāo)遮擋而形成，另一部分是由于θ>π/2，即α<β-π/2時(shí)電磁波鏡面反射導(dǎo)致。

圖1 幾何關(guān)系Fig.1 Geometric relationship

1.2 基于信號(hào)特征的傳統(tǒng)逐像素檢測(cè)法

1.2.1 基于幅度閾值分割方法

1.2.2 基于局部頻率估計(jì)方法

文獻(xiàn)[8]指出InSAR系統(tǒng)形成的圖像之間存在頻率偏移Δf，并給出了具體表達(dá)式

(2)

式中，f0為載波的頻率(單位Hz)；R為雷達(dá)到目標(biāo)點(diǎn)的斜距(單位m)；k為常數(shù)，對(duì)于雙星收發(fā)分置雷達(dá)系統(tǒng)取2；B為基線長(zhǎng)度(單位m)；γ為基線傾角(單位rad)。由式(2)可以看出，當(dāng)β-α<0或β-α>π/2時(shí)Δf取值為負(fù)，此時(shí)對(duì)應(yīng)疊掩和陰影區(qū)域。

本文采用文獻(xiàn)[8]中頻譜分析方法估計(jì)局部頻率，在指定窗口內(nèi)做二維FFT求干涉相位圖頻譜，頻譜峰值出現(xiàn)在二維局部頻率對(duì)應(yīng)的位置，通過(guò)判斷頻率符號(hào)確定可能的疊掩和陰影區(qū)域。

1.2.3 基于CFAR檢測(cè)方法

均值類(lèi)CFAR算法是目前主流的CFAR檢測(cè)算法，其實(shí)現(xiàn)思路是以被檢測(cè)單元為中心的周?chē)鷨卧鳛閰⒖即埃瑢⒋翱趦?nèi)的所有數(shù)據(jù)平均作為背景功率估值，將估計(jì)的背景功率與系數(shù)(門(mén)限因子)相乘得到?jīng)Q策閾值T，根據(jù)不同的背景功率和門(mén)限因子。文獻(xiàn)[10]采用SO-CFAR算法進(jìn)行檢測(cè)，將參考窗分為對(duì)稱(chēng)的兩個(gè)窗口，以較小的背景雜波數(shù)據(jù)來(lái)估計(jì)背景功率，其背景功率為

(3)

(4)

式中，εso為門(mén)限因子。SO-CFAR算法適用于多個(gè)相鄰的目標(biāo)檢測(cè)，目標(biāo)周?chē)摼实停请s波邊緣區(qū)域虛警概率提升。

2 SAR圖像疊掩和陰影檢測(cè)算法實(shí)現(xiàn)

2.1 SLIC超像素分割算法的改進(jìn)

SAR圖像反映的是場(chǎng)景反射電磁波的能力，因此SAR圖像是強(qiáng)度圖像，傳統(tǒng)的SLIC算法適用于RGB圖像的分割，鑒于SAR圖像像素?cái)?shù)量通常達(dá)到106～108數(shù)量級(jí)，需要在保證分割精度的同時(shí)盡可能減少傳統(tǒng)SLIC算法運(yùn)算量以達(dá)到工程上快速分割的目的，本節(jié)對(duì)SLIC算法進(jìn)行改進(jìn)，并給出了改進(jìn)后的SLIC具體實(shí)現(xiàn)流程。預(yù)處理時(shí)需要將SAR圖像轉(zhuǎn)換成取值范圍為[0,1]的灰度圖像，然后放大100倍作為圖像的亮度信息，得到亮度圖L。

2.1.1 初始化聚類(lèi)中心

圖2 初始化種子點(diǎn)Fig.2 Initializing the seed point

為避免將超像素放置在圖像邊緣或噪聲點(diǎn)上影響聚類(lèi)性能，本文將聚類(lèi)中心移動(dòng)到與八連通域中的最低梯度位置相對(duì)應(yīng)的種子位置。梯度度量方法為

(5)

2.1.2 像素點(diǎn)的聚類(lèi)

SLIC算法通過(guò)檢測(cè)每個(gè)像素點(diǎn)周?chē)?S×2S范圍內(nèi)聚類(lèi)中心之間的相似度，將相似度最高的聚類(lèi)中心的標(biāo)簽賦予該像素。SLIC算法定義了最大的空間距離Ns=S和最大的顏色差異Nc來(lái)對(duì)空間距離和顏色距離進(jìn)行差異整合，由于Nc的值不容易直接定義，所以SLIC算法將Nc取值為常量m∈[1,40]，最終給出了像素間的相似性公式為

d(i,Ck)=

(6)

d(i,Ck)描述了像素i和聚類(lèi)中心Ck之間的相似度?？梢詫⑹?6)表示為向量的形式，特征向量f的加權(quán)L2范數(shù)形式為

(7)

式中，對(duì)角矩陣W=diag(1,λ,λ)，λ=m2/S2；f=[l,x,y]T。則d(i,Ck)可寫(xiě)成

(8)

(9)

進(jìn)一步通過(guò)近似的方法將d(i,Ck)近似到ds(i,Ck)減少了計(jì)算量，將聚類(lèi)問(wèn)題轉(zhuǎn)化為argminds(i,Ck)。

2.1.3 迭代更新

初始化標(biāo)簽圖像l={l(i)=-1|i∈[1,N]}，距離圖像D={d(i)=∞|i∈[1,N]}，循環(huán)k∈[1,K]計(jì)算Ck鄰域2S×2S范圍內(nèi)ds(i,Ck)，若ds(i,Ck)

2.1.4 后處理增強(qiáng)連通性

經(jīng)過(guò)上述迭代優(yōu)化可能出現(xiàn)多連通情況、超像素尺寸過(guò)小以及出現(xiàn)孤立超像素等，這些情況可以通過(guò)增強(qiáng)連通性解決。本文具體處理步驟如下：

(1) 為確保不相鄰區(qū)域的標(biāo)記是不同的，將具有不同區(qū)域但相同標(biāo)簽的重新分配標(biāo)簽。

(2) 對(duì)標(biāo)簽圖像進(jìn)行形態(tài)學(xué)開(kāi)操作處理，去除孤立的小點(diǎn)和毛刺。

(3) 上述操作完成后，某些區(qū)域可能會(huì)被關(guān)聯(lián)，分成多個(gè)區(qū)域或被其他區(qū)域吸收，為確保不相鄰區(qū)域的標(biāo)記是不同的，需要對(duì)區(qū)域重新編號(hào)，以便他們從1開(kāi)始順序增加。

上述步驟可能導(dǎo)致超像素個(gè)數(shù)增多或減少，最終輸出標(biāo)簽圖像l={l(i)|i∈[1,N]}，其中l(wèi)(i)表示SAR圖像上每個(gè)像素對(duì)應(yīng)的超像素編號(hào)，對(duì)應(yīng)的超像素特征結(jié)構(gòu)可表示為

(10)

2.2 基于多信息融合的超像素檢測(cè)算法

基于改進(jìn)后的SLIC超像素分割算法，結(jié)合InSAR相干系數(shù)、傳統(tǒng)閾值分割和局部頻率估計(jì)方法，本文提出了一種疊掩陰影判別算法。具體流程如圖3所示。

圖3 算法流程Fig.3 Algorithm flowchart

對(duì)于給定輸入為R×A=N的SAR圖像。首先，采用均值濾波的方式去除SAR圖像中的斑點(diǎn)噪聲；然后，對(duì)SAR圖像中的輪廓特征進(jìn)行精細(xì)化分割，獲得標(biāo)簽圖像l={l(i)|i∈[1,N]}和超像素特征結(jié)構(gòu)SPCk,k∈[1,K′]參數(shù)。初始化疊掩、陰影和疑似疊掩陰影的標(biāo)簽圖像Il、Is和Ils。計(jì)算SAR圖像平均強(qiáng)度Lavr，遍歷k∈[1,K′]，對(duì)每個(gè)超像素做局部頻率估計(jì)得到疊掩和部分陰影區(qū)的區(qū)域集合Ils，通過(guò)強(qiáng)度和相干系數(shù)閾值分割法在區(qū)域Ils中識(shí)別疊掩區(qū)Il，則可得到部分陰影區(qū)Is=Ils-Il。估計(jì)出部分陰影區(qū)的平均強(qiáng)度avr(L(Is))，將2×avr(L(Is))和相干系數(shù)小于0.6作為最終判定陰影區(qū)的閾值進(jìn)行陰影區(qū)域的檢測(cè)。陰影區(qū)域通常沒(méi)有回波能量，因此相干性較差，根據(jù)處理過(guò)程中的經(jīng)驗(yàn)，一般認(rèn)為相干系數(shù)小于0.6的區(qū)域包含陰影區(qū)域。在相干系數(shù)差的區(qū)域中通過(guò)設(shè)置閾值提取陰影區(qū)域。閾值分割得到的陰影區(qū)域中同樣包含回波能量弱的區(qū)域，這是由于SAR圖像每個(gè)像素值反映的是電磁波與目標(biāo)點(diǎn)作用后被接收機(jī)接收的回波能量值，鑒于散射系數(shù)受到入射角影響較大，對(duì)于入射角大于70°的面目標(biāo)而言，回波能量相對(duì)較弱，使得在數(shù)據(jù)分布特性方面與陰影區(qū)域相似，造成這部分區(qū)域相干系數(shù)較差。對(duì)于高程反演而言，相干系數(shù)較差的區(qū)域也會(huì)影響相位生成精度，因此回波能量弱的區(qū)域也應(yīng)該被有效檢測(cè)。

2.3 基于回波模擬技術(shù)的算法有效性檢驗(yàn)

為驗(yàn)證本文算法的正確性，需要明確疊掩和陰影的實(shí)際位置，由1.1節(jié)的分析可知，對(duì)于給定的觀測(cè)場(chǎng)景數(shù)字高程模型(DEM)，通過(guò)坡度和下視角的關(guān)系可判斷疊掩區(qū)域，通過(guò)射線追蹤的方法可判斷陰影區(qū)域。因此通過(guò)準(zhǔn)確地模擬回波錄取過(guò)程以獲取真實(shí)場(chǎng)景的SAR數(shù)據(jù)，可以檢驗(yàn)本文算法的有效性。SAR回波過(guò)程可被描述為發(fā)射信號(hào)與系統(tǒng)沖激響應(yīng)的卷積[29]，將其擴(kuò)展到雙基，具體可表達(dá)為

SRx(τ,f,t)=STx(τ,f)?τhc(τ,t)+n0(τ,f,t)

(11)

式中，τ為快時(shí)間參數(shù)(單位s)；t為慢時(shí)間參數(shù)(單位s)；f為頻率參數(shù)(單位Hz)；接收信號(hào)SRx(τ,f,t)被描述為式(11)中STx(τ,f)與通道脈沖響應(yīng)函數(shù)hc(τ,t)的卷積以及0均值高斯白噪聲n0(τ,f,t)的形式。在快時(shí)間的相干窗口內(nèi)假設(shè)發(fā)射機(jī)位置PosTx(t)和接收機(jī)位置PosRx(t)相對(duì)時(shí)間是不變的。此時(shí)目標(biāo)點(diǎn)的通道脈沖響應(yīng)可被寫(xiě)為

(12)

式中，δ(τ-R(t,r)/c)為雙基距離延時(shí)的δ函數(shù)；σ(t,r)為在采樣時(shí)刻t對(duì)應(yīng)的WGS-84坐標(biāo)系位置r=(x，y，z)處目標(biāo)點(diǎn)的反射率函數(shù)，反射率函數(shù)是由觀測(cè)幾何角、目標(biāo)點(diǎn)材質(zhì)等特性決定的；雙基距離函數(shù)R(t,r)為

R(t,r)=|PosTx(t)-r|+|PosRx(t)-r|

(13)

在每個(gè)方位向慢時(shí)間采樣時(shí)刻，4 dB波束寬度內(nèi)對(duì)應(yīng)的照射區(qū)域內(nèi)的脈沖響應(yīng)hc(τ,t)可被描述為

(14)

式中，W(t,r)是由雷達(dá)發(fā)射功率、收發(fā)天線增益以及雷達(dá)工作過(guò)程中由熱損失等引起的系統(tǒng)損耗等因素共同決定的參數(shù)；St為慢時(shí)刻t時(shí)4 dB波束照射到的場(chǎng)景目標(biāo)集合，去除其中的陰影區(qū)域。根據(jù)雙站SAR雷達(dá)方程[29]可得到W(t,r)的具體表達(dá)為

(15)

式中，Pt為發(fā)射信號(hào)的功率；Gt(t,r)為慢時(shí)間t接收位置r處目標(biāo)的發(fā)射天線方向圖增益；Gr(t,r)為接收天線方向圖增益。

SAR成像過(guò)程實(shí)質(zhì)上是重構(gòu)目標(biāo)點(diǎn)反射率函數(shù)的過(guò)程，因此SAR數(shù)據(jù)仿真的準(zhǔn)確性取決于反射率函數(shù)計(jì)算的準(zhǔn)確性。文獻(xiàn)[24]將AIEM模型與實(shí)測(cè)數(shù)據(jù)進(jìn)行對(duì)比，證明了AIEM模型在計(jì)算隨機(jī)粗糙面雙站散射系數(shù)的可靠性。AIEM模型考慮基爾霍夫項(xiàng)、交叉項(xiàng)和補(bǔ)償項(xiàng)，可以簡(jiǎn)化寫(xiě)成最終的簡(jiǎn)單形式，散射系數(shù)計(jì)算公式為

(16)

(17)

式中，K=sqrt((ksx-kx)+(ksy-ky)2)；ρ為相關(guān)長(zhǎng)度(單位m)。

在AIEM模型的推導(dǎo)坐標(biāo)系下φ=π，因此只要給出δ、ρ、θ、θs、φs這5個(gè)參數(shù)以及介電常數(shù)即可準(zhǔn)確計(jì)算目標(biāo)點(diǎn)的散射系數(shù)。通常對(duì)于目標(biāo)場(chǎng)景采用三角面元建模，每個(gè)場(chǎng)景面元都有一個(gè)確定的觀測(cè)幾何角θ、θs、φs，若用A表示三角面元的面積，則每個(gè)面元接收電磁波的有效面積為Acosθ。由于SAR圖像存在相干斑噪聲，因此模擬的SAR圖像也必須具有斑點(diǎn)噪聲特性。文獻(xiàn)[19]指出斑點(diǎn)噪聲被建模為實(shí)部和虛部服從0均值高斯分布的乘性噪聲。因此最終的反射率函數(shù)表達(dá)為

(18)

式中，N(0,1)表示服從均值為0方差為1的高斯分布特性的數(shù)。

根據(jù)上述建模分析，最終的雙站SAR回波方程可被描述為

δ(τ-R(t,r)/c)dr+n0(τ,f,t)

(19)

3 試驗(yàn)結(jié)果及分析

3.1 改進(jìn)SLIC算法性能分析

一方面，通過(guò)將疊掩和陰影區(qū)域檢測(cè)之后的結(jié)果作為掩膜圖屏蔽疊掩和陰影區(qū)域進(jìn)行干涉相位解纏，可以大大縮短解纏時(shí)間并提高解纏精度。對(duì)檢測(cè)出的疊掩和陰影區(qū)域進(jìn)行標(biāo)注，以便后續(xù)通過(guò)升降軌融合方法或已有的DEM對(duì)疊掩和陰影區(qū)域進(jìn)行修正，生成高精度的DEM產(chǎn)品。另一方面，TH-2致力于構(gòu)建全球DEM庫(kù)，因此必須盡可能提高疊掩和陰影的檢測(cè)效率，而超像素分割過(guò)程作為疊掩和陰影檢測(cè)的重要步驟，需要提高其運(yùn)算效率和分割性能。為了證明本文提出的改進(jìn)SLIC算法的有效性，對(duì)TH-2獲取的SAR數(shù)據(jù)集進(jìn)行了性能和運(yùn)算效率分析，并將本文算法與目前已被開(kāi)發(fā)能夠用于SRA圖像分割的方法(傳統(tǒng)SLIC[18]、MR方法[20]及PDF-based SLIC方法[21])進(jìn)行對(duì)比。其中，MR方法是集成在商業(yè)軟件eCognition v8.9中的方法，在該版本中MR方法并未使用多核并行計(jì)算技術(shù)進(jìn)行加速。為對(duì)比4種算法的運(yùn)算效率，其余方法均采用C++編程語(yǔ)言實(shí)現(xiàn)且不進(jìn)行并行加速。

為了評(píng)估不同分割方法的性能，本文采用3個(gè)指標(biāo)來(lái)評(píng)估分割結(jié)果，分別為邊界契合度(boundary fit，BF)、平均變異系數(shù)(mean coefficient of variation，MCV)和運(yùn)算效率。

BF參數(shù)用來(lái)衡量超像素分割邊界與真實(shí)邊界的一致性，定義為

(20)

式中，ES∩T用來(lái)描述通過(guò)算法分割的邊界和地物真實(shí)邊界的重合像素?cái)?shù)；ET表示地物真實(shí)邊界的像素?cái)?shù)。因此BF的值域?yàn)閇0,1]之間，BF的值越接近于1表示超像素塊與圖像邊緣的一致性更好。本文中BF參數(shù)選取局部分割區(qū)域進(jìn)行統(tǒng)計(jì)求解。

通常不同分割算法分割出的超像素塊之間差異較大，因此每個(gè)超像素塊的平均值略有差異，采用方差難評(píng)估多種算法之間的性能，CV參數(shù)可以消除由于不同超像素塊的平均水平不同對(duì)離散程度比較的影響。MCV參數(shù)定義為

(21)

如圖4所示，試驗(yàn)數(shù)據(jù)選自TH-2號(hào)2020年6月獲取的哥倫比亞區(qū)域?qū)崪y(cè)數(shù)據(jù)，大小為18 214×23 566像素，截取其中大小為3000×3000像素的SAR圖像作為輸入。在MR方法中，尺度參數(shù)設(shè)置為50，形狀參數(shù)設(shè)置為0.4，顏色參數(shù)設(shè)置為0.5。在SLIC方法中，初始超像素個(gè)數(shù)K取值均為10 000，顏色和空間重要性度量參數(shù)m取值15，其余所需參數(shù)均參考對(duì)應(yīng)文獻(xiàn)。由圖4結(jié)果可以看出，SAR圖像中通過(guò)分割得到的超像素區(qū)域中的顏色是均勻的，即屬于同一類(lèi)的像素在回波能量上具有相似的表現(xiàn)。從細(xì)節(jié)上可以看出不同算法在性能上有一定的差異，主要體現(xiàn)在邊界的定位以及分割效率上。對(duì)圖4中矩形的區(qū)域進(jìn)行評(píng)估，由圖4中黃色框和表1中BF和MCV參數(shù)結(jié)果可以看出MR方法和本文算法都能較好地分割出不同區(qū)域的邊界，而PDF-based SLIC方法和傳統(tǒng)SLIC方法分割性能相對(duì)較差。本次試驗(yàn)中采用MR方法進(jìn)行分割用時(shí)18.28 s，本文算法用時(shí)1.99 s，結(jié)果表明本文算法在更短的時(shí)間內(nèi)分割出的超像素更加均勻和緊湊，這為后續(xù)疊掩和陰影區(qū)域檢測(cè)提供了重要的保障。

表1 性能比較

圖4 TH-2數(shù)據(jù)超像素分割結(jié)果Fig.4 Superpixel segmentation results for the TH-2

在TH-2號(hào)2020年6月獲取的哥倫比亞區(qū)域數(shù)據(jù)集中截取[500,1000,2000,…,10 000]2數(shù)據(jù)量大小的圖像作為輸入，在每次試驗(yàn)中超像素個(gè)數(shù)設(shè)置為圖像行數(shù)的10倍，顏色和空間重要性度量參數(shù)m固定為15。具體結(jié)果如圖5所示。

圖5 運(yùn)算效率Fig.5 Operating efficiency

隨著輸入數(shù)據(jù)量的增大，MR方法需要根據(jù)相異性準(zhǔn)則對(duì)空間上非相鄰的區(qū)域進(jìn)行合并使得運(yùn)算量增大，PDF-based SLIC需要估計(jì)像素周?chē)母怕拭芏群瘮?shù)進(jìn)行聚類(lèi)使得計(jì)算量增加，由于本文算法在計(jì)算上的有效近似以及后處理的優(yōu)化使得計(jì)算效率明顯優(yōu)于其他方法。

3.2 仿真驗(yàn)證

本文利用STK軟件生成觀測(cè)河北赤城精度檢測(cè)場(chǎng)區(qū)域的軌道數(shù)據(jù)，通過(guò)構(gòu)造500 m基線得到發(fā)射機(jī)和接收機(jī)對(duì)應(yīng)的天線相位中心文件，其余具體仿真參數(shù)見(jiàn)表2。

表2 仿真參數(shù)

觀測(cè)區(qū)域的DEM見(jiàn)圖6(a)，本文采用西安電子科技大學(xué)InSAR課題組研發(fā)的RDSpace軟件對(duì)仿真回波數(shù)據(jù)進(jìn)行成像，本文給出發(fā)射星的成像結(jié)果見(jiàn)圖6(b)，對(duì)比圖6(d)成像后地理編碼的結(jié)果和圖6(a)可以證明本文回波仿真方法的有效性。

圖6 回波仿真SAR成像結(jié)果Fig.6 SAR results of echo simulation

對(duì)于3500×3588像素大小的SAR圖像，分別設(shè)置初始超像素個(gè)數(shù)K=[3,10,3]×10 000，顏色和空間重要性度量參數(shù)m=[15,15,5]。為方便呈現(xiàn)分割后的結(jié)果，本文選取圖7(a)紅色框出來(lái)的區(qū)域進(jìn)行分析，圖7(a)為紅色框中的SAR圖像，3種分割參數(shù)均實(shí)現(xiàn)了對(duì)SAR圖像的分割，并很好地保證了SAR圖像中目標(biāo)的輪廓信息。對(duì)比圖7(b)、(c)可以看出當(dāng)m取值相同時(shí)，超像素個(gè)數(shù)越多分割的結(jié)果越細(xì)，分割完成后的超像素個(gè)數(shù)接近初始的超像素個(gè)數(shù)，初始超像素越多適應(yīng)目標(biāo)邊界的能力越強(qiáng)。對(duì)比圖7(b)、(d)可以看出當(dāng)初始超像素取值相同時(shí)，m取值越小超像素形狀和尺寸越不規(guī)則，但是分割出的超像素更加均勻和緊湊，邊界粘附性更好。從圖7(d)結(jié)果可以看出，生成的超像素個(gè)數(shù)比輸入的超像素多，這是因?yàn)閰?shù)m取值越小顏色權(quán)值越大，此時(shí)算法更傾向于像素強(qiáng)度值上的鄰近性，使得在初次分割時(shí)大量相同標(biāo)簽的區(qū)域并未接壤，在后續(xù)處理中為確保不相鄰區(qū)域的標(biāo)記是不同的，對(duì)具有不同區(qū)域但相同標(biāo)簽的重新分配標(biāo)簽，進(jìn)而使得超像素個(gè)數(shù)增多。在疊掩和陰影檢測(cè)時(shí)，更關(guān)注不同類(lèi)別區(qū)域地物能量的變化，因此通常m取值較小，為保證結(jié)果的可靠性，一般m取值范圍為1～40。

超像素形狀和尺寸的規(guī)則性并不能提高疊掩和陰影檢測(cè)正確性，對(duì)于疊掩和陰影的檢測(cè)通常更關(guān)心目標(biāo)能量上的變化。通常希望分割的結(jié)果能夠更緊密地黏附到目標(biāo)邊界，因此設(shè)置初始超像素個(gè)數(shù)K=30 000，顏色和空間相對(duì)重要性度量指標(biāo)m=5進(jìn)行疊掩和陰影的檢測(cè)。并將本文算法與文獻(xiàn)[5—6，10]以及根據(jù)幾何關(guān)系判斷的真實(shí)疊掩和陰影區(qū)域進(jìn)行對(duì)比驗(yàn)證，詳細(xì)結(jié)果見(jiàn)圖8。圖8(a)、(e)是根據(jù)DEM和軌道信息幾何關(guān)系準(zhǔn)確判斷的疊掩和陰影位置，由于SAR圖像反映的是目標(biāo)點(diǎn)的反射率函數(shù)被脈沖響應(yīng)函數(shù)調(diào)制后的點(diǎn)散布函數(shù)形式，在該試驗(yàn)中SAR圖像一個(gè)像素上的點(diǎn)是5 m×4 m范圍內(nèi)的DEM點(diǎn)的疊加，因此首先需要將DEM反投影到SAR圖像上再進(jìn)行插值，才能確保高程點(diǎn)與SAR像素對(duì)應(yīng)?？梢钥闯鯟FAR檢測(cè)方法出現(xiàn)很多漏檢區(qū)域和無(wú)效區(qū)域，疊掩的檢測(cè)正確率不足50%。本文算法幾乎檢測(cè)出全部的疊掩區(qū)域，未檢測(cè)的部分大多是離散的孤立點(diǎn)區(qū)域，這些區(qū)域通常只占1～2個(gè)像素在形態(tài)學(xué)處理時(shí)被去除，對(duì)后續(xù)InSAR處理影響可忽略。對(duì)于陰影區(qū)域而言，陰影區(qū)域多為噪聲，回波能量較弱，因此陰影區(qū)域幾乎都可以根據(jù)強(qiáng)度值檢測(cè)出來(lái)，由圖8(h)中可以看出，檢測(cè)出的陰影區(qū)域比真實(shí)陰影區(qū)域較多，這是因?yàn)檫@些區(qū)域普遍集中在背坡面，具有較大的入射角，幾乎沒(méi)有回波信號(hào)，導(dǎo)致在圖像上呈現(xiàn)與陰影區(qū)域相似的特征。

圖7 SLIC超像素分割結(jié)果Fig.7 The results of SLIC superpixel segmentation

圖8 不同方法檢測(cè)結(jié)果Fig.8 The results of different methods

為驗(yàn)證檢測(cè)結(jié)果的有效性，定義疊掩和陰影檢測(cè)正確率為檢測(cè)出的真實(shí)疊掩和陰影數(shù)占準(zhǔn)確的疊掩陰影總數(shù)比例，定義疊掩檢錯(cuò)率為檢測(cè)錯(cuò)誤的總像素占檢測(cè)出的總像素的比例，為驗(yàn)證本文算法對(duì)DEM產(chǎn)品誤差精度的提升效果，通過(guò)計(jì)算未掩膜處理時(shí)DEM生成的高程誤差均值以及經(jīng)過(guò)疊掩和陰影區(qū)域掩膜處理后的高程誤差均值，并計(jì)算相位解纏速度提升指標(biāo)和DEM產(chǎn)品誤差精度提升指標(biāo)。具體指標(biāo)結(jié)果見(jiàn)表3。

表3 疊掩陰影檢測(cè)正確率

由表3可以看出，CFAR檢測(cè)方法對(duì)于復(fù)雜地形的SAR圖像表現(xiàn)效果較差，相干系數(shù)聯(lián)合幅度閾值檢測(cè)的方法檢測(cè)正確率相對(duì)于本文算法較低，具體是因?yàn)橄喔上禂?shù)是由周?chē)畔⒐烙?jì)得到的，幅度閾值是按照像素點(diǎn)去檢測(cè)的，所以實(shí)際的相干系數(shù)和像素點(diǎn)的幅度會(huì)有不一致性導(dǎo)致正確率降低，錯(cuò)檢率增加，這一現(xiàn)象在邊界處更為明顯。通過(guò)將本文算法檢測(cè)出的疊掩和陰影區(qū)域作為掩膜圖，在相位解纏過(guò)程中進(jìn)行掩膜使得相位解纏時(shí)間從16.64 s降低到13.13 s，相位解纏速度提升17.31%，生成的DEM產(chǎn)品高程誤差均值從2.32 m降低到1.91 m，DEM產(chǎn)品誤差精度提升17.58%，優(yōu)于其他方法。圖9給出了仿真數(shù)據(jù)InSAR處理后生成的DEM產(chǎn)品高程誤差值，誤差分析試驗(yàn)區(qū)域?yàn)閳D9中紅色線段位置，圖9(b)為未掩膜的DEM生成誤差，圖9(c)給出了掩膜后的DEM生成誤差，DEM生成誤差是通過(guò)InSAR處理得到的DEM產(chǎn)品和原始回波仿真輸入的觀測(cè)場(chǎng)景DEM做差得到。圖9中對(duì)比紅色方框中的結(jié)果可以明顯看出，由于解纏錯(cuò)誤導(dǎo)致誤差較大的區(qū)域，通過(guò)被有效地掩膜使得高程反演誤差降低，證明本文算法的有效性。

圖9 DEM產(chǎn)品高程誤差Fig.9 Elevation errors of DEM product

3.3 天繪二號(hào)數(shù)據(jù)試驗(yàn)

將本文所提算法應(yīng)用在天繪二號(hào)衛(wèi)星在2019年9月獲取的河北赤城精度檢測(cè)場(chǎng)SAR圖像中，對(duì)于6046×6971像素大小的SAR圖像，設(shè)置初始超像素個(gè)數(shù)K=30 000，顏色和空間重要性度量參數(shù)m=5。考慮到原始接收到的SAR數(shù)據(jù)中由于地物信息復(fù)雜存在強(qiáng)散射點(diǎn)數(shù)據(jù)，導(dǎo)致如果直接對(duì)SAR圖像量化后顯示，造成SAR圖像顯示較暗問(wèn)題，若不進(jìn)行處理會(huì)使得后續(xù)陰影檢測(cè)失效，因此需要去除SAR數(shù)據(jù)中后5%強(qiáng)度值數(shù)據(jù)再進(jìn)行后續(xù)試驗(yàn)。

如圖10所示，圖10(a)、(c)為本文算法檢測(cè)出的疊掩區(qū)域以及將疊掩區(qū)域邊界標(biāo)記到SAR圖像上結(jié)果，圖10(b)、(d)為本文算法檢測(cè)出的陰影區(qū)域以及將陰影區(qū)域邊界標(biāo)記到SAR圖像上結(jié)果。為證明本文算法的有效性，選取圖10中綠色框中的區(qū)域的SAR圖像以及對(duì)應(yīng)的相干系數(shù)圖，分別采用本文算法、CFAR檢測(cè)方法和相干系數(shù)聯(lián)合幅度閾值檢測(cè)方法進(jìn)行試驗(yàn)。試驗(yàn)結(jié)果如圖11所示，從中可以看出相干系數(shù)聯(lián)合幅度閾值檢測(cè)的方法檢測(cè)出的疊掩區(qū)域較真實(shí)疊掩區(qū)域更大，無(wú)法適應(yīng)疊掩的邊界，檢測(cè)出的陰影區(qū)域較小，且疊掩和陰影漏檢嚴(yán)重。相干系數(shù)聯(lián)合幅度閾值檢測(cè)方法受到SAR圖像固有的斑點(diǎn)噪聲影響檢測(cè)結(jié)果離散，檢測(cè)性能較差。實(shí)際SAR圖像相比于仿真的SAR圖像而言斑點(diǎn)噪聲更嚴(yán)重，使得相干系數(shù)聯(lián)合幅度閾值檢測(cè)方法失效。而本文算法檢測(cè)的疊掩和陰影區(qū)域邊界一致性較好，檢測(cè)性能更優(yōu)。

圖10 河北赤城精度檢測(cè)場(chǎng)SAR圖像試驗(yàn)結(jié)果Fig.10 The results of the SAR Image of Hebei Chicheng calibration field

圖11 不同方法檢測(cè)結(jié)果Fig.11 The results of different methods

為進(jìn)一步驗(yàn)證本文算法應(yīng)用于實(shí)測(cè)數(shù)據(jù)的有效性，將檢測(cè)出的疊掩和陰影區(qū)域在相位解纏過(guò)程中進(jìn)行掩膜。用于計(jì)算高程誤差的輔助DEM選自天繪中心測(cè)量得到的5 m精度的DEM數(shù)據(jù)，未掩膜時(shí)相位解纏用時(shí)23.77 s，DEM生成的高程誤差均值為3.79 m，掩膜后相位解纏用時(shí)16.93 s，DEM生成的高程誤差均值為2.74 m，相位解纏速度提升了28.78%，DEM產(chǎn)品精度提升了27.7%。選取圖12中紅色線段區(qū)域進(jìn)行分析，其對(duì)應(yīng)的經(jīng)緯度范圍從(115.623 683°E,41.210 212°N)到(115.666 799°E,41.179 376°N)，圖12(a)為未掩膜時(shí)DEM生成后的高程誤差值，圖12(b)為掩膜后DEM生成后的高程誤差值。為便于對(duì)比掩膜前后高程誤差變化結(jié)果，采用紅色框標(biāo)記同一區(qū)域，可以看出經(jīng)過(guò)掩膜處理后，圖12中紅色框區(qū)域內(nèi)的高程誤差明顯降低，證明了本文提出的基于超像素分割與多信息融合的疊掩和陰影區(qū)域檢測(cè)方法的有效性。

圖12 TH-2實(shí)測(cè)數(shù)據(jù)DEM產(chǎn)品高程誤差Fig.12 Elevation errors of DEM products based on TH-2 measured data

4 結(jié) 論

試驗(yàn)結(jié)果表明，本文提出的基于多信息融合的超像素檢測(cè)算法能夠快速準(zhǔn)確地提取SAR圖像中的疊掩和陰影區(qū)域，在疊掩區(qū)域識(shí)別上具有顯著的優(yōu)勢(shì)，在陰影區(qū)域識(shí)別上，由于低回波能量區(qū)域和實(shí)際地物遮擋形成的陰影區(qū)域在數(shù)據(jù)特征上相似，在沒(méi)有實(shí)際精確DEM輔助的情況下無(wú)法準(zhǔn)確區(qū)分低回波能量區(qū)域和實(shí)際地物遮擋形成的陰影區(qū)域，但是在大部分實(shí)際需求中這些區(qū)域都會(huì)對(duì)后續(xù)處理造成影響，因此陰影區(qū)域的準(zhǔn)確區(qū)分并不是必要的。本文所提算法檢測(cè)的結(jié)果可作為相位解纏環(huán)節(jié)中的掩膜產(chǎn)品，提升相位解纏的速度和DEM產(chǎn)品的精度。同時(shí)提取的超像素信息還可用于SAR目標(biāo)識(shí)別、升降軌融合等領(lǐng)域。本文提出的基于AIEM模型的頻域回波仿真算法，實(shí)現(xiàn)了星載全極化雙站SAR的原始回波數(shù)據(jù)仿真，模擬的回波能夠反映SAR圖像的疊掩、陰影、相干斑噪聲等物理特征，回波仿真結(jié)果可用于成像算法驗(yàn)證和系統(tǒng)設(shè)計(jì)參數(shù)驗(yàn)證。