PS-InSAR和SBAS-InSAR技術對昆明主城區地面沉降監測的對比分析

范 軍，左小清，李 濤，陳乾福

(1.昆明理工大學 國土資源工程學院，云南 昆明 650093；2.國家測繪地理信息局衛星測繪應用中心，北京 100830)

地面沉降是一種十分嚴重的地災現象，與人們的生活息息相關。目前根據文獻顯示， 2001年姜朝松等人利用1979—1998年水準數據對昆明地面沉降的發展過程及其特征進行研究[1]；2003年孟國濤利用1987—1998年四期水準數據對昆明南市區地面沉降進行研究[2]；2004年薛傳東[3]對昆明市區地面沉降進行深入研究。精密水準測量方法精度滿足要求，但是存在工作力度大、周期長、成本高、效率低和覆蓋范圍小等缺點，同時很難滿足時間、空間尺度的需要。2007年高照中結合“3S”技術建立自動監測系統對昆明市地面進行高精度沉降監測[4]，此方法對比傳統水準測量，GPS則具有覆蓋范圍廣、精度高、速度快等優點，但作業周期長、測量結果不連續難以高效監測與分析整個城市區域性沉降。

合成孔徑雷達差分干涉測量(DInSAR)近年來發展起來的一種監測地表沉降的有效手段，能夠對地表進行全天候、大范圍、亞厘米級精度的監測，但是傳統DInSAR易受時空失相干與大氣效應的影響，較難完成高精度長時間間隔的地表監測，特別是永久散射體雷達干涉測量技術(PS-InSAR)和短基線干涉測量技術(SBAS-InSAR)的提出，能夠充分有效地利用多景雷達影像，利用時序分析技術有效解決傳統DInSAR技術中一些問題，提高DInSAR地面監測的精度[5-6]。目前，PS-InSAR和 SBAS-InSAR技術已經廣泛應用在城市地面沉降監測。2016年尹振興選取2007—2010年SAR數據基于SBAS-InSAR技術對昆明地面沉降進行監測[7]，在前人研究基礎上找出多處新的沉降漏斗。目前，采用實時SAR數據探究昆明地區沉降情況較少，也很少有人應用PS-InSAR技術監測昆明地區沉降。因此，本文選取了2014—2017年29景升軌Sentinel-1A數據基于PS和SBAS兩種技術對昆明市地面沉降進行監測并作對比分析。

1 PS-InSAR方法

2001年，Ferretti等人首次提出并驗證了PS-InSAR技術[8]。即在一組長時間序列相干SAR影像內一般選取高相干性且較為穩定的永久散射體(PS)，這些PS點受時空相干的影響小，分析這些PS點的相位組成，利用各相位分量的時頻特性，分離大氣相位，從而獲得精確的地表形變信息。近年來，利用PS-InSAR方法監測地面沉降已取得良好的研究成果，在北京、天津、上海等地獲得的沉降結果與外業實測結果一致且能達到較到精度，證明該方法的應用潛力。

PS-InSAR技術基本原理主要是獲取高質量PS點，通過InSAR技術分析以獲取PS點高精度的時間序列， 本文主要利用振幅離差法初選 PS 備選點，再對選取的PS點進行相位分析，解算時間相干系數γx。Ferretti等人在2001年提出振幅離差法[8]，主要利用振幅離差和相位標準偏差的關系建立提取初始相干點閾值，振幅離差指數

(1)

其中，σA，μA為振幅的標準偏差和均值，σφ，DA分別表示相位標準偏差和振幅離差，在這里通常選取0.25作為振幅離差閾值，小于0.25的點，認為相位離散度較小，可選為PS備選點，在得到PS備選點后，分析其相位穩定性來篩選最佳的PS點基本模型，式(2)表示第i幅干涉條紋圖中第x個PS備選點的干涉相位為

φx,i=wrap{φD,x,i+φA,x,i+Δφh,x,i+φN,x,i}.

(2)

其中，φD為雷達視線方向的形變相位，φA為大氣延遲相位，φh為高程相位，φN為噪聲相位，形變相位、大氣相位為空間相關相位，高程相位大部分為相關相位，小部分為非相關相位這里不考慮，只有噪聲相位為空間不相關相位，采用Hooper提出的低通濾波法來獲取前面3個分量的相位和，根據式(2)對PS備選點直接進行相位濾波求解噪聲相位，為此采用式(3)時間相干系數閾值作為評價最終PS點提取的標準[9]。

(3)

2 SBAS-InSAR方法

SBAS方法是由Berardino和Lanari等研究人員于2002年正式提出的一種InSAR時間序列分析方法[11]，該方法通過多景SAR數據自由組合所有短基線干涉對，得到一系列短基線差分干涉圖，這些差分干涉圖能夠有效地去除基線失相關現象，不僅連接所有差分中短基線集，同時連接較大空間基線中分開的孤立SAR數據集，來提高觀測數據時間采樣率，在求解形變速率時，SBAS方法采用Delaunay MCF方法進行3D解纏，Goldstein方法對干涉圖進行濾波處理，最后奇異值分解(SVD)算法來獲得高相干點的沉降速率和時間序列，在此基礎上Hooper在2008年提出了新的SBAS方法，通過提取SDFP點求解形變，SDFP點提取方法同PS點一樣[9，12]。

SBAS-InSAR技術基本原理就是將同一地區多景SAR影像進行兩兩組合配準得到短基線干涉對，采用SVD法生成相干點平均形變速率圖和時間序列，要根據短基線原則生成差分干涉圖，得到干涉相位。假設本文選取同一地區N+1幅SAR影像，按照時間順序(t0，t1,…，tN)排列，選取其中一幅影像為主影像，其他SAR影像配準到這幅主影像上，在短基線干涉組合條件限制下，生成M幅差分干涉圖，且M滿足：

(4)

對于任意干涉圖j在tA與tB時間內獲取的影像進行干涉生成(tB>tA)，方位向坐標x和距離向坐標r的像素的干涉相位可以表示為

δφj(x,r)=φB(x,r)-φA(x,r)≈

(5)

式中：j∈(1,2，…,M)；λ為雷達信號的波長；φB(x,r)和φA(x,r)分別為tB與tA時刻SAR影像的相位值；d(tB,x,r)和d(tA,x,r)分別為tB與tA時刻相對于d(t0,x,r)=0的雷達視線方向累積形變量。

為了更好的表示地面沉降時間序列，可將式(5)相位值表示為兩個獲取時間段的平均相位速率vj與該時間段之間的乘積

(6)

則第j幅干涉圖的相位值可表示為

(7)

即第j幅干涉圖的相位值表示為各時間段平均速率在主、從影像時間間隔上的積分，可表示為新的矩陣方程

Bv=δφ.

(8)

式中：B是一個M×N的系數矩陣，由于SBAS采用多個主影像時空基線獲取干涉對的方法，因此矩陣B產生秩虧，運用SVD法對矩陣B進行解算，得到矩陣B的廣義逆矩陣，進而得到速率的最小范數解，最后在各個時間段對速率積分就可得到各個時間段的累積沉降量[13]。

3 實例分析

3.1 研究區概況及實驗數據

昆明市區坐落于昆明晚新生代斷陷盆地內，廣泛發育第四紀松散沉積層，以湖沼相粉砂和軟粘土為主，夾多層淤泥和草煤層等松軟沉積物[3]。近年來，由于昆明市城市與軌道交通建設的飛速發展，主城區發生大規模地面沉降現象，沉降范圍日益擴大，新的沉降中心不斷產生，為此，本文選擇升降Sentinel-1A(S1A)數據的昆明主城區為研究區域，對實驗數據進行裁剪得到本文研究范圍如圖1所示，中心經緯度為25°00′N，102°43′E，覆蓋面積約為400 km2，研究區主要覆蓋官渡區，五華區，盤龍區，西山區等部分區域。

圖1 研究區地理位置

本文從歐空局共獲取2014-12-06—2017-02-05間的29景升軌(Track-26)間Level-0級原始S1A數據，S1A衛星在2014年發射，入射角為39.5°，C波段(波長5.6 cm)，分辨率5 m×20 m，觀測模式為干涉寬幅(IW)，極化方式為單極化(VV)，SAR影像大小為7800行×2130列。

3.2 數據處理

3.2.1 PS干涉處理

首先對29景SAR影像經過配準干涉去平等步驟進行差分干涉處理，生成28個差分干涉對，在去除地形相位時，引入了90 m分辨率的SRTM DEM數據。以2016-05-05為主影像，如圖2所示，空間基線分布最大116 m，最小3 m，平均為34 m，可以看出基線分布較小較穩定，高程模糊度小，DEM誤差對形變的影響比較小，相干性也會比較高，故地形誤差引起的誤差相位將非常小，時間基線最大值為540 d，最小24 d。

本文采用振幅離差指數閾值與時間相干系數閾值相結合的方法[8-9]，首先將幅度離差指數閾值設置為0.65，作為初始PS候選點，將像素采樣到30 m格網空間，濾波窗口大小設為15×15，從而使其滿足空間相關條件，通過設置相干性系數γx的迭代收斂條件，進一步篩選最終得到43 441個PS點，研究區域總面積約400 km2，平均每平方公里約150個PS點，本文研究昆明主城區，PS點密度較高，PS點大多分布在建筑物、橋梁、裸露的巖石等處，且水泥路面、滇池、綠色植被覆蓋綠化區域基本沒有 PS 點，可見識別出PS 點質量較好。

3.2.2 小基線干涉處理

SBAS方法同樣選取2016-05-05為超級主影像，并對所有影像進行4×1(距離向×方位向)的多視處理。試驗中設置時間基線和空間基線閾值分別為550天和±250 m，進行多主影像干涉對自由組合，剔除質量較差的干涉對，最終選取361個較高質量的差分干涉對，如圖3所示黃色代表超級主影像，綠色代表SAR影像，線段代表干涉對。其中最大時間基線540 d，最小時間基線12 d，平均時間基線為235 d，最大空間基線為227.718 m，最小空間基線為3.028 m，平均空間基線為55.560 m。由于所有干涉圖的空間基線分布較小、相干性較穩定，所構集合網絡更強，保證后續計算結果的可靠性。

圖2 PS干涉處理時空基線分布

4 結果分析

4.1 沉降結果分析

經過上述PS和SBAS數據處理，得到2014-12—2017-02昆明市區主城區地面年沉降速率圖(底圖為昆明市Eoogle地圖)如圖4所示，本次沉降量選取3個相對較大的區域沉降分別記為A、B、C，獲取各區域地面沉降年平均速率(沿雷達視線方向)，紅色代表遠離衛星，表示地面下沉，藍色代表接近衛星，表示地表上升[14]。

A區域位于官渡區，由PS和SBAS監測結果如圖5所示，沉降漏斗主要位于滇池北岸以及東北岸，分布范圍廣泛，出現多處沉降漏斗，其中沉降最嚴重為彼岸小區-楓林盛景小區-星宇園小區-子君欣景小區-羊甫村-義路村-廣衛村-螺螄灣北等多處沉降漏斗連成一體，構成研究區域內最大的漏斗沉降區域，范圍有明顯加劇的趨勢，其中沉降漏斗形成面積已達到5.63 km2，兩種方法最大年沉降速率分別為-39.580 mm/a和-37.405 mm/a，最大沉降量達到85 mm。昆明地區的地面沉降主要以軟土沉降為主，地面受到建筑用地施壓造成軟粘土中的水被擠出，土層被壓密實，地基下沉；其次，這些沉降區域位于昆磨高速公路和東環城高速公路交界的廣衛立交橋的西北處附近，此處有正在建設的昆明地鐵4號線，其中螺螄灣北為4號地鐵站，年沉降速率分別為-32.510和-31.801 mm/a，這些區域受到地鐵隧道開挖和相應的工程降水措施會直接擾動軟土，同時地面上下車荷載交通流量引起的土體變形和地層移動可能影響地表結構，共同觸發了該路段以及其周邊區域的沉降[15]。小板橋位于云大西路和昆磨高速公路與廣昆高速公路交界的鳴泉村立交橋西南處附近，年沉降速率為-19.360 mm/a和-20.174 mm/a，沉降原因不僅受到地質因素的影響，同時還受到相對集中的居民地和交通流量荷載的影響，造成該區域沉降。曉東村、雨龍村、小板橋、竹園村和陳旗營等沉降中心靠近地鐵1號線已連成一片，雨龍村在沉降區域表現最為嚴重，年沉降速率分別為-17.870和-17.922 mm/a。

圖4 研究區年平均沉降速率圖

圖5 A區域年平均沉降速率圖

B區域位于官渡區，為昆明老城區，靠近滇池，由監測結果如圖6所示，福保路附近一帶且出現了多處沉降漏斗，沉降區域有明顯向市中心擴張的趨勢，從南向北分別為丁姚村，韓家村，龔家村，金家村，六甲村，葉家村，羅家村，五甲河，龍馬村，羅衙村和陳家營村等沉降區域，其中六甲村出現了多處沉降漏斗，韓家村-金家村-龔家村等沉降區域連成一片并產生了多處新的沉降中心，構成了研究區域較為嚴重的漏斗沉降區，最大年沉降速率分別達到-27.240 mm/a與-25.809 mm/a，最大累積沉降量達到-59 mm，該地區沉降主要因為韓家村，金家村，龔家村等沉降區域靠近滇池，滇池湖水面是盆地最低的徑流基準面，大量地下水向滇池方向徑流和排泄，其次滇池附近地區主要存在大量軟土、飽水粉砂層，并夾有多層淤泥、草煤層等松軟沉積物，這些地層形成時間短、固結度低、孔隙率高、含水量大，導致滇池周邊地區地基太軟，隨著軟土下沉而出現明顯的沉降現象[2]。此外，由于B區域位于老城區，從理論上分析，一般可能是鄰近地區有人類工程活動影響，如深基坑開挖和地下水超采、城市建筑加載等。

圖6 B區域年平均沉降速率

C區域位于西山區，由監測結果圖7可知，沉降漏斗出現在春苑小區和棕樹營兩個區域，最大年沉降速率分別為-21.780 mm/a和-20.944 mm/a，而其他大部分區域無明顯升降變化，則認為是地面相對穩定。

圖7 C區域年平均沉降速率

4.2 對比分析

4.2.1 整體沉降情況分析

由于缺少同期測量數據，將利用PS-InSAR技術所得的PS點和利用SBAS-InSAR技術得到SDFP點進行統計分析，其中PS點為43 441個，SDFP點為61 860個，由圖8可以看出PS點與SDFP點在各個數值上的分布趨勢非常相似，說明兩種方法得到的年沉降速率相符，但是在這兩種方法選擇備選點時由于算法局限性等因素得到的結果還是有些細微的差別[16]。

圖8 年平均沉降速率分布

4.2.2 相關性分析

圖9 PS與SBAS年平均沉降速率相關性圖

為了進一步對兩種方法精度分析，本文選取A區域，篩選972個同名(具有相同的地理坐標)相干性較高PS點和SDFP點，根據篩選得到的PS點和SDFP點上的沉降速率進行對比得到相關性分析結果如圖9所示，可以看出PS點和SDFP點的沉降速率分布比較一致，兩者線性相關性R達到了0.893以上，說明兩種方法得到的年沉降速率具有高度的一致性，因此利用PS-InSAR和SBAS-InSAR方法監測地面沉降的可行性。

4.2.3 剖面圖分析

結合沉降剖面圖10—圖11對A區域和B區域進行分析，剖面線的方向為S-N方向，通過兩種方法結果顯示，A區域最大年沉降速率達到-39.580 mm/a和-37.405 mm/a， B區域最大年沉降速率達到-27.240 mm/a與-25.809 mm/a，如下圖所示，PS-InSAR和SBAS-InSAR兩種技術監測得到的剖面圖沉降趨勢基本一致。

4.2.4 時序分析

選擇A區域和B區域中最大沉降漏斗作時序分析，相對于起始時間2014-12-06的累積沉降量的變化情況，如圖12—圖13所示，兩種技術得到A區域最大漏斗累積沉降量分別為85 mm和81 mm，B區域最大漏斗累積沉降量分別為59 mm和55 mm，兩種技術得到沉降量差值均為4 mm。

綜上所述，證明PS-InSAR和SBAS-InSAR兩種技術監測結果的準確性、可靠性和一致性。

圖10 A區域沉降漏斗剖面

圖11 B區域沉降漏斗剖面

圖12 A區域最大沉降漏斗累積沉降量

圖13 B區域最大沉降漏斗累積沉降量

5 結束語

由于未收集到同期的實測水準數據，通過使用PS-InSAR和SBAS-InSAR技術分別對昆明主城區的地面沉降進行監測，得到了昆明主城區沿LOS方向年平均沉降速率圖，通過對比分析兩種技術的沉降趨勢基本吻合。試驗證明：研究區共出現3處較大的形變區，其中A區域螺螄灣北，小板橋南一帶地面沉降基本呈線性沉降規律，年平均沉降速率達到-39.580 mm/a，最大累積沉降量達到85 mm。B區域靠近滇池一帶，福保路周邊多處區域出現沉降漏斗連成一體，年平均沉降速率達到-27.240 mm/a，最大累積沉降量達到59 mm，這些嚴重的沉降現象不僅與工程地質環境密切有關，近幾年來昆明市城市和軌道交通建設的飛快發展，居民區的興建與密集道路網增大地面荷載力，地鐵施工造成地下降水量增大，導致廣泛分布的粘松散型軟土層很容易受到外力作用產生地面沉降，這成為昆明主城區地面沉降的主要原因。

參考文獻：

[1] 姜朝松,邵德晟,樊友心,等. 昆明市地面沉降發展過程及其特征[J]. 地震研究,2001(1):55-60.

[2] 孟國濤.昆明南市區地面沉降研究[D].昆明：昆明理工大學,2003.

[3] 薛傳東,劉星,李保珠,等.昆明市區地面沉降的機理分析[J].中國地質災害與防治學報,2004(3):51-58.

[4] 高照忠,羅志清,魏海霞. 3S技術在昆明市地面沉降監測中的應用[J]. 地礦測繪,2007(1):27-29.

[5] CASU F, ELEFANTE S, IMPERATORE P, et al. SBAS-DInSAR Parallel Processing for Deformation Time-Series Computation[J]. IEEE Journal of Selected Topics in Applied Earth Observations & Remote Sensing, 2014, 7(8):3285-3296.

[6] MARSELLA M. Analysis of ground deformation using SBAS-DInSAR technique applied to COSMO-SkyMed images, the test case of Roma urban area[C] SAR Image Analysis, Modeling, and Techniques XII. SAR Image Analysis, Modeling, and Techniques XII, 2012:448-453.

[7] 尹振興,鐘麗云,許兵,等.基于SBAS-InSAR的昆明地面沉降監測研究[J].地礦測繪,2016,(04):1-5.

[8] FERRETTI A, PRATI C, ROCCA F. Permanent scatterers in SAR interferometry[J]. IEEE Transactions on Geoscience and Remote Sensing, 2001, 39(1):8-20.

[9] HOOPER A, SEGALL P, ZEBKER H. Persistent scatterer interferometric synthetic aperture radar for crustal deformation analysis, with application to Volcán Alcedo, Galápagos[J]. Journal of Geophysical Research Atmospheres, 2007, 112(B7):B07407.

[10] 李濤. 基于點面散射體的多時相雷達干涉模型與形變探測方法[D].成都：西南交通大學,2014.

[11] BERARDINO P, FORNARO G, LANARI R,et al.A new Algorithm for Surface Deformation Monitoring based on Small Baseline Differential SAR Interferograms[C]. IEEE Transactions on Geoscience and Remote Sensing 40(11), 2375-2383 (2002).

[12] HOOPER A. A multi‐temporal InSAR method incorporating both persistent scatterer and small baseline approaches[J]. Geophysical Research Letters, 2008, 35(16):96-106.

[13] 廖明生，王騰. 時間序列InSAR技術與應用[M].北京：科學出版社，2014.

[14] 周志偉,鄢子平,劉蘇,等.永久散射體與短基線雷達干涉測量在城市地表形變中的應用[J].武漢大學學報(信息科學版),2011(8):928-931.

[15] 秦曉瓊,楊夢詩,廖明生,等.應用PSInSAR技術分析上海道路網沉降時空特性[J].武漢大學學報(信息科學版),2017(2):170-177.

[16] 劉欣,商安榮,賈勇帥,等. PS-InSAR和SBAS-InSAR在城市地表沉降監測中的應用對比[J]. 全球定位系統,2016(2):101-105.