相位分塊與擬合法結(jié)合的InSAR相位解纏算法

馬靚婷，盧小平，余振寶

(河南理工大學(xué) 自然資源部礦山時空信息與生態(tài)修復(fù)重點實驗室，河南 焦作 454000)

0 引言

合成孔徑雷達(dá)干涉技術(shù)(interferometric synthetic aperture radar，InSAR)是在主動遙感基礎(chǔ)下發(fā)展的對地觀測技術(shù)，其利用同一目標(biāo)區(qū)域的SAR復(fù)圖像對共軛相乘得到干涉圖，根據(jù)干涉圖的相位值，計算出2次成像過程中信號產(chǎn)生的路程差，從而獲取監(jiān)測區(qū)域高精度的三維地形信息和微小的地形形變信息[1-2]。因全天時、全天候、方便快捷等優(yōu)點，InSAR技術(shù)被廣泛用來監(jiān)測地面沉降、地裂縫、火山活動、地震形變等[3]。

相位解纏是InSAR處理的重要環(huán)節(jié)，其結(jié)果的優(yōu)劣直接影響到地形測量的精度。受限于合成孔徑雷達(dá)的成像和處理方式，直接利用影像獲得的干涉圖一般含有較大的噪聲，局部相位殘差點的增多會形成不可靠數(shù)據(jù)斑塊，使該區(qū)域相位解纏出現(xiàn)漏解或錯解，導(dǎo)致InSAR圖像恢復(fù)失敗，從而影響到形變提取的精度。因此，提高相位解纏精度是InSAR處理中提高形變精度的重要環(huán)節(jié)[4]。

相位解纏是通過在解纏路徑上進行積分從而還原真實目標(biāo)信息。當(dāng)干擾因素少、相位質(zhì)量高時，能很好地還原真實相位信息；當(dāng)干擾因素多時，誤差會通過積分進行積累與傳播，得到的相位數(shù)據(jù)與真實數(shù)據(jù)會有較大的差異?，F(xiàn)有相位解纏方法主要分為路徑跟蹤解纏法[5-6]、最小范數(shù)解纏法[7]和網(wǎng)格規(guī)劃解纏法[8]3大類。路徑跟蹤法是通過設(shè)置合適的積分路徑，將誤差限制在一定區(qū)域內(nèi)，防止相位誤差全局傳遞，其涵蓋了經(jīng)典的Goldstein枝切法、質(zhì)量圖引導(dǎo)法、最小范數(shù)法等。枝切法是利用殘差點的連接得到枝切線，最后沿著枝切線進行積分得到解纏結(jié)果，但枝切法易出現(xiàn)孤島現(xiàn)象[5]。質(zhì)量圖引導(dǎo)法是在質(zhì)量圖的引導(dǎo)下確定積分路徑，這種算法對干涉圖質(zhì)量要求較高[9-10]。最小范數(shù)法是將相位解纏轉(zhuǎn)換成數(shù)學(xué)上的最小范數(shù)問題，其常用的是最小二乘法，但這種方法穿過殘差點會造成誤差的全局傳遞[11]。網(wǎng)格規(guī)劃法是將相位解纏問題轉(zhuǎn)化為求解費用流的網(wǎng)絡(luò)優(yōu)化問題，主要有最小費用流和統(tǒng)計費用流等，但這種方法噪聲會沿著積分路徑傳遞，使得解纏結(jié)果不理想[12-13]。

針對當(dāng)前相位解纏算法在相位圖像存在嚴(yán)重噪聲時解纏效果較差的問題，本文提出了基于相位分區(qū)與擬合法結(jié)合的InSAR干涉圖相位解纏算法。該算法：首先對纏繞相位進行分塊，將相位在相同區(qū)間的相鄰像素進行合并；塊內(nèi)像素屬于同一包裹次數(shù)，通過線性擬合求取合適的補償系數(shù)K進行塊間解纏；最后利用曲面擬合的方法進行殘差塊解纏，合并所有解纏塊得到最終的解纏結(jié)果。

1 基本理論

1.1 相位解纏繞

在進行InSAR數(shù)據(jù)處理時，將主輔影像共軛相乘并取相位信息即可得到復(fù)干涉條紋圖，但通過共軛相乘得到的相位差，與真實相位差相差2Kπ，即真實相位與纏繞相位的關(guān)系如式(1)所示。

Φ=φ+2Kπ

(1)

式中：Φ表示真實相位；φ表示纏繞相位；K表示補償系數(shù)，其值為整數(shù)。

相位解纏是從纏繞相位φ確定補償系數(shù)K值，進而估計真實相位Φ的過程。相位解纏必須兼顧一致性和精確性。一致性是指解纏后的相位數(shù)據(jù)矩陣中任意兩點間的相位差與其路徑無關(guān)；精確性是指解纏后的相位數(shù)據(jù)能真實地恢復(fù)原始相位信息[14]。

1.2 區(qū)域生長算法

區(qū)域生長算法的基本思想是將具有相似性質(zhì)的像素合并到一起。對每一個區(qū)域要先指定一個種子點作為生長的起點，然后將種子點周圍鄰域的像素點和種子點進行對比，將具有相似性質(zhì)的像素合并起來繼續(xù)向外生長，直到?jīng)]有滿足條件的像素被包括進來為止[15]。

2 本文算法

InSAR干涉圖雖然經(jīng)過了濾波處理，但干涉圖中仍存在噪聲，噪聲使得纏繞相位出現(xiàn)殘差點和低相干區(qū)域，使干涉圖解纏結(jié)果不全面或解纏誤差較大?；诜謮K與合并策略的相位解纏算法是一類高效的方法[16-18]，該類方法把整個纏繞相位圖分成若干塊并執(zhí)行塊間的相位解纏，合并所有塊可得到最終的解纏結(jié)果?，F(xiàn)有的解纏算法無法將誤差和孤島現(xiàn)象同時避免，因此本文提出針對InSAR干涉圖的分區(qū)與擬合法結(jié)合的相位解纏算法。

2.1 相位分區(qū)

本文將相位區(qū)域擴張法(prelude)中的相位分區(qū)用于InSAR干涉圖相位解纏中，當(dāng)分區(qū)間隔較大時相位塊內(nèi)會出現(xiàn)纏繞的現(xiàn)象，分區(qū)間隔較小時相位塊較為零碎，不利于殘差塊的識別，影響算法的效率和精度。因此，本文選擇以π/3為區(qū)間長度，對InSAR干涉圖纏繞的相位矩陣實行分區(qū)，將相位信息在同一個區(qū)域的相鄰像素點合并為像素塊。纏繞相位的相位區(qū)間都在(-π，π]間，根據(jù)選定的區(qū)間長度，將相位分為(-π，-2π/3]，(-2π/3，-π/3]，(-π/3，0]，(0，π/3]，(π/3，2π/3]，(2π/3，π]6個區(qū)間。遍歷整個相位矩陣，對其中相位屬于各個分區(qū)的像素點進行標(biāo)記，如圖1(a)所示，將屬于同一分區(qū)且相連的像素合并為塊，得到若干相位處于同一區(qū)間的像素塊，如圖1(b)所示。合并像素采用四鄰域法，即像素周圍的4個與其相鄰的同區(qū)像素可以進行合并，如圖1(c)所示，塊內(nèi)像素?fù)碛型话鼣?shù)K。設(shè)置閾值為50，像素塊大于等于50為正常塊，像素塊小于50為殘差塊。

圖1 相位解纏過程示意圖

2.2 擬合法解纏

在理想狀態(tài)下干涉圖不存在噪聲，相位梯度小于π，選擇一個點作為起始點，可直接進行積分解纏，但現(xiàn)實中噪聲、地形起伏和相干性較低等現(xiàn)象給相位解纏帶來了困難。利用線性擬合求取補償系數(shù)K值的方法對噪聲不敏感，在區(qū)域生長解纏繞中利用臨近的已解纏繞像素的相位信息求取K值可快速得到相鄰像素的真實相位信息，且能夠適應(yīng)較大的噪聲水平和相位的快速變化等情況[19]。

擬合法相位解纏是基于相鄰像素真實相位變換不大于π，使用相位分區(qū)方法對相位圖像進行分塊，塊內(nèi)相鄰像素的相位變化小于給定相位間隔，即相鄰像素之間的相位差小于π/3且塊內(nèi)像素?fù)碛型话鼣?shù)K，將塊間相位解纏問題轉(zhuǎn)化為線性擬合法求K值。選取塊間相鄰的2列像素，如圖1(d)所示，其2列像素的相位關(guān)系為式(2)所示。

Y=X+2Kπ

(2)

式中：X表示已解纏像素塊的相位值；Y表示未解纏像素塊的相位值；K為補償系數(shù)，其值為整數(shù)。

根據(jù)式(2)可知k為

(3)

式中：k為線性擬合系數(shù)，認(rèn)為其最接近的整數(shù)值為整數(shù)補償系數(shù)K；Xi表示已解纏像素的相位值；Yi表示未解纏像素的相位值。

纏繞相位塊之間的相位解纏采用區(qū)域生長方式進行。把塊作為一個解纏處理的基本單元，塊內(nèi)像素?fù)碛邢嗤陌螖?shù)K，選擇第一個塊為起始塊，距離起始塊中心位置最近的塊為生長塊。生長塊的解纏繞即是找到一個最佳的相位包裹次數(shù)K，選擇起始塊與生長塊相臨近的2組像素，進行線性擬合，找到最優(yōu)的補償系數(shù)進行相位解纏，并合并解纏后的2塊，之后進行下一個塊的解纏，直至所有正常的塊完成解纏與合并。

2.3 殘差塊曲面擬合

當(dāng)所有的正常塊完成相位解纏后，利用曲面擬合的方法對殘差塊進行解纏。根據(jù)解纏后相位是連續(xù)曲面的原理，通過殘差塊周圍像素的真實相位對殘差塊進行解纏，將殘差塊內(nèi)像素點10×10窗口內(nèi)所有經(jīng)過解纏的像素進行最小二乘曲面擬合，將其結(jié)果作為殘差像素的真實相位值，選擇第一個像素為起點利用區(qū)域生長算法將殘差塊內(nèi)的像素解纏完畢。

2.4 本文算法描述

本文所提相位解纏算法主要分3個部分，即進行相位分區(qū)，將像素點分為若干個塊；進行塊間的相位解纏；殘差像素塊解纏。其算法的具體實現(xiàn)有以下4步。

1)將輸入的InSAR干涉圖纏繞的相位矩陣實行相位分區(qū)，將相位信息在同一個區(qū)間且相鄰的像素點合并為像素塊，將像素小于50的像素塊標(biāo)記為殘差塊。

2)實現(xiàn)相鄰塊間的相位解纏。利用線性擬合法求取補償系數(shù)K進行相鄰塊間的相位解纏，根據(jù)區(qū)域生長算法將正常塊解纏完畢。

3)當(dāng)正常像素塊全部合并完成后，對標(biāo)記的殘差像素塊進行曲面擬合解纏。

4)合并所有解纏結(jié)果，完成整個區(qū)域的相位解纏工作。

3 實驗結(jié)果與分析

本文對模擬數(shù)據(jù)進行實驗并用實測數(shù)據(jù)驗證。針對InSAR相位的特性，從解纏的準(zhǔn)確度和算法運行的時間對解纏結(jié)果進行評價，除主觀視覺評價外，用均方根誤差和算法運行時間對實驗進行客觀評價。

3.1 模擬實驗

本實驗為MATLAB仿真的地形圖，用基于雷達(dá)傳感器參數(shù)和軌道數(shù)據(jù)的方法模擬InSAR干涉圖[20]。首先模擬無噪相位圖，如圖2(a)所示，加入噪聲后進行相位纏繞形成纏繞相位，如圖2(c)所示,并將其中50像素×50像素大小的范圍加重噪聲，以驗證算法在高強度噪聲下的效果。

圖2 模擬數(shù)據(jù)示意圖

將本文算法與Goldstein枝切法、質(zhì)量圖引導(dǎo)的路徑追蹤法、四向剪切最小二乘法解纏結(jié)果進行比較與分析。圖3(a)～圖3(l)為不同解纏方法結(jié)果。從目視效果看，Goldstein枝切法解纏結(jié)果在高強度噪聲區(qū)域出現(xiàn)了孤島現(xiàn)象，無論是曲面圖和俯視圖都存在明顯的錯誤。質(zhì)量圖引導(dǎo)的路徑追蹤法解纏結(jié)果在噪聲較弱、相位質(zhì)量好的區(qū)域解纏效果較好，在噪聲密集區(qū)域解纏結(jié)果不連續(xù)。四向剪切最小二乘法解纏結(jié)果相位曲面圖較為連續(xù)，但據(jù)其誤差直方圖可知解纏結(jié)果誤差較大。本文算法解纏結(jié)果與原曲面圖相位較為一致，在噪聲密集區(qū)域依然有好的效果且誤差值較小。

圖3 模擬數(shù)據(jù)實驗結(jié)果

表1是對模擬數(shù)據(jù)解纏結(jié)果的定量比較。從中可看出，四向剪切最小二乘法運算速度最快，但誤差全局傳遞造成其精度最低。質(zhì)量圖引導(dǎo)法均方根誤差小于四向剪切最小二乘法，但其運行時間大于四向剪切最小二乘法。本文算法的解纏結(jié)果均方根誤差最小，相較于四向剪切最小二乘法減少了44%，相較于質(zhì)量圖引導(dǎo)法減少了38%，且運行時間與質(zhì)量圖引導(dǎo)法相當(dāng)。

表1 不同方法的模擬干涉圖解纏結(jié)果定量比較

3.2 實測數(shù)據(jù)驗證

為驗證本文算法的有效性，利用山西平朔地區(qū)2011年12月17日和2012年2月27日的RadarSat-2實測數(shù)據(jù)進行驗證。對其數(shù)據(jù)進行配準(zhǔn)、干涉等處理得到其真實相位干涉數(shù)據(jù)，截取其中500像素×500像素大小進行實驗處理，如圖4(a)所示。

圖4(b)～圖4(h)顯示了不同方法的解纏結(jié)果。從目視效果看，Goldstein枝切法出現(xiàn)了孤島現(xiàn)象，存在明顯的錯誤。質(zhì)量圖引導(dǎo)法在低質(zhì)量區(qū)域解纏效果不佳，其解纏結(jié)果俯視圖上出現(xiàn)多個解纏效果不連續(xù)的區(qū)域。四向剪切最小二乘法解纏結(jié)果較為連續(xù)，但誤差較大。本文算法結(jié)果較為光滑避免了很多尖峰毛刺現(xiàn)象。

圖4 實測數(shù)據(jù)實驗結(jié)果

表2是對實測數(shù)據(jù)解纏結(jié)果的定量比較?？梢钥闯觯珿oldstein枝切法運算時間最長且解纏失敗，四向剪切最小二乘法運算時間最短但其均方根誤差最大。質(zhì)量圖引導(dǎo)法運算時間較長，但其精度高于四向剪切最小二乘法。本文算法解纏結(jié)果精度最高且運算時間適中。

表2 不同方法的實測干涉圖解纏結(jié)果定量比較

4 結(jié)束語

本文提出了一種相位分區(qū)與擬合法結(jié)合的相位解纏算法，該算法首先將獲得的相位圖像分為多個塊，塊內(nèi)相位都在給定的相位區(qū)間內(nèi)，把像素個數(shù)小于給定閾值的塊歸類為殘余像素；然后利用區(qū)域生長的擬合方法，依次進行塊與塊之間相位解纏繞和殘余像素相位解纏繞，合并后得到最終的解纏結(jié)果。通過對模擬數(shù)據(jù)進行實驗并用實測數(shù)據(jù)驗證，實驗結(jié)果表明，該算法無論是目視效果還是定量指標(biāo)分析均優(yōu)于其他算法，減小了誤差傳播的范圍，提高了相位解纏的精度，但該算法受相位連續(xù)性約束面對迭掩區(qū)域和復(fù)雜地形相位跳變問題時仍存在較大的改進空間。