時序InSAR數(shù)據(jù)在大范圍地面沉降監(jiān)測中的應用研究

胡 守 超

(山東電力工程咨詢院有限公司，山東 濟南 250013)

1 概述

合成孔徑雷達差分干涉測量(DInSAR)技術(shù)能夠快速、高精度地提取大范圍的地表形變信息，具有強大的技術(shù)優(yōu)勢和應用潛力。然而，常規(guī)DInSAR技術(shù)由于受時間/空間失相干和大氣效應等因素的影響，在實際應用中具有一定的局限性，尤其對于地面沉降這種緩慢累積的形變而言，形變信號往往被淹沒在噪聲之中，如何將其準確可靠地提取出來是目前研究的重點問題。

InSAR時間序列處理方法利用長時間序列SAR數(shù)據(jù)集，從中尋找具有穩(wěn)定相位特性的點，通過對這些相干點相位信息時間序列的分析和建模，將形變信號和噪聲分離，從而獲得可靠的形變信息時間序列。InSAR時間序列處理方法將InSAR面狀測量轉(zhuǎn)換為點測量，以犧牲空間分辨率為代價，換取高精度、可靠的結(jié)果。本文采用基于最小二乘的相干點分析方法并將其應用于大范圍地面沉降監(jiān)測，選擇13景ALOS PALSAR數(shù)據(jù)進行實驗，并利用實測的水準數(shù)據(jù)進行對比分析。

2 InSAR時間序列處理方法

2.1 多主影像干涉組合

為了使數(shù)據(jù)得到充分利用，采用多主影像的組合方式。在(t0,t1,…,tM-1)時刻獲得同一區(qū)域M幅SAR影像，以t0時刻的SAR影像為參考影像，其余影像與參考影像配準并重采樣。根據(jù)基線條件，選擇不同的主影像形成N個干涉組合M/2≤N≤M(M-1)/2，經(jīng)干涉處理、去除平地效應和去除地形相位，得到N幅差分干涉圖。

2.2 大氣延遲改正模型

在沒有嚴重地質(zhì)災害的情況下，城市地面沉降一般是一個比較緩慢且穩(wěn)定的過程，兩景SAR影像之間的地表形變只是一個微小量，而大氣延遲相位和軌道誤差等各種誤差導致的相位在干涉相位中卻占有較大的比重。要獲取高精度的地表形變信號，需要采取措施提高干涉相位的信噪比。考慮到大氣延遲相位在空間上具有強相關(guān)性，可以采用多項式擬合的方法對其進行擬合，同時為了盡可能地減少擬合模型對形變信號的影響，本文采用雙線性模型消除大氣延遲相位，獲得每幅干涉圖大氣延遲相位的初始估計，并分別將其從差分干涉相位中去除。

2.3 最小二乘方法求解線性形變

消除大氣延遲初始相位后，將差分干涉相位組成觀測方程：

y=Bx+V

(1)

其中，y=[δφ1，δφ2，…,δφN]T;B為系數(shù)陣；V為殘差矩陣。利用最小二乘，并設觀測值的權(quán)為P，可得參數(shù)x的最佳無偏估計：

(2)

3 實驗區(qū)域及數(shù)據(jù)

本文選擇上海市主城區(qū)為實驗區(qū)，上海地面沉降始現(xiàn)于1921年，具有顯著的階段性，經(jīng)歷了微量回彈、微量沉降、沉降加速和沉降減緩四個階段。目前上海市年均沉降量在10 mm以內(nèi)。

實驗數(shù)據(jù)采用日本先進陸地觀測衛(wèi)星(ALOS)搭載的相控陣L波段合成孔徑雷達(PALSAR)數(shù)據(jù)，收集了2007年2月～2009年10月上海地區(qū)的13景PALSAR影像。選擇中心城區(qū)為實驗區(qū)域，整個區(qū)域約180 km2。實驗區(qū)內(nèi)有9個水準點，水準測量結(jié)果將用于InSAR結(jié)果的校正和檢驗。

4 實驗結(jié)果分析

選擇垂直基線小于1 km的影像，形成27組干涉組合。差分干涉基于ROI_PAC軟件進行，采用二軌差分方式，利用SRTM的DEM數(shù)據(jù)去除地形相位。由于數(shù)據(jù)包含兩種模式，F(xiàn)BD和FBS模式，因此在進行干涉之前首先對FBD數(shù)據(jù)進行重采樣，使其與FBS數(shù)據(jù)具有相同的采樣率。經(jīng)干涉、去平、去除地形相位和濾波等步驟，并對濾波后的差分干涉圖進行地理編碼，得到實驗區(qū)域?qū)?7幅差分干涉圖，如圖1所示。

從圖1中可以看出，多數(shù)差分干涉圖的相位變化比較緩和，沒有條紋產(chǎn)生，此時可以認為大氣在空間上沒有劇烈變化，利用多項式擬合可以將其影響減弱至一定的程度。對于071125-080527,080110-080527,080225-080527 3幅干涉圖，局部區(qū)域的相位出現(xiàn)顯著的跳躍，而三者的成像間隔均小于1年，地表形變量很小，不可能導致如此大的相位梯度，因此，相位的不連續(xù)應該是大氣延遲在空間上分布不均引起的，同時，3幅干涉圖對應同一景副影像，說明主要是受副影像大氣效應的影響。對于這種空間上存在較大差異的干涉圖，進行多項式擬合時必然會引入較大的誤差，為了保證結(jié)果的可靠性，將該3幅干涉圖剔除，利用其余24幅干涉圖進行求解。

基于相干系數(shù)閾值篩選出一定數(shù)目的高相干點。利用24幅干涉圖對應的相干圖，采用平均相干系數(shù)和最小相干系數(shù)雙閾值選擇高相干點，設置平均相干系數(shù)閾值為0.45，最小相干系數(shù)閾值為0.25，得到初始相干點的個數(shù)為47 415，平均每平方千米接近300個點。

由于InSAR所得的相位值為相對量而非絕對量，因此需要將24景干涉圖統(tǒng)一到相同的基準。本文的做法是從相干點集中選擇一個點作為參考點，參考點須滿足以下兩個條件：1)參考點的相位與所有相干點相位的平均值相差不大；2)在影像中心附近。這樣一方面最大程度保證參考點上不存在粗差，另一方面所有點與參考點差分也可以消除一部分公共誤差。

利用多項式擬合方法，依次估計24幅差分干涉圖對應的大氣延遲相位，并將其從干涉圖中去除。

干涉對的垂直基線小于1 km，根據(jù)前文分析，DEM誤差引起的相位在1.2 rad以內(nèi)。干涉組合中最長時間基線為782 d，而近年來上海市區(qū)年均沉降速率在1 cm/年以內(nèi)，最大年均沉降速率不超過3 cm/年，投影到視線方向不足2 cm/年，因此，地表形變所引起的相位應在2 rad以內(nèi)。綜上所述，線性模型項(形變相位和DEM誤差相位)不會超過半個周期。因此，當誤差不是很大的時候，可以認為干涉相位中不含有整周模糊度，而含有整周模糊度的點被視為粗差點，利用粗差探測的方法予以剔除，即：

|δφx-mean(δφ)|>Tφ

(3)

其中，mean(δφ)為某一幅干涉圖中所有相干點相位的均值；Tφ為閾值，可取為3倍中誤差或先驗值。經(jīng)過一次多項式擬合，多數(shù)干涉圖的均方根在0.5以內(nèi)，而少數(shù)干涉圖的均方根接近1 rad，本文選定2 rad作為粗差探測的閾值。利用剔除粗差點的相干點集，重新進行多項式擬合，估計大氣延遲相位并去除，再次剔除粗差點，以上過程重復執(zhí)行，直至沒有粗差點。最終剩余21 488個相干點。利用剔除粗差點后的相干點集，基于最小二乘方法求解各點相對參考點的線性形變速率，如圖2所示。

所有相干點的平均年沉降速率在-10 mm/年左右，平均中誤差為2 mm/年，最大沉降速率值約-30 mm/年。利用水準數(shù)據(jù)對校正后的結(jié)果進行檢驗，二者結(jié)果對比見圖3，其中有三個點的差值超過7 mm，其余點的差值均在4 mm以內(nèi)，具有較好的一致性。

5 結(jié)語

本文基于最小二乘的相干點分析方法，利用13景ALOS PALSAR數(shù)據(jù)，反演了實驗區(qū)的線性形變速度，并用實驗區(qū)內(nèi)的水準數(shù)據(jù)進行了驗證。實驗結(jié)果表明，該方法具有一定的有效性。

本文利用多項式擬合消除空間相關(guān)誤差(如大氣效應)的影響，但該方法的有效范圍有限，同時多項式擬合對形變量也會造成一定的影響，因此實驗結(jié)果中仍包含一定的不確定因素。

參考文獻：

[1] 范景輝，郭華東，郭小方，等.基于相干目標的干涉圖疊加方法監(jiān)測天津地區(qū)地面沉降[J].遙感學報，2008，12 (1):111-118.

[2] 方 正，武健強，趙建康.長江三角洲地區(qū)地面沉降監(jiān)測[J].上海地質(zhì)，2003(2)：1-4.

[3] 葛大慶，王 艷，郭小方，等.基于相干點目標的多基線D-InSAR技術(shù)與地表形變監(jiān)測[J].遙感學報，2007，11(4):574-510.

[4] 劉國祥，陳 強，丁曉利.基于雷達干涉永久散射體網(wǎng)絡探測地表形變的算法與實驗結(jié)果[J].測繪學報，2007，36(1):13-18.

[5] 王 艷，廖明生，李德仁，等.利用長時間序列相干目標獲取地面沉降場[J].地球物理學報，2007,50(2):598-604.

[6] 吳宏安.基于改進最小二乘DInSAR技術(shù)反演大區(qū)域地面沉降研究[D].北京：中國科學院遙感應用研究所，2009.