劉萬(wàn)利 張秋昭 胡 江

（1.徐州工程學(xué)院數(shù)學(xué)與統(tǒng)計(jì)學(xué)院，江蘇 徐州 221018；2.中國(guó)礦業(yè)大學(xué)環(huán)境與測(cè)繪學(xué)院，江蘇 徐州 221116）

二維相位解纏是合成孔徑雷達(dá)干涉測(cè)量技術(shù)（Interferometric Synthetic Aperture Radar，InSAR）中的一個(gè)關(guān)鍵環(huán)節(jié)，解纏的成功與否直接影響到目標(biāo)高程信息的準(zhǔn)確性，因此一直是InSAR技術(shù)研究的熱點(diǎn)和難點(diǎn)［1-2］。近年來(lái)，許多學(xué)者對(duì)相位解纏進(jìn)行了深入研究，主要可分為兩大類：一類是基于積分的路徑跟蹤算法，另一類是最小范數(shù)算法［3］。質(zhì)量圖在這兩類算法中起著非常重要的作用，直接影響其解纏精度［4］。尤其對(duì)于路徑跟蹤算法，一個(gè)好的質(zhì)量指標(biāo)通常可以避免或減少積分路徑穿過(guò)殘差點(diǎn)，可以使得高質(zhì)量區(qū)域獲得比較高的精度［5］。

基于此，許多學(xué)者就如何選擇合適的質(zhì)量指標(biāo)展開(kāi)了相關(guān)研究［6-14］。GHIGLIA和PRITT在文獻(xiàn)［12］中比較了相干系數(shù)（Correlation Coefficient，CC）、偽相干系數(shù)（Pseudo Correlation Coefficient，PCC）、相位導(dǎo)數(shù)變化（Phase Derivative Variation，PDV）、最大相位梯度（Maximum Phase Gradient，MPG）4種質(zhì)量指標(biāo)，并指出相干系數(shù)是最有效的一個(gè)。OSMANOGLU等［13］描述了6種質(zhì)量指標(biāo)，試驗(yàn)結(jié)果顯示Fisher距離的性能總體上優(yōu)于相干系數(shù)、相位導(dǎo)數(shù)變化、結(jié)合枝切法的相位導(dǎo)數(shù)變化、二階導(dǎo)數(shù)可靠度及線掃描指標(biāo)等5種質(zhì)量指標(biāo)。LIU等［15］提出一種基于灰度共生矩陣（Gray Level Co-occurrence Matrix，GLCM）的熵差質(zhì)量指標(biāo)，該指標(biāo)對(duì)分級(jí)級(jí)數(shù)要求比較高，需要事先確定級(jí)數(shù)值。劉萬(wàn)利等［5］對(duì)幾種常見(jiàn)的質(zhì)量指標(biāo)做了比較深入的分析研究，并提出了一種改進(jìn)的Fisher距離指標(biāo)，試驗(yàn)結(jié)果表明，F(xiàn)isher和改進(jìn)的Fisher信息質(zhì)量指標(biāo)總體上比較穩(wěn)定。

綜上分析，基于相干系數(shù)和相位導(dǎo)數(shù)變化的質(zhì)量指標(biāo)簡(jiǎn)單易行，相干系數(shù)強(qiáng)調(diào)時(shí)間相干性，相位導(dǎo)數(shù)變化強(qiáng)調(diào)空間相似性［5］；并且上述質(zhì)量圖指標(biāo)大多單方面偏重時(shí)間相關(guān)、空間相關(guān)、噪聲等某一個(gè)指標(biāo)。InSAR技術(shù)在監(jiān)測(cè)山區(qū)礦區(qū)地表形變時(shí)，獲得的干涉圖通常具有噪聲大、條紋密且復(fù)雜等特點(diǎn)，單一質(zhì)量圖指標(biāo)用于山區(qū)礦區(qū)等復(fù)雜區(qū)域干涉圖相位解纏時(shí)，適應(yīng)性較差。為此，本研究提出一種融合相位導(dǎo)數(shù)變化和相關(guān)系數(shù)兩種質(zhì)量指標(biāo)的新的質(zhì)量圖指標(biāo)。該指標(biāo)包含了空間相關(guān)性、時(shí)間相關(guān)性等更多的質(zhì)量信息，因此更適用于地物覆蓋復(fù)雜地區(qū)、條紋比較密集干涉圖的相位解纏（尤其是山區(qū)煤礦區(qū)沉降監(jiān)測(cè)干涉圖）。采用模擬數(shù)據(jù)和實(shí)測(cè)數(shù)據(jù)試驗(yàn)驗(yàn)證了所提質(zhì)量指標(biāo)的可靠性和適用性。

1 相位解纏理論

基于濾波的相位解纏方法在展開(kāi)密條紋高噪聲干涉圖相位時(shí)具有獨(dú)特的優(yōu)勢(shì)［5，15］，本研究采用容積卡爾曼濾波相位解纏方法［16］進(jìn)行分析。

1.1 相位解纏模型

容積卡爾曼濾波相位解纏模型在繼承卡爾曼濾波相位解纏模型的基礎(chǔ)上，采用容積卡爾曼濾波器處理非線性方程，其狀態(tài)方程和觀測(cè)方程為

式中，y(k)和x(k)分別為像元k的觀測(cè)值和解纏相位；φ(k)為相位梯度估計(jì)值，通過(guò)最大似然估計(jì)器獲得；h[?]表示觀測(cè)函數(shù)；v(k)為k像元的觀測(cè)誤差向量；w(k)為梯度估計(jì)誤差，被認(rèn)為高斯白噪聲且滿足如下條件：

式中，E[?]表示期望；D[?]表示方差；為梯度估計(jì)誤差方差，其計(jì)算參考文獻(xiàn)［5］。v1(k)和v2(k)分別為復(fù)觀測(cè)值實(shí)部和虛部的觀測(cè)誤差，被認(rèn)為是高斯白噪聲且滿足如下條件：

式中，R(k)為觀測(cè)誤差協(xié)方差陣；，SkSNR為像元k的信噪比。

容積卡爾曼濾波相位解纏是沿著某一特定路徑進(jìn)行的，而路徑的選擇是基于某一質(zhì)量指標(biāo)從高質(zhì)量像元到低質(zhì)量像元。質(zhì)量指標(biāo)選擇得當(dāng)與否直接關(guān)系到解纏效果。

1.2 相位解纏步驟

本研究相位解纏實(shí)施流程如圖1所示。具體步驟為：①根據(jù)某一質(zhì)量指標(biāo)生成質(zhì)量圖；②選取質(zhì)量圖中的一個(gè)高質(zhì)量點(diǎn)A1作為解纏始點(diǎn)；③檢測(cè)像素點(diǎn)A1的4個(gè)鄰接像素點(diǎn)；④利用容積卡爾曼濾波解纏模型估計(jì)4個(gè)鄰接像素點(diǎn)中最高質(zhì)量的像素點(diǎn)A2；⑤選取像素點(diǎn)A1、A2鄰接的質(zhì)量最高的未解纏像素點(diǎn)A3并解纏。重復(fù)以上步驟，直到相位圖中的全部像素點(diǎn)都被解纏完畢。

2 基于相干系數(shù)與相位導(dǎo)數(shù)變化的新質(zhì)量指標(biāo)

相干系數(shù)和相位導(dǎo)數(shù)變化是兩種較常用的衡量干涉圖像質(zhì)量的質(zhì)量指標(biāo)，分別介紹如下。

2.1 相干系數(shù)

相干系數(shù)是最直接、最常用的衡量［12］指標(biāo)，計(jì)算公式為

式中，s1=c1+n1，s2=c2+n2，c1、c2為圖像的信號(hào)部分，n1、n2為圖像的噪聲部分；[?]*為共軛運(yùn)算。

對(duì)于實(shí)測(cè)干涉圖，CC是綜合多種失相干因素得出的結(jié)果，表征的是兩幅參與干涉的影像的相干性。通常情況下，反射能力強(qiáng)的地物對(duì)應(yīng)像元的CC值較高。本質(zhì)上是用相位的一階導(dǎo)數(shù)來(lái)描述相位質(zhì)量，對(duì)相位梯度比較敏感。因此，對(duì)于地形陡峭的區(qū)域，通過(guò)計(jì)算相位梯度，會(huì)將其標(biāo)識(shí)為低質(zhì)量數(shù)據(jù)。

2.2 相位導(dǎo)數(shù)變化

對(duì)于二維干涉圖數(shù)據(jù)，某像元的PDV值計(jì)算公式為［12］

由式（5）可以看出，PDV是以當(dāng)前像元為中心的一定窗口內(nèi)的像元在水平方向和豎直方向上梯度方差大小之和，反映的是目標(biāo)窗口內(nèi)梯度數(shù)值的離散程度，是梯度空間相似性好壞的評(píng)價(jià)指標(biāo)。本質(zhì)上是用相位的二階導(dǎo)數(shù)來(lái)描述相位質(zhì)量，對(duì)梯度變化比較敏感。因此在一些相干性較低的同質(zhì)地物覆蓋區(qū)，該指標(biāo)會(huì)表現(xiàn)出高質(zhì)量結(jié)果。

2.3 PDV和CC融合質(zhì)量指標(biāo)

由2.1及2.2節(jié)可以看出CC可以很好地反映像元的個(gè)體特征，但不能很好地反映空間相似性，對(duì)于一些存在孤立高相干點(diǎn)的干涉圖，容易造成解纏路徑穿過(guò)低質(zhì)量區(qū)，去解纏那些高質(zhì)量點(diǎn)然后再穿過(guò)低質(zhì)量區(qū)去尋找次高質(zhì)量點(diǎn)，從而造成誤差過(guò)早累積。而PDV會(huì)將一些空間上勻質(zhì)但相干性不高的區(qū)域誤判為高質(zhì)量區(qū)?；诖耍狙芯刻岢隽艘环N融合PDV和CC的質(zhì)量指標(biāo)。運(yùn)算步驟如下：

（1）數(shù)據(jù)標(biāo)準(zhǔn)化。由式（4）和式（5）獲得的CC和PDV值量綱和數(shù)量級(jí)不同，首先對(duì)其進(jìn)行標(biāo)準(zhǔn)化處理可得：

由式（7）可知，PDV_CC值越小，質(zhì)量越高，該指標(biāo)融合了像元相位的空間相似性和時(shí)間相干性兩個(gè)特征，是一種信息比較全面的質(zhì)量指標(biāo)。

3 試驗(yàn)與分析

我國(guó)多數(shù)煤礦區(qū)地表下沉速度快、植被覆蓋多，尤其是山區(qū)礦區(qū)地形起伏大，因此InSAR干涉圖通常呈現(xiàn)出條紋密集復(fù)雜、噪聲大、殘差點(diǎn)多等特點(diǎn)。如果質(zhì)量指標(biāo)選用不當(dāng)，會(huì)造成解纏路徑穿過(guò)大量殘差點(diǎn)，進(jìn)而導(dǎo)致沉降監(jiān)測(cè)精度不理想。本研究采用一組多山場(chǎng)景模擬數(shù)據(jù)（干涉圖具有條紋密集、噪聲大的特點(diǎn)）和礦區(qū)實(shí)測(cè)數(shù)據(jù)對(duì)所提出的質(zhì)量指標(biāo)進(jìn)行試驗(yàn)分析。為詳細(xì)比較試驗(yàn)結(jié)果，整個(gè)相位解纏過(guò)程只有質(zhì)量指標(biāo)不同，解纏算法采用容積卡爾曼濾波相位解纏方法。采用解纏誤差d和偏差χ2兩個(gè)指標(biāo)進(jìn)行定量評(píng)價(jià)，兩者計(jì)算公式為

式中，φu為解纏相位值；φr為相位先驗(yàn)值；σ為標(biāo)準(zhǔn)差；N為樣本數(shù)。

3.1 模擬數(shù)據(jù)試驗(yàn)

多山場(chǎng)景地形模擬數(shù)據(jù)如圖2（a）所示。為了詳細(xì)分析地形起伏對(duì)相位解纏的影響，模擬數(shù)據(jù)的垂直基線設(shè)為150 m，像元個(gè)數(shù)為250×250。數(shù)據(jù)具有地勢(shì)起伏大、陡峭度較高等特點(diǎn)。此場(chǎng)景的模擬成像加噪干涉圖如圖2（c）所示。本試驗(yàn)的含噪干涉圖體現(xiàn)的是幾何失相干程度，幾何失相干嚴(yán)重的區(qū)域噪聲強(qiáng)度大，由干涉圖可以看出多個(gè)區(qū)域條紋比較密集。

圖3給出了模擬數(shù)據(jù)的3種質(zhì)量圖。從大致的紋理來(lái)看，CC質(zhì)量圖與其余兩種質(zhì)量圖差別較大，高低質(zhì)量區(qū)區(qū)分比較明顯。這是因?yàn)槟M數(shù)據(jù)CC值的計(jì)算僅體現(xiàn)了幾何失相干，條紋密集的區(qū)域質(zhì)量低，稀疏的區(qū)域質(zhì)量高。PDV_CC與PDV的大致紋理比較接近，但顏色有所區(qū)別，前者顏色介于CC和PDV之間。這是由于PDV_CC是融合了CC和PDV兩種質(zhì)量指標(biāo)所致。另外，由圖3還可以看出在條紋稀疏區(qū)，PDV_CC和PDV也有少量低質(zhì)量散點(diǎn)，是因?yàn)楹敫缮鎴D的噪聲添加具有隨機(jī)特性，而這兩者的計(jì)算都是基于含噪干涉相位的。

圖4給出了3種質(zhì)量指標(biāo)引導(dǎo)解纏路徑的相位解纏結(jié)果圖。相對(duì)于真實(shí)相位圖（圖2（b）），可以看出CC引導(dǎo)的解纏失?。▓D4（a）），PDV（圖4（b））和PDV_CC（圖4（c））解纏結(jié)果總體保持了原樣，但在條紋非常密集的區(qū)域解纏結(jié)果表現(xiàn)出了不連續(xù)現(xiàn)象，即解纏相位誤差較大。這是因?yàn)檫@些區(qū)域噪聲強(qiáng)度大、圖像質(zhì)量低、條紋密集且變化較快，進(jìn)而導(dǎo)致條紋梯度的估計(jì)誤差較大。進(jìn)一步觀察可以發(fā)現(xiàn)，在條紋密集的區(qū)域，PDV_CC較PDV解纏結(jié)果更連續(xù)光滑（主要體現(xiàn)在橢圓和矩形區(qū)域）。說(shuō)明兼顧空間相似性和時(shí)間相干性的PDV_CC質(zhì)量指標(biāo)優(yōu)于僅考慮空間相似性的PDV質(zhì)量指標(biāo)。

為了進(jìn)一步討論所提質(zhì)量指標(biāo)的可靠性，圖5給出了3種指標(biāo)的解纏誤差圖及其對(duì)應(yīng)的誤差直方圖。從誤差圖可以明顯看出CC解纏誤差圖右側(cè)部分出現(xiàn)了大面積的大誤差區(qū)域，PDV和PDV_CC解纏誤差圖只是在條紋密集噪聲較大的區(qū)域存在較大誤差，這說(shuō)明CC指標(biāo)最不可靠。對(duì)比圖5（b）和圖5（c），可以看出PDV_CC方法解纏圖中橢圓和矩形區(qū)域的不連續(xù)區(qū)域明顯小于PDV方法，這說(shuō)明PDV_CC質(zhì)量指標(biāo)優(yōu)于PDV。從誤差直方圖也可以看出PDV_CC的誤差分布圖形中心峰（落在0附近部分）高度最高且寬度最窄，說(shuō)明該指標(biāo)引導(dǎo)的解纏結(jié)果誤差方差最小。因此PDV_CC的解纏結(jié)果優(yōu)于其他兩種指標(biāo)。

表1給出了模擬數(shù)據(jù)3種質(zhì)量指標(biāo)性能的定量評(píng)價(jià)指標(biāo)取值，分別為誤差均值、誤差方差以及偏差均值。從解纏誤差和偏差表達(dá)式可知，以上3個(gè)定量評(píng)價(jià)指標(biāo)的值越接近0說(shuō)明解纏結(jié)果越好。因此從表1可以看出CC指標(biāo)最不可靠，PDV_CC優(yōu)于CC，是一種比較魯棒可靠的質(zhì)量指標(biāo)。

分析偏差值隨解纏路徑的變化情況可以很好地解釋誤差從哪里產(chǎn)生、誤差隨路徑如何傳播等問(wèn)題。一個(gè)好的質(zhì)量指標(biāo)應(yīng)該是較大的偏差值出現(xiàn)得越晚越好，這樣可以避免誤差較早累積。模擬數(shù)據(jù)解纏結(jié)果的偏差值隨解纏路徑的走向曲線圖如圖6所示?？梢钥闯觯珻C的偏差值最大，PDV次之，PDV_CC最?。籆C方法的大誤差出現(xiàn)最早（在3 700個(gè)像元左右），PDV次之，PDV_CC方法的大誤差出現(xiàn)最晚（在56 000個(gè)像元左右）。這說(shuō)明PDV_CC方法出現(xiàn)誤差累積現(xiàn)象較晚，有更多的像元解纏結(jié)果較可靠。

3.2 礦區(qū)實(shí)測(cè)數(shù)據(jù)試驗(yàn)

研究區(qū)為神東礦區(qū)大柳塔煤礦52304工作面，像元個(gè)數(shù)為300×250，區(qū)域內(nèi)地形略有緩坡，有荒草沙石覆蓋。植被類型有低矮灌木、沙生植被、荒草等。干涉圖由2013年4月2日及2013年4月24日兩景TerraSAR-X 影像生成（圖7（a）），外部DEM采用課題組通過(guò)收集的航拍地形圖資料并結(jié)合GPS控制點(diǎn)生成的Relief-DEM。干涉圖的生成采用Gamma軟件。由于地表覆蓋導(dǎo)致的失相干較嚴(yán)重，故干涉圖采用的是2視處理結(jié)果。3種質(zhì)量指標(biāo)對(duì)應(yīng)的質(zhì)量圖見(jiàn)圖7（b）至圖7（d）。由于試驗(yàn)采用的是沉降監(jiān)測(cè)干涉圖，整個(gè)區(qū)域的真實(shí)沉降值難以獲得，因此解纏結(jié)果主要從連續(xù)性角度進(jìn)行分析。

由圖7（b）至圖7（d）可以看出，CC質(zhì)量圖的高質(zhì)量區(qū)比較小，而PDV和PDV_CC顯現(xiàn)出大面積高質(zhì)量區(qū)，這是因?yàn)楹髢蓚€(gè)指標(biāo)考慮了空間相似性，而研究區(qū)地面覆蓋空間相似性較高。PDV_CC與PDV質(zhì)量圖相似度較高，但比后者顏色稍淡，與公式所體現(xiàn)的特征相符。由圖7（e）至圖7（g）解纏相位圖可以看出，CC質(zhì)量指標(biāo)解纏結(jié)果在形變區(qū)右下角出現(xiàn)了較大的相位突變，相位連續(xù)性最差。在此區(qū)域，PDV和PDV_CC兩種方法解纏效果相當(dāng)，均優(yōu)于CC解纏方法。結(jié)合干涉圖和質(zhì)量圖可以看出，此區(qū)域噪聲大并且失相干比較嚴(yán)重，說(shuō)明在低質(zhì)量區(qū)PDV和PDV_CC指標(biāo)魯棒性較好。在形變區(qū)右上方，信噪比較高，CC和PDV_CC取得了比較連續(xù)的解纏結(jié)果，而PDV方法出現(xiàn)了相位跳變。這說(shuō)明PDV指標(biāo)在引導(dǎo)解纏路徑時(shí)出現(xiàn)了從低質(zhì)量到高質(zhì)量的情況，導(dǎo)致誤差提前大量累積從而造成高質(zhì)量區(qū)解纏結(jié)果不理想。從3種質(zhì)量指標(biāo)引導(dǎo)的解纏順序圖（圖7（h）至圖7（j），尤其是圖7（j）橢圓區(qū)域）可以看出，解纏結(jié)果的不同是由于質(zhì)量指標(biāo)引導(dǎo)的路徑不同導(dǎo)致的誤差累積不同所致。由此試驗(yàn)可以看出本研究提出的PDV_CC質(zhì)量指標(biāo)魯棒性較好。

綜合以上試驗(yàn)結(jié)果可以看出，無(wú)論是模擬數(shù)據(jù)還是實(shí)測(cè)數(shù)據(jù)，PDV_CC都取得了優(yōu)于PDV和CC的解纏結(jié)果。這說(shuō)明融合指標(biāo)可靠性好，是一種比較適合引導(dǎo)復(fù)雜地物覆蓋、沉降較快工作面干涉圖相位解纏的指標(biāo)。此外，對(duì)于模擬數(shù)據(jù)，CC指標(biāo)解纏結(jié)果極其不可靠，是因?yàn)槟M的地形在小范圍內(nèi)起伏變化較快，存在高質(zhì)量孤島；CC指標(biāo)引導(dǎo)的路徑出現(xiàn)了從低質(zhì)量到高質(zhì)量情形，從而造成了誤差累積傳播。對(duì)于實(shí)測(cè)數(shù)據(jù)，在沉陷區(qū)PDV_CC和PDV取得了相似的解纏結(jié)果，較符合實(shí)際沉降情況，但在非沉陷區(qū)PDV出現(xiàn)了解纏相位不連續(xù)現(xiàn)象，而CC和PDV_CC結(jié)果相似且較符合實(shí)際情況。兩種試驗(yàn)結(jié)果有所差別，主要是因?yàn)槟M數(shù)據(jù)的噪聲產(chǎn)生原因比較單一，主要由幾何失相干引起，而實(shí)測(cè)數(shù)據(jù)噪聲原因比較復(fù)雜。

4 結(jié)論

（1）提出了一種融合PDV和CC的新質(zhì)量指標(biāo)，該指標(biāo)不僅考慮了時(shí)間相干性還兼顧了空間相似性，更適用于噪聲大、條紋密集的復(fù)雜干涉圖相位解纏。

（2）高山地形模擬數(shù)據(jù)和礦區(qū)實(shí)測(cè)數(shù)據(jù)試驗(yàn)均表明，新的質(zhì)量指標(biāo)與單一質(zhì)量指標(biāo)相比，具有較好的魯棒性和實(shí)用性。

（3）相位梯度估計(jì)噪聲與質(zhì)量指標(biāo)的相互影響機(jī)制是下一步研究的方向。