馮瀚， 趙峰， 汪云甲， 閆世勇， 彭鍇， 王騰，張念斌， 徐東彪

1. 中國礦業(yè)大學(xué) 環(huán)境與測繪學(xué)院, 徐州 221116;

2. 國土環(huán)境與災(zāi)害監(jiān)測國家測繪地理信息局重點實驗室, 徐州 221116;

3. 黃河勘測規(guī)劃設(shè)計研究院有限公司, 鄭州 450003

1 引 言

時 序InSAR 技 術(shù)PSI （Persistent scatterer Interferometry）在地表變形監(jiān)測上具有監(jiān)測范圍廣、監(jiān)測精度高等優(yōu)勢，目前已被廣泛用于大區(qū)域地表形變探測（郭華東和張露，2019）。PSI 技術(shù)最早由Ferretti等（2001）提出，該技術(shù)作為差分干涉合成孔徑雷達DInSAR（Differential Interferometric Synthetic Aperture Radar）的延伸，為減小差分干涉中的失相干及大氣延遲影響，對多幅SAR 影像中的永久散射體PS（Persistent Scatterers）進行時序分析，其適用范圍和形變監(jiān)測精度均優(yōu)于傳統(tǒng)DInSAR 技術(shù)（朱建軍 等，2019）。之后，研究學(xué)者相繼推出SBAS（Berardino 等，2002）、干涉點目標分析（Werner 等，2003a 和2003b）、StaMPS（Hooper 等，2004）、時域相干點目標分析（Zhang等，2011）等時序InSAR 技術(shù)，并在各類變形監(jiān)測中進行了應(yīng)用。

但是，時序InSAR 技術(shù)在人工建筑物較少的非城市區(qū)域，可獲取的PS 監(jiān)測點密度低，難以準確反映地表變形，對相干性較低的分布式散射體區(qū)域的形變探測存在不足（Cao等，2015；Crosetto等，2016；林琿 等，2017）。因此，國內(nèi)外研究學(xué)者基于經(jīng)典PSI 算法提出了更高級的分布式DSI（Distributed Scatterer Interferometry） 技術(shù)，如聯(lián)合永久散射體（PS）與分布式散射體（DS）的SqueeSAR 技術(shù)（Ferretti 等，2011）、SqueeSAR 的各類變體技術(shù)（Parizzi 和Brcic，2011；Wang 等，2012）以及CAESAR 技術(shù)（Fornaro 等，2015），很好地彌補了傳統(tǒng)PSI 對DS 目標監(jiān)測的不足，推動了PSI的進一步發(fā)展。但是上述方法只使用了單極化數(shù)據(jù)，對多極化SAR 數(shù)據(jù)信息挖掘研究相對較少。然而，各類地表散射體對不同極化的響應(yīng)不同，利用不同極化通道的豐富極化信息可提高PSI技術(shù)地表形變監(jiān)測能力（Lee和Pottier，2009）。長期以來，長時序極化數(shù)據(jù)獲取難度大，基于多極化數(shù)據(jù)進行時序地表形變監(jiān)測較為困難。但是，隨著近年來具備多極化SAR 數(shù)據(jù)衛(wèi)星傳感器的相繼發(fā)射，如Radarsat-2、TerraSAR-X、ALOS-2、Sentinel-1 等，已經(jīng)具備了將PSI 技術(shù)擴展到多極化的現(xiàn)實基礎(chǔ)。2009 年，Pipia 等（2009）基于全極化地基SAR 數(shù)據(jù)，最早提出了極化時序InSAR的PolPSI（Polarimetric PSI）并首次成功將極化優(yōu)化算法應(yīng)用于PSI技術(shù)；此后許多學(xué)者利用星載極化SAR 數(shù)據(jù)開展了大量的研究工作，獲得了較好的研究成果（Navarro-Sanchez 等，2010；Navarro-Sanchez 和Lopez-Sanchez，2014；Iglesias 等，2014；Zhao和Mallorqui，2019a和2019b）。

Sentinel-1 數(shù)據(jù)作為一種可免費獲取的多極化SAR 數(shù)據(jù)，現(xiàn)有相關(guān)研究對其極化信息的利用不足，已有的研究案例只是簡單對比兩種極化數(shù)據(jù)的PSI 監(jiān)測效果（熊佳誠 等，2019），或者是基于單一相位質(zhì)量指標（例如DA或相干系數(shù)）對所有像元進行統(tǒng)一的干涉相位優(yōu)化（Shamshiri 等，2018；Azadnejad 等，2020），并未針對PS 與DS 目標不同特性進行自適應(yīng)的極化優(yōu)化，Sentinel-1 極化數(shù)據(jù)在地表形變中的應(yīng)用潛力還需進一步挖掘和研究。

有鑒于此，為充分利用Sentinel-1 數(shù)據(jù)的雙極化信息，本文借鑒DSI技術(shù)和相關(guān)研究思想，針對PS 與DS 目標各自的極化特性，提出了一種基于雙極化Sentinel-1 數(shù)據(jù)的PS 與DS 目標自適應(yīng)極化優(yōu)化的PolPSI 方法。以上海浦東國際機場為研究區(qū)，對比分析本文方法較傳統(tǒng)PSI技術(shù)優(yōu)勢，并對相關(guān)結(jié)果進行討論。

2 研究區(qū)概況及數(shù)據(jù)介紹

浦東國際機場為4F 級民用機場，是中國3 大門戶復(fù)合樞紐之一。該機場以南北向姿態(tài)坐落于長江三角洲沖積平原，地處上海市浦東新區(qū)，距上海市中心約30 km。據(jù)資料顯示，浦東機場旅客吞吐量居全國第二，僅次于首都機場，目前已擁有4 條運行跑道、兩個航站樓以及在建的第五跑道，其中第四跑道于2015 年建成并正式投入使用（劉冬明，2018）。該機場所在區(qū)域土層軟弱，機場由西向東地質(zhì)條件差異大：其西邊場地成陸時間早，東邊場地為圍涂駐塘后形成的灘地，地質(zhì)條件較為脆弱復(fù)雜。圖1為研究區(qū)域示意圖：浦東機場位于左圖中白框區(qū)域；右圖為機場的Sentinel-1強度影像時序平均圖，圖中黃色圓形與紅色方形標志分別標示典型PS、DS 目標位置，將用于后續(xù)結(jié)果分析。

圖1 研究區(qū)域及SAR強度影像圖Fig. 1 Study area and SAR intensity image

本文選取34 景雙極化Sentinel-1 數(shù)據(jù)進行PolPSI時序形變監(jiān)測。圖2為Sentinel-1影像時空基線分布圖，監(jiān)測時間從2016年5月至2017年12月，以2017年5月22日為主影像。

圖2 SAR影像時空基線分布圖Fig. 2 SAR images’ spatial and temporal baseline distribution

3 方法原理

本文借鑒DS-InSAR 思想，提出了PolPSI 技術(shù)。不同于常見DS-InSAR 技術(shù)對DS 目標基于其同質(zhì)像元平均來進行優(yōu)化，本研究利用最小均方差MMSE（Minimum Mean Square Error）濾波方法（Lee 等，1999）對DS目標進行處理。該方法是一種同質(zhì)像元加權(quán)平均技術(shù)，可更好保持地物細節(jié)信息。提出PolPSI技術(shù)整體方法流程如圖3所示。

首先通過統(tǒng)計檢驗識別同質(zhì)像元（Jiang 等，2017），并設(shè)置同質(zhì)像元閾值劃分DS 與PS 目標（Ferretti 等，2011）。在劃分目標過程中，對DS 目標再次使用DA準則（DA設(shè)為0.25）進行剔除，剔除的像元劃歸為PS（圖3）。

圖3 PolPSI方法流程圖Fig. 3 Flowchart of PolPSI

之后，針對不同目標，為保持PS 目標的空間分辨率，使用DA準則對PS目標進行極化優(yōu)化；對DS目標，為最大程度提升其相位質(zhì)量，先進行DS 目標MMSE濾波，之后使用相干性準則進行極化優(yōu)化。

最后，基于生成的優(yōu)化差分干涉圖，通過TPC方法選取高質(zhì)量像元（Zhao 和Mallorqui，2019c）進行PSI 處理（Hooper 等，2004），得到PolPSI 形變監(jiān)測結(jié)果。本文TPC閾值選為0.92（對應(yīng)相位標準差約為25°）。

3.1 DS目標MMSE自適應(yīng)濾波

DS 目標MMSE 自適應(yīng)濾波主要包括同質(zhì)像元識別和基于同質(zhì)像元識別結(jié)果的極化相干矩陣MMSE 濾波兩個步驟。本研究使用蔣彌（Jiang 等，2017）提出的HTCI 方法，以15×15 的窗口識別同質(zhì)像元，結(jié)果如圖4（b）所示：結(jié)合Google Earth影像分析（圖4（a）），識別結(jié)果在機場跑道、綠植草地等地區(qū)同質(zhì)像元個數(shù)較多，在建筑物等區(qū)域的同質(zhì)像元個數(shù)較少，符合實際情況，識別效果較好。

圖4 浦東機場影像Fig. 4 Pudong airport image

MMSE濾波能很好地保留影像的極化特征和邊緣清晰度（Lee 等，1999；Lee 和Pottier，2009）。該濾波方法基于識別的同質(zhì)像元進行自適應(yīng)加權(quán)平均濾波，較常見DSI方法中使用的同質(zhì)像元平均濾波具有更好的干涉圖濾波效果。

式（1）—（3）給出了VV 極化與VH 極化哨兵數(shù)據(jù)的Pauli 基向量ki與極化相干矩陣T4矩陣的構(gòu)建方法（Zhao和Mallorqui，2019b）：

式（1）—（3）中，S 代表SLC 影像數(shù)據(jù)，下標i、j表示兩個不同時刻，上標T表示矩陣的轉(zhuǎn)置運算，H表示矩陣的共軛轉(zhuǎn)置。極化相干矩陣T4的濾波權(quán)重b按式（5）計算，并通過式（6）完成MMSE 濾波（Lee等，1999）。

式（4）—（6）中，tr表示矩陣的跡，下標SHPS表示同質(zhì)像元，下標mean 代表窗口內(nèi)同質(zhì)像元的平均，coef 值取0.51（Lee 等，1999；Lee 和Pottier，2009）。影像濾波前后假彩色顯示如圖5 所示，圖中用于合成為假彩色影像的R、G、B 三個通道分別 為R =|SVV|2，G =|SVH|2，B = |SVV+SVH|2，其 中第二行對應(yīng)于第一行白框區(qū)域的放大顯示，黃色圓形標示PS1目標，紅色方形標示DS1目標，后續(xù)討論了兩種典型目標的極化優(yōu)化空間參數(shù)搜索結(jié)果。從圖5中可看到MMSE 濾波較好保持了強散射體的空間分辨率，如城市區(qū)域強散射體的旁瓣信號細節(jié)在濾波后被保留下來；濾波后同質(zhì)像元區(qū)域斑點噪聲明顯降低，同質(zhì)區(qū)更加平滑。

圖5 假彩色影像Fig. 5 The false color image

需要說明的是，利用MMSE 自適應(yīng)濾波處理之后的極化相干矩陣Tfilter，得到優(yōu)化后VV 極化差分干涉圖（干涉圖中PS 目標保持不變，DS 目標基于同質(zhì)像元進行自適應(yīng)濾波處理），進一步利用常見DSI 時序分析技術(shù)可得地表形變監(jiān)測結(jié)果。該DSI 形變監(jiān)測結(jié)果（MMSE 方法結(jié)果）將用于后續(xù)同本文提出PolPSI方法結(jié)果對比分析。

3.2 干涉相位的時序極化優(yōu)化

時序極化優(yōu)化通過構(gòu)建極化空間ω，遍歷搜索使像元質(zhì)量最優(yōu)的空間投影方向，隨后將影像重投至該方向，獲取最優(yōu)極化干涉相位。對雙極化哨兵影像，空間投影矢量ω構(gòu)建如下（Zhao 和Mallorqui，2019b）：

振幅離差為SAR 強度信息的統(tǒng)計，相干性可以用于評價干涉質(zhì)量。通常利用振幅離差來識別永久散射體（Ferretti 等，2001），用相干性評價干涉圖相位質(zhì)量（Lee和Pottier，2009）。因此，搜尋最優(yōu)極化也就是尋找使振幅離差最小或相干性最高的空間投影方向，并求解出α與ψ參數(shù)。極化數(shù)據(jù)的相干性（Shamshiri 等，2018）公式如式（8）和式（9）所示：

式（8）和（9）中k表示第k幅干涉圖，Nintf表示干涉圖總數(shù)量，γ表示干涉圖之間的相干性，Ωij與Tii為式（3）極化相干矩陣中的分量。

極化數(shù)據(jù)的振幅離差DA（Shamshiri 等，2018）公式如式（10）和式（11）：

式中，N表示第N幅SLC 影像，ki表示式（1）里的Pauli 基向量。需指出PS 目標優(yōu)化是基于SLC的，即最優(yōu)投影空間ω是使SLC像元質(zhì)量最好的矢量方向，后續(xù)需要基于優(yōu)化SLC 重新構(gòu)建差分干涉圖。

圖6 展示了PS1 目標（圖5 中黃色圓形標志）與典型DS1目標（圖5中紅色方形標志）的極化空間參數(shù)搜索結(jié)果。可發(fā)現(xiàn)PS與DS目標的極化參數(shù)在極化空間中呈平滑變化，并在某一位置取到最優(yōu)（圖中白色星號所示），說明優(yōu)化方法的有效性。

圖6 典型散射體空間參數(shù)搜索結(jié)果Fig. 6 Parameter searching of typical scatters

4 結(jié)果與討論

為分析本文提出的PolPSI 技術(shù)性能，將其與傳統(tǒng)單極化PSI 以及基于MMSE 濾波的DSI 方法（只利用VV 極化干涉圖，既MMSE 方法）結(jié)果進行比較分析。此外，從極化散射特性角度，對本文PolPSI優(yōu)化結(jié)果進行討論。

4.1 差分干涉圖結(jié)果分析

選取時間基線最長的第一幅差分干涉圖進行質(zhì)量分析。顯示原始VV 極化干涉圖、MMSE 濾波后VV 極化干涉圖與時序極化優(yōu)化干涉圖如圖7 所示。圖7中可以看出濾波、極化優(yōu)化后的干涉圖條紋比原始干涉圖更為清晰，時序極化優(yōu)化子圖7放大的機場航站樓（黃色虛線圈出）條紋最平滑，輪廓最清晰，說明時序極化優(yōu)化較好地保持了分辨率、相位優(yōu)化質(zhì)量最好；進一步選取相位梯度均值、相位標準偏差、殘差點數(shù)目均值作為相位優(yōu)化的評價指標，這3種評價指標越小，對應(yīng)干涉圖噪聲水平越低、優(yōu)化效果越好。如表1 所示，在3 項指標上，MMSE 濾波相比于原始干涉圖分別下降了10.76%、49.46%、22.81%，時序極化優(yōu)化（表中黑色字體加粗顯示）相比原始干涉圖分別下降了19.97%、57.47%、41.75%，優(yōu)化效果最好。

表 1 干涉圖質(zhì)量評價結(jié)果Table 1 Interferogram quality evaluation results

圖7 差分干涉圖Fig. 7 Differential Interferogram phase

此外，計算時序平均相位標準偏差如圖8 所示，圖中黃色圓形與紅色方形標志分別標示典型PS、DS 目標位置，后續(xù)統(tǒng)計散射體優(yōu)化參數(shù)表格（表2）。圖8 中可看出本文方法在機場跑道區(qū)域（紅色實線表示區(qū)域）以及建筑物區(qū)域取得了更好的優(yōu)化效果，其時序平均相位標準偏差最??；同時，統(tǒng)計不同方法時序平均相位標準偏差直方圖如圖9所示，可明顯看出兩種優(yōu)化方法都降低了時序干涉圖的相位標準偏差，說明兩種方法對時序干涉圖都起到了優(yōu)化效果。其中，本文的極化優(yōu)化方法效果最好。

圖8 時序平均相位標準偏差圖Fig. 8 Mean phase standard deviation of Interferograms

圖9 不同方法時序平均相位標準偏差直方圖Fig. 9 Standard deviation histograms of different methods

4.2 極化散射特性討論

為比較典型散射體最優(yōu)極化優(yōu)化效果，分別統(tǒng)計了3個典型PS目標（圖1中紅色圓形標示）與DS 目標（圖1 中紅色方形標示）的相位標準偏差與優(yōu)化空間參數(shù)如表2所示。從表中可以看出極化優(yōu)化最大程度降低了兩種典型散射體的相位標準差。根據(jù)Azadnejad 等（2020）的研究，極化空間參數(shù)α 值趨于0°對應(yīng)二面散射體、垂直極散射體等，VV 極化的貢獻度較大；α趨于90°對應(yīng)隨機體散射，VH極化能有較好響應(yīng)；而對于空間參數(shù)ψ，趨于±180°時表示目標對垂直極化響應(yīng)較好，越趨于0°則表示水平極化具有較好的響應(yīng)。表2中，PS散射體α 值趨于0°，ψ值趨于-180°，DS 則表現(xiàn)相反，兩種典型目標滿足此規(guī)律。

表 2 典型散射體優(yōu)化參數(shù)統(tǒng)計Table 2 Statistics of optimization parameters for typical scatterers

為進一步說明最優(yōu)極化與地表散射體之間的聯(lián) 系（Cloude 和Papathanassiou，1998；Navarro-Sanchez 等，2010），計算所有PS 與DS 目標極化優(yōu)化前后時序相位標準偏差均值，優(yōu)化前分別為0.51、1.73，優(yōu)化后為0.33、0.81；同時，統(tǒng)計研究區(qū)內(nèi)所有PS與DS目標最優(yōu)空間參數(shù)，繪制統(tǒng)計直方圖如圖10所示。

圖10 最優(yōu)參數(shù)統(tǒng)計直方圖Fig. 10 The optimal parameter histograms

從圖10（a）和10（c）所有PS 目標的統(tǒng)計圖中可以看出，其α 最優(yōu)值基本集中在0°附近，ψ最優(yōu)值集中分布在-180°附近，表明PS 目標對VV 極化的響應(yīng)更好，因此在使用哨兵數(shù)據(jù)進行城市區(qū)域PS 目標形變監(jiān)測時，利用VV 單通道可得到較好的監(jiān)測效果；對于DS 目標而言（圖10（b）和10（d）），α 最優(yōu)值從0°到90°均有覆蓋，ψ最優(yōu)值從-180°到180°近似均勻分布，表明VH 極化對于DS 目標優(yōu)化有一定的貢獻度，從側(cè)面說明VH 通道中蘊含一定的DS相位信息。

4.3 形變監(jiān)測結(jié)果

VV 極化PSI、MMSE 濾波DSI、時序PolPSI 這3 種方法時序平均形變結(jié)果如圖11 所示。傳統(tǒng)PSI算法獲得29961 個高質(zhì)量監(jiān)測點；基于MMSE 濾波的DSI 方法得到46511 個高質(zhì)量監(jiān)測點，相比傳統(tǒng)PSI 增加了55%；本文PolPSI 獲取到60829 個高質(zhì)量監(jiān)測點，相對傳統(tǒng)PSI 則增長近103.0%，相對MMSE的DSI方法增加30.8%。圖11 中東北角沿海岸邊的線性目標為修建的沿海堤壩，在針對該沿海堤壩的變形監(jiān)測上，傳統(tǒng)PSI算法僅能選取到零星的高質(zhì)量監(jiān)測點；基于MMSE 濾波的DSI 結(jié)果對堤壩上的沉降有所體現(xiàn)；本文PolPSI 方法結(jié)果，該線性堤壩上點位有明顯增加。

圖11 地表形變監(jiān)測結(jié)果（白框區(qū)域放大圖位于最右下角）Fig. 11 Ground deformation velocity maps obtained （The close-up area of white rectangle are locating in the bottom right corner）

除了沿海堤壩工程，選取第四跑道南北兩端的典型區(qū)域進行放大展示：傳統(tǒng)PSI方法在跑道同質(zhì)區(qū)基本沒有監(jiān)測出沉降信息；而基于MMSE 的DSI技術(shù)能監(jiān)測出一些形變點，但是密度不足，不能很好反映跑道上的形變情況，甚至在跑道南端區(qū)域沒能探測出形變信息（圖11（B2））；本文PolPSI在跑道上的點位數(shù)量明顯增多，沉降信息更為突顯，跑道南部區(qū)域（白框C1，C2 區(qū)域）探測出了小范圍的沉降（圖11（C2））；此外，PolPSI方法在其它區(qū)域的表現(xiàn)也更好。本文研究區(qū)監(jiān)測結(jié)果與前人利用高分辨率SAR 數(shù)據(jù)研究結(jié)果吻合（廖明生 等，2020）。

為進一步說明3種方法形變結(jié)果可靠性，統(tǒng)計了不同方法中共同監(jiān)測點的形變速率差異直方圖，如圖12 所示。從圖12 中可看出，本文PolPSI 方法與MMSE 濾波的DSI 方法和傳統(tǒng)PSI 的沉降速率差異都很小，基本在±4 mm 之內(nèi)，說明PSI、MMSE濾波的DSI 和本文PolPSI 這3 種方法監(jiān)測結(jié)果差異小，方法可靠。

圖12 不同監(jiān)測結(jié)果之間共同監(jiān)測點形變速率差異直方圖Fig. 12 Histogram of deformation rate difference of corresponding-points between different monitoring results

選取3 種方法穩(wěn)定和形變區(qū)域的共同監(jiān)測點P1 與P2（圖11（a）中紅色標示），繪制3 種方法得到時序形變，如圖13 所示。兩個特征點3 種方法的時序形變一致性很高，進一步說明了本文PolPSI監(jiān)測方法的可靠及有效性。

圖13 P1點（a）和P2點（b）時序形變圖Fig. 13 Time series deformation of points P1 （a） and P2 （b）

5 結(jié) 論

為充分挖掘Sentinel-1 雙極化數(shù)據(jù)信息用于地表形變監(jiān)測，本文提出一種基于雙極化Sentinel-1數(shù)據(jù)的PolPSI 地表形變監(jiān)測技術(shù)方法。以上海浦東國際機場為例，對比分析本文PolPSI 方法的監(jiān)測效果，得到以下結(jié)論：（1）本文PolPSI 技術(shù)通過對PS進行DA準則的極化優(yōu)化，對DS進行自適應(yīng)濾波與極化優(yōu)化，能很好地保持PS 目標分辨率，提升DS 目標相位質(zhì)量，在降低干涉圖相位噪聲的同時較好保留地物細節(jié)；（2）本文PolPSI 技術(shù)能顯著提升高質(zhì)量監(jiān)測點密度，相比于VV 極化傳統(tǒng)PSI 方法，監(jiān)測點增加103.0%，相比基于MMSE 自適應(yīng)濾波的DSI 方法，監(jiān)測點增加了30.8%；在一些機場跑道區(qū)域，由于哨兵數(shù)據(jù)分辨率較低，傳統(tǒng)方法未能得到形變監(jiān)測點，而本文PolPSI 技術(shù)能探測出這些區(qū)域形變；（3）對研究區(qū)SAR 像元最優(yōu)極化參數(shù)的統(tǒng)計分析顯示，PS目標對VV極化響應(yīng)較好，VH 極化對DS 目標干涉相位極化優(yōu)化有較大貢獻。