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

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

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

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

1 PS-InSAR方法

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

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

(1)

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

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

(2)

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

(3)

2 SBAS-InSAR方法

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

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

(4)

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

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

(5)

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

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

(6)

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

(7)

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

Bv=δφ.

(8)

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

3 實例分析

3.1 研究區(qū)概況及實驗數(shù)據(jù)

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

圖1 研究區(qū)地理位置

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

3.2 數(shù)據(jù)處理

3.2.1 PS干涉處理

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

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

3.2.2 小基線干涉處理

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

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

4 結(jié)果分析

4.1 沉降結(jié)果分析

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

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

圖4 研究區(qū)年平均沉降速率圖

圖5 A區(qū)域年平均沉降速率圖

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

圖6 B區(qū)域年平均沉降速率

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

圖7 C區(qū)域年平均沉降速率

4.2 對比分析

4.2.1 整體沉降情況分析

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

圖8 年平均沉降速率分布

4.2.2 相關(guān)性分析

圖9 PS與SBAS年平均沉降速率相關(guān)性圖

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

4.2.3 剖面圖分析

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

4.2.4 時序分析

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

綜上所述，證明PS-InSAR和SBAS-InSAR兩種技術(shù)監(jiān)測結(jié)果的準(zhǔn)確性、可靠性和一致性。

圖10 A區(qū)域沉降漏斗剖面

圖11 B區(qū)域沉降漏斗剖面

圖12 A區(qū)域最大沉降漏斗累積沉降量

圖13 B區(qū)域最大沉降漏斗累積沉降量

5 結(jié)束語

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

參考文獻(xiàn)：

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

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

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

[4] 高照忠,羅志清,魏海霞. 3S技術(shù)在昆明市地面沉降監(jiān)測中的應(yīng)用[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的昆明地面沉降監(jiān)測研究[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á)干涉模型與形變探測方法[D].成都：西南交通大學(xué),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技術(shù)與應(yīng)用[M].北京：科學(xué)出版社，2014.

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

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

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