利用PALSAR數據提取地面沉降短基線集方法研究

王成

(重慶市勘測院，重慶 400020)

利用PALSAR數據提取地面沉降短基線集方法研究

王成?

(重慶市勘測院，重慶 400020)

短基線集DInSAR方法可以很好地克服數據相干性差的困難，本文在分析短基線集方法的基礎上，對大氣相位和軌道誤差相位的影響進行處理得到較為精確的形變相位。采用天津實驗區從2007年1月17日～2009年10月25日的17景PALSAR單視幅圖像提取地面相對沉降信息，將DInSAR與水準測量數據進行融合得到外部精度為3.8mm的沉降結果，同時確定了天津地面的沉降漏斗的分布和絕對形變量。

合成孔徑雷達干涉測量；短基線集；地面沉降；軌道誤差；大氣相位誤差

1 引 言

合成孔徑雷達差分干涉測量技術(Differential SAR Interferometry，DInSAR)是在InSAR技術基礎上進一步發展而來的，它通過比較地表形變前后干涉圖之間的相位差來監測地表目標的位移[1]。具有全天時、全天候、廣覆蓋、高分辨率及高地表垂直形變靈敏度等優點。但在長時間微小地表形變監測中，時空去相干和大氣相位延遲的影響制約了它的發展[2]。

Berardino[3]等提出的短基線集算法，可以很好地解決DInSAR長時間序列干涉紋圖中缺乏大面積連續高相干區域的問題。它是在常規DInSAR基礎上限制時空基線的大小后自由組合干涉對，利用短基線干涉圖提取高相干區域，從而有效地進行形變反演。

本文收集了天津地區自2007年1月17日～2009年10月25日的16景PALSAR單視復影像，基于GAMMA軟件平臺進行常規二軌差分處理得到136組干涉像對。由于是對長時間地面形變進行監測，因此引入了短基線集方法，分別選取其中干涉質量好的垂直基線較短的78組干涉像對采用短基線集DInSAR方法進行形變反演，通過基線精化處理和擬合最佳軌道面的方法來削弱軌道誤差影響，采用高斯低通濾波的方法來減弱大氣相位的影響，通過阻尼最小二乘算法分離平均線性形變速率和DEM殘余誤差的改正量，得到最終的垂直形變量。

最后利用短基線集DInSAR方法所得到的垂直形變量與7個已知水準點同時段的沉降量進行融合，使得短基線DInSAR方法監測得到的相對形變量轉化為絕對形變量。PALSAR數據監測結果的精度為3.8 mm，高相干點位密度為666個/km2。將監測結果繪制成沉降等值線并疊加到Google Earth地圖上，由此分析確定天津地區的地面沉降漏斗中心位置分布及其絕對沉降量。

2 實驗區概況

地面沉降是許多沿海城市面臨的頑癥，天津地區的地面沉降現象在中國具有代表性，其主要原因是地下水過量開采。由于地下流體資源過量開采和構造沉降活動，在寶坻斷裂和薊運河斷裂以南約7 300 km2的廣大地區，出現不同程度的地面沉降現象，南部和濱海地區尤為明顯，并且與河北省沉降區連成一片，在這一大面積的沉降范圍內，又形成了市區、塘沽、漢沽、大港及海河下游工業區形成了5個二級沉降中心。最近幾年隨著區縣、鄉鎮經濟的發展，出現了武清、津南、靜海等新的沉降中心，其發展速度比當年老的沉降區域還要快，局部沉降速率超過80 mm/y[4]。

常規地面沉降監測技術主要有水準測量和GPS測量方法，但二者的空間點位稀疏、周期長的缺點大大限制其對區域地面沉降發展趨勢的研究。DInSAR技術具有大范圍、全天候、高精度、高分辨率地測量地表三維空間的微小變化，在地表變形的監測方面顯示出前所未有的優越性。

3 短基線集DInSAR方法與原理

設有按時間序列t0，……，tN獲取的N+1幅單視復數SAR影像，將它們以任意影像為主影像進行配準，設定垂直基線閾值(一般小于400 m)，垂直基線小于該閾值的影像作為一個基線集，并將其劃分為L組，每組內的影像進行差分干涉處理，最終可得到M幅差分干涉圖，且有:

現以t0時刻為參考時刻，任意時刻ti(i=1，…，N)相對于參考時刻的差分相位φ(ti)為未知數，數據處理獲取的差分干涉相位δφ(tk)(i=1，…，M)為觀測量，不考慮大氣延遲相位和數字高程模型(DEM)誤差的影響，對第k幅(k=1，……，M)差分干涉相位圖中的任意像元(x，r)可以組成如下方程:

上式中λ為雷達波長，d(ti，x，r)和d(tj，x，r)分別為ti和tj時刻像元(x，r)相對于參考時刻t0的雷達視線方向(line of sight，LOS)地表形變。將式中的φ(ti，x，r)變為相鄰影像像元沿LOS方向的平均相位變化速率v，有:

其中IE和IS分別為按時間順序排列的主影像序列和從影像序列，滿足IEj＞ISj，?j=1，……，M化簡即有:

其中B是M×N矩陣，B為奇異矩陣，可采用阻尼最小二乘算法進行求解(Levenberg，1944)，可得到線性形變速率為:

式中γ是阻尼因子，阻尼因子γ越大，分辨就越低；另一方面，當γ增大時分母就越大，從而減小了模型參數的誤差[5]，阻尼因子γ對反演結果有一定影響，本文對阻尼因子進行實驗分析選取了一個較為合適的數值，阻尼因子為0.25進行阻尼最小二乘的反演計算。

將式(3)帶入式(2)得:

4 數據預處理與誤差消除

4.1 實驗區數據準備與預處理

針對天津實驗區，選取了16景ALOSPALSAR從2007年1月17日～2009年10月25日高分辨率模式單視復影像，所獲取為升軌時期的數據。將其中任意兩幅圖進行組合形成120個像對，使用GAMMA軟件利用常規二軌差分干涉測量方法進行處理，對于PALSAR還需要進一步精化基線來削弱線性相位的影響。利用美國航天飛機雷達地形測繪任務(SRTM)3″分辨率(90 m)DEM去除地形相位和地平相位，利用最小費用流方法對干涉相位進行相位解纏，從而得到差分干涉圖。

處理得到的差分相位由形變相位、殘余軌道誤差相位、DEM殘余誤差相位和大氣相位等組成。由In-SAR原理可知垂直基線長度與高程模糊度成反比，基線越短，高程模糊度越大，地形對相位的影響也就越小，所以需要選擇時空基線較短的干涉像對組合。通過實驗分析，最終選取了4 000 m以下的干涉質量較好的78組干涉像對，采用短基線集方法進行形變反演分析，對以上誤差相位進行分離并予以削弱，得到最終的形變相位。所選取的干涉對的時空基線分布如圖1所示。

圖1 PALSAR時空基線分布圖

4.2 衛星軌道誤差的消除

SAR衛星軌道誤差是指用來計算的衛星星歷中存在誤差引起的基線計算不精確而導致后續的干涉相位誤差。這里首先采取精化基線的方法對軌道誤差進行一定的改正，如圖2所示。為了進一步削弱殘余軌道誤差的影響，然后通過線性擬合的方法來對解纏后的干涉相位進行擬合，得到最佳擬合軌道面，即殘余的軌道誤差，并將其從干涉圖中去除。按照如下方程來計算每個干涉圖的最佳的軌道面[6]:

式中，[x，y]是每個像素點的在像素坐標系下的坐標，u，v是斜率參數，w是截距，如圖3所示。

圖2 PALSAR基線精化處理

圖3 PALSAR軌道誤差改正

4.3 消除DEM殘余誤差

此外，在實際處理過程中需要考慮DEM殘余誤差的影響，則式(5)表示為:

Dv+C·ε=△φ(8)

式中cT=[(4π/λ)(B⊥1/r sinθ)，…，(4π/λ)(B⊥N/ r sinθ)]，λ為雷達波長，B⊥為垂直基線，r是衛星到地面間的距離，θ是衛星視角。通過阻尼最小二乘算法計算出DEM殘余誤差的改正量ε，其大小在[-10，10]m范圍內，結果如圖4所示。同時可得到沿雷達視線方向的線性形變速率v乘以時間就得到雷達視線方向的形變量再除以cosθ，從而轉化為垂直形變量。

圖4 PALSAR DEM殘余誤差改正量

4.4 大氣相位影響

短基線集方法是采用線性模型，為了進一步消去殘余相位中的大氣相位和非線性形變相位的影響，可以在線性模型的基礎上進行適當的濾波處理[7]。根據大氣相位空域的低頻特性，本文采用在空間上實施高斯低通濾波處理，移除“形變信號”中上的高頻部分。如圖5所示，是對其中一組數據進行濾波的前后對比圖。

圖5 PALSAR濾波前后對比圖

5 結果分析

5.1 獲取天津地面相對形變量

此次形變反演過程中相干閾值設置為0.67，得到高相干點個數為599 728，影像面積約為900 km2，即點密度為666個/km2。DInSAR測量結果是相對于基準點而言的，該基準點是研究區內所有其他點形變量的參考基準，所以需要選擇形變較為穩定的點作為基準點。同時選擇干涉質量較好的2007年7月20日主影像為沉降基準時刻，將另外的16景影像與之差分，從而得到16幅天津實驗區地面相對形變量變化圖，如圖6所示。圖中“+”表示的參考基準點，是研究區內相對穩定的點，負值表示沉降。通過地理編碼，將沉降結果由制圖坐標系轉換到WGS84坐標下，可以獲取實驗區內自2007年01月17日～2009年10月25日的累積相對形變量圖，如圖7所示。圖中“★”表示參考基準點，“?”表示累積沉降量最大的點，我們定義為A點，其相對形變量不超過15 cm，將其相對于參考基準點的沉降趨勢表示出來，如圖8所示。

圖6 天津實驗區地面相對形變量變化圖

圖7 2007年01月17日～2009年10月25日天津地面相對形變量累積圖

圖8 A點相對于基準點的沉降趨勢圖

5.2 DInSAR與水準數據的融合

通過短基線集DInSAR方法獲得的是大范圍相對形變量，而通過水準測量可以獲得高精度的局部范圍的絕對沉降信息。為了獲取大范圍的高精度地面沉降信息，則可以將二者進行數據融合。由于DInSAR測量的基準與水準測量基準不一致，所以需要進行基準的統一。通過數據分析，兩者監測結果關系可采用下列函數表達[8]，令:

式中y為水準測量結果，x為雷達干涉測量結果。

基于2007.01～2009.10的水準測量結果和2007.01.17～2009.10.25的DInSAR結果(由于水準測量的周期較長，所以不能保證與DInSAR數據的時間段完全重合)，對水準與DInSAR結果進行了融合。首先按線性關系對兩者進行了擬合分析，剔除偏差較大的點后，再用7個點分布均勻的數據對兩者重新擬合，得到新的函數關系結果如圖9所示。

圖9 DInSAR數據與水準數據的線性關系圖

用于擬合的點的校正后的DInSAR結果與水準數據的標準差為3.8 mm。則可利用水準與DInSAR之間的函數關系對DInSAR形變結果進行校正，可以得到整個區域的經過水準校正的地面沉降成果，從而獲取天津實驗區的從2007年01月17日～2009年10月25日的絕對形變量，結果如圖10所示。5.3 實驗區沉降結果分析

圖10 絕對累積形變量

將經過水準數據修正之后的InSAR數據疊加到GoogleEarth地圖上，各相干點目標沉降累積形變量的大小以不同的顏色分級表示，結果顯示DInSAR數據點主要分布在天津市中心的范圍，越靠近郊區的范圍的相干點數據越為稀疏(主要為農田或魚塘)，如圖11所示。將沉降結果繪制成沉降等值線疊加到GoogleE-arth地圖上，如圖12所示。

從圖5中沉降趨勢分析，自2007年1月～2009年10月天津市主城區沉降速率較為穩定，累積沉降量約為4 cm左右。沉降速率從中心到郊區逐漸變大，說明主城區的地面沉降已經得到了有效的遏制，而郊區的沉降速率仍然偏大，且發現局部地區存在明顯的沉降漏斗。

圖11 天津實驗區DInSAR相干點數據疊加到Google Earth圖

圖12 DInSAR沉降等值線疊加到Google Earth圖

從圖11和圖12可以看出沉降漏斗主要分布在武清區、北辰區以及西青區。研究區域內最大(負值代表下沉)沉降漏斗中心位于西青區南河鎮附近該沉降漏斗中心的累積沉降量約為26.6 cm和大寺鎮附近該沉降漏斗中心的累積形變量達24.6 cm。武清區的沉降漏斗位于王慶坨鎮附近，累積沉降量達18.5 cm。還有位于北辰區的西南方向的上河頭鎮附近的沉降漏斗，漏斗中心的沉降量約為21.6 cm。郊區沉降較大的原因是近幾年天津工業快速向郊區轉移和發展，郊區地下水抽取日益嚴重，導致地面沉降加速[9，10]。

6 結 語

本文針對天津地區地面沉降監測試驗主要做了以下幾項工作，首先是利用16景PALSAR單幅視影像進行常規二通差分干涉測量方法處理，采用精軌數據去除較大的軌道誤差，利用SRTM DEM去除地形相位和地平相位的影響，利用最小費用流方進行相位解纏，從而得到差分干涉圖。然后根據短基線集方法的原理，選取78個干涉對進行時間序列分析，采用擬合最佳軌道面的方法消除殘余軌道誤差的影響，并利用阻尼最小二乘算法求解DEM殘余誤差改正量以及所需要得到的平均線性形變速率。為了減弱大氣相位的影響，這里做了空間低通濾波的處理，得到天津地區2007年1月～2009年10月的相對累積形變。最后為了獲得絕對形變量，與7個水準點測量數據進行了融合，獲得外部精度為3.8 mm的沉降結果，并確定天津地區沉降漏斗的分布及其絕對沉降量。

基于上述工作得到了以下幾點認識:

(1)短基線集方法算法較PS方法簡單，本文的特點在于在該算法基礎上對軌道誤差和大氣相位的影響做了一定的分析和削弱；

(2)在利用阻尼最小二乘算法進行反演過程中，選擇較為合適的阻尼因子是關鍵，一般選擇在0.01～1之間；

(3)在水準數據與DInSAR數據融合過程中，由于水準點處沒有高相干點數據，其水準點的DInSAR數據是通過插值得到，對最終結果有一定影響；

(4)利用水準數據可以將DInSAR的相對沉降信息轉為絕對沉降信息，同時兼顧水準測量高精度的優點和DInSAR大面積、分辨率高的優點，能夠在地面沉降監測的實際應用工作中發揮一定作用。

[1] 周春霞.星載SAR干涉測量技術及其在南極冰貌地形研究中的應用[D].武漢:武漢大學，2004.

[2] 王超，張紅，劉智.星載合成孔徑雷達干涉測量[M].北京:科學出版社，2002.

[3] Berardino，P.Fornaro，G.Lanari，R.et al.A new algorithm for surface deformationmonitoring based on smallbaseline differential interferograms[J].IEEE Transactions on Geoscience and Remote Sensing，Vol.40，No.11，2375～2383，2002.

[4] 董克剛，王威，于強等.天津地面沉降區土水比論述[J].水文地質工程地質，2008，19(3):54～59.

[5] 崔希璋，於宗儔，陶本藻等.廣義測量平差(第二版)[M].武漢:武漢大學出版社，2009.

[6] 湯益先.基于相干目標的DInSAR方法及其地表沉降應用研究[D].北京:中國科學研究院，2006.

[7] 張紅，王超，吳濤.基于相干目標的DInSAR方法研究[M].北京:科學出版社，2008.

[8] 張詩玉.干涉圖大氣效應校正方法研究及InSAR在天津地面沉降監測中的應用[D].武漢:武漢大學，2009.

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

[10] 張詩玉，李淘，夏耶.基于InSAR技術的城市地面沉降災害監測研究[J].武漢大學學報·信息科學版.2008，33 (8):850～853.

Study on Using PALSAR Data Extract Land Subsidence w ith Short Baseline Subset M ethod

Wang Cheng
(Chongqing Survey Insitute，Chongqing 400020，China)

The short baseline subset DInSAR technology can be very good to overcome the difficulties of poor coherence data，this paper is based on the analysis of shortbaseline subsetmethod toweak the atmospheric phase and orbital phase error，then gainmore accurate deformation phase.This paper collected 17 ALOSPALSAR single complex images since the January 17，2007 to October 2009，25which cover the Tianjin area，Intergration of DInSAR and leveling data to get external precision as the settlement results of 3.8mm.Finally，Tianjin ground subsidence funnel distribution and absolute variables are determined.

DInSAR；short baseline subset；groud subsidence；orbit error；atmospheric phase error

1672-8262(2013)05-81-06

P237

A

2013—03—18

王成(1988—)，女，工程師，主要從事航空攝影測量與遙感生產工作。

國家科技支撐計劃(2011BAH12B07-03)