運城市地表沉降的時序InSAR監測與分析

牛寶勝，楊德宏，成飛飛，高霞霞，石建楊

(昆明理工大學國土資源工程學院，云南 昆明 650093)

1 引 言

隨著地下水開采及大型建筑物的密集建造，城市地面沉降問題日益突出，尤其是在經濟發達的城市，可能導致嚴重的財產損失和基礎下部建筑的塌陷。在我國北方多數城市和沿海許多地區，都面臨著地面沉降的問題[1，2]。地表沉降是一種由人為或自然原因引起的一種具有緩變性、不均勻及不可逆性的地質災害[3]，由地表形變引發的地質災害已經成為我國一大熱點問題。傳統的地面沉降監測手段主要依賴于精密水準測量和全球定位系統(GPS)等技術[4]，雖然測量精度高，但是水準測量存在作業周期長、耗費大、監測區域小、復雜地形人力無法到達等缺點；全球定位系統測量在高程精度方面，不能滿足大范圍形變監測的需要[5～7]，難以滿足城市地表形變精細化防控監測的要求。

合成孔徑雷達干涉測量技術(D-InSAR)是近二十年發展起來的極具潛力的微波遙感新技術，與常規大地測量形變監測技術相比，InSAR具有全天時、全天候、覆蓋范圍廣、成本低、安全和觀測連續等特點，能夠獲取高分辨率、高精度的地形及地表形變信息[8，9]。為克服傳統D-InSAR的時空失相關、大氣延遲和軌道誤差等影響[10]，由Ferretti[11]等2001年提出的PS-InSAR技術和Berandino[12]等2002年提出的SBAS-InSAR技術等新技術應運而生，被廣泛應用到有關形變監測的多個領域，包括由地下水開采和建筑施工等引起的城市地面沉降等方面。

地面沉降嚴重威脅著城市的發展，但目前針對運城市的地面沉降研究較少，難以對地面沉降及其時空演化特征進行詳細表征。因此，為了完善和豐富運城市地面沉降的時空變化監測，本文從歐空局獲得2019-1～2021-3時間段內的30景Sentel-1A升軌數據，綜合利用PS-InSAR與SBAS-InSAR兩種方法對運城市的沉降進行監測并交叉驗證，分析運城地面沉降現狀和空間分布特點，為城市地面沉降防控及科學治理提供重要的依據。

2 研究區概況及數據來源

2.1 研究區概況

運城市(圖1底圖使用空間分辨率為30 m的SRTM-1DEM數據)位于山西省南端，與陜西、河南兩省隔黃河而相望，北與臨汾市毗連，東與晉城接壤。總面積 14 182.78 km2，介于110°15′～112°04′E，34°35′～35°49′N之間。運城全年受季風活動影響，屬暖溫帶大陸性季風氣候，全市氣溫平均14.3℃，平均年降雨量 525.2 mm，區域內地形復雜，相對高差明顯，具有平原、山地、丘陵、盆地、臺地等多種地貌類型，地處山西斷陷盆地中地裂縫分布范圍最廣、發育規模最大的區域，地面沉降的形成和發展與區內數條隱伏斷裂有關[13]，且本地水資源匱乏，年際變化大，且分布不均勻，常年農田灌溉、居民用水等一直靠抽取地下水作為供水水源，已造成地下水位大幅下降，引起地面沉降等引發了一系列地質災害問題。

圖1 研究區域位置示意圖

2.2 數據源

為了研究該區域的地表沉降，主要選取了歐空局(ESA)于2014年發射的Sentine-1A衛星的30景干涉寬幅(interferometric wide swath，IW)升軌SAR影像，極化方式為VV，重訪周期 12 d，影像入射角為39.02°，時間跨度為2019年1月～2021年3月，共計769天，影像具體參數如表1所示。同時采用定位精度優于 5 m的AUX_POEORB精密軌道星歷參數用于主從影像間的空間基線估計，選取SRTM-1的 30 m分辨率的DEM用于地形相位去除與地理編碼。

研究區Sentine-1A影像參數 表1

3 數據處理方法

3.1 PS-InSAR技術

PS-InSAR技術是利用多景(最少25景)同一地區的SAR影像，通過統計分析時間序列上幅度和相位信息的穩定性，探測不受時間、空間基線去相關影響的穩定點(PS)目標。將得到的PS點進行時序建模分析及線性參數求解，通過時間、空間濾波，提取PS點的形變量，進一步對各PS點在各段時間內的實際總變形量進行空間規則格網插值，即得到研究區內在各時間段內的形變場，具體流程如圖2所示。

圖2 PS數據處理流程

考慮到時間基線、空間基線、多普勒中心的因素影響[14]，選取2020年3月8日成像的SAR影像為單一主影像，其余29景為輔助影像，共形成了29幅干涉像對，設置基線閾值500%，干涉對的時間、空間基線信息如圖3所示，其中黃點表示主影像位置。從圖3可以看出，空間基線主要分布在 -135 m～130 m之間，最長空間基線為 130.141 m，最短空間基線為 7.609 m。

圖3 PS時空基線分布

在PS點的選取上，采用振幅離差指數的方法對研究區內的穩定散射目標點進行識別，相干系數閾值設置為0.75，振幅離差閾值為0.25，最后探測出PS點 647 067個。其表達式為：

(1)

將獲取的PS點利用Delaunay不規則三角網建立PS候選點之間的連接關系，使用帶權的最小二乘方法進行相位解纏得到每個點的線性形變速率和DEM的絕對誤差值，通過時間域和空間域濾波分離殘余相位的非線性形變和大氣相位，將線性形變和非線性形變疊加即可得到每個PS點上的形變時間序列。實驗獲得的形變結果如圖4(a)，其中不同的顏色代表不同的沉降速率，負值表示沉降，正值表示抬升。在ArcGIS軟件上采用克里金(Kriging)插值獲取運城市區的地表沉降分布，為了保證插值的準確性，選取了干涉點較多、沉降速率較快的城區進行插值[15]，插值結果如圖4(b)，不同的顏色代表不同的沉降區域，顏色越深，沉降速率越大。

圖4 PS監測地面沉降結果

從圖4(a)、(b)可以看出，運城市南部區域有明顯的沉降，其他地區形變在空間上分布不均勻，且發生的沉降在東西方向上呈帶狀分布，形成了永濟-絳縣、稷山-新絳縣的帶狀沉降區域。2019-1到2021-2監測時間段內的沉降速率在 -28.50 mm/y～20.47 mm/y之間，最大沉降速率為 28.50 mm/y。運城市中心城區沉降明顯，出現了多個沉降漏斗，運城關公機場周圍及空港開發區、運城學院等地沉降最為嚴重，沉降速率介于 -28.50 mm/y～ -7.96 mm/y之間，平均沉降速率為 23.15 mm/y。稷山縣城區地表穩定，沒有沉降發生，但其西社鎮北部出現了以薛家莊村、麻古垛村、劉家莊村為中心的沉降漏斗，最大沉降速率為 27.50 mm/y；新絳縣城區域整體穩定，轄區的北張鎮出現大幅度沉降，平均沉降速率為 19.10 mm/y；絳縣城區出現小幅沉降，轄區柳泉村-橫水鎮區域呈帶狀沉降，平均沉降速率為 6.44 mm/y；夏縣整體呈現不均勻沉降，沉降較快的區域包括尉郭鄉、禹王鄉、牛家凹村、井溝村、西陰村等。

3.2 SBAS-InSAR技術

SBAS技術是對相干目標(coherent target，CT)進行相位分析來獲取時序形變，通過選擇合適的空間基線和時間基線閾值組成差分干涉對，利用具有較短時-空基線的影像對產生干涉圖提高相干性，選取相干目標點利用線性相位變化模型進行建模和解算，并通過時空濾波去除大氣延遲，在減少DInSAR處理中的失相關影響及高程、大氣誤差的同時獲取地表的形變時間序列[12]，主要技術流程圖如圖5所示。

圖5 SBAS數據處理流程

同樣以20200308影像為主影像，為了保證有足夠的干涉影像對，且限制時間和空間失相干的影響，選擇臨界基線閾值為45%和時間臨界基線閾值為365天進行短基線集干涉影像對組合。干涉處理時在方位向和距離向對影像進行了 4∶1的多視處理以削弱相位噪聲的影響。相干系數的閾值設為0.2，采用Goldstein方法進行濾波，最終生成的干涉基線分布如圖6所示，共生成284個干涉像對，最大臨界基線 257.362 m，最小臨界基線 0.030 m。

圖6 SBAS時空基線分布

相位解纏過程中，由于用于探測的緩慢失相關濾波相位像素點(slowly-decorrelating filtered phase，SDFP)點[16]是離散的，在相位解纏前，將滿足相干系數閾值條件的像元作為節點生成不規則三角網(Delaunay)，鑒于研究區域存在植被和水體，基于這些離散點的格網利用最小費用網絡流法(MCF)進行相位解纏。相位解纏完成后，將PS-InSAR獲取的GCP文件轉換到SAR坐標下，用于精確估算軌道參數，去除殘余相位和軌道誤差，即軌道精煉和重去平。將前面得到的相干點建立線性模型，采用奇異值分解法(SVD)估算形變速率和高程系數。

圖7 從殘余相位中分出大氣相位

圖8 剩余相位分離出非線性相位

圖9 基于SBAS-InSAR方法地面沉降監測結果

4 結果對比及分析

4.1 交叉驗證

由于缺少實測同期水準監測數據，本文通過比較PS(圖7)和SBAS(圖10)得到形變速率結果進行精度驗證。不難看出，兩種結果所呈現的沉降范圍及量級具有較好的一致性，主要沉降區分布在運城市區、稷山縣、夏縣、絳縣區域。區別在于SBAS獲得點的密度大于PS的結果，個別區域的沉降速率存在較大偏差。

兩種方法反演得到的地面沉降速率統計指標如表2所示，其中，由PS與SBAS兩種方法獲得的地面沉降速率平均值分別為 -3.02 mm每年(mm/y)和 -0.54 mm/y，相差 2.08 mm，兩種方法的反演結果很接近，且SBAS方法得到的相干點數是PS的15倍，結合圖7與圖9，可以看出SBAS一定程度上提高了相干點在非城區的密度。PS與SBAS方法反演地面沉降的速率分別是 -28.5 mm/y、-88.99 mm/y，分別位于稷山縣的劉家莊村和馬首官莊村。SBAS方法的標準差為5.93，稍小于PS方法的標準差6.70，表明SBAS方法反演的沉降速率分布更集中。

兩種方法沉降速率監測結果對比 表2

上面的分析可以看出，兩種方法得到的部分點的沉降速率有較大偏差，可能受個別因素的影響(如房屋倒塌等)，因此對其沉降速率分布進行統計，如圖10所示，PS和SBAS方法得到的相干點速率在 -10 mm/y～0 mm/y之間的點各占相干點總數的43.12%、40.18%，說明運城市地表沉降相對穩定，這與運城市對地下水的壓采有關，表明治理成果顯著。沉降速率大于 10 mm/y分別占16.49%、5.10%，表明小部分區域發生較大沉降，且已經形成沉降漏斗。整體來看，出現沉降的區域分別占整個運城市的59.61%、49.09%，結合圖7與圖10，沉降區域主要分布在運城市區、夏縣、稷山縣等地。

圖10 地面沉降速率數值分布

在沉降形成的漏斗區域，本文選取6對典型點進行分析，如圖11所示，在運城市區選取了運城學院、軍屯村、圣惠北路附近三個點，在絳縣選取了橫水鎮的典型點，在夏縣選取了馬村的典型點，在新絳縣選取了西行莊村的典型點，繪制了在2019年～2021年期間的時序沉降折線圖。從圖中可以看出，在研究期間，兩種數據處理方法得到的6個點的沉降趨勢具有較高的一致性，圖11(a)、(e)中PS和SDFP(SBAS)點的沉降量基本保持一致，運城學院在2020年6月18號以后SDFP點沉降加劇，在2021年后PS點的沉降超過SDFP點。圖11(b)、(c)、(d)、(f)中PS和SDFP點的沉降量保持一致，受到了共模系統誤差的影響，SDFP點沉降量整個大于PS點。總體來看，它們的沉降趨勢一直呈持續加劇趨勢。

圖11 特征點沉降序列

為進一步確定監測結果的精度，對PS和SBAS進行內符合精度驗證，在運城市城區選取兩種方法中經緯度相同的點735個同名點，比較兩種方法的年均沉降速率。以PS點沉降速率為橫軸，SBAS點沉降速率為縱軸，得到的散點圖如圖12所示，由圖可知，兩者的相關系數R2達到了 0.887 2，進一步驗證了PS和SBAS的沉降速率具有較高的一致性。

圖12 PS與SBAS年均沉降速率關系圖

4.2 原因分析

(1)構造活動：研究區域地處山西地塹系西南端的一個半封閉斷陷盆地(運城盆地)，盆地南、東、北三面環山，是一個強烈的沉降盆地，內部形成多個斷裂帶(如圖13所示)，結合PS與SBAS的沉降速率圖，可以發現稷山縣與新絳縣的交界帶附、夏縣的尉郭鄉至鹽湖區陶村鎮等地的沉降是由斷裂活動引起的，說明斷裂帶的出現往往伴隨著地面沉降，兩者相輔相成。

圖13 運城市地裂縫分布

(2)地下水開采：在運城盆地內，多數地面沉降與地下水的超采有關[17]，地下水的超采導致地下水水位下降，從而引起可壓縮土層的壓縮，進而引起地面沉降。運城屬于嚴重缺水地區，中深層地下水是運城市工農業生產和城鎮生活用水的主要水源，開采量占總用水量的75%～80%。從圖14可以看出，自20世紀80年代以來，運城盆地地下水開始大量開采，持續的開采地下水致使水位不斷下降[18]，產生了沉降漏斗，2002年漏斗中心水位降到了最低點，2007年完成關井壓采121眼后，2008年水位有所回升。截至2015年底，地下水超采仍占全市總面積的28%[20]。受區域超采的影響，漏斗區域水位仍呈下降趨勢，史雙雙[19]監測2012年運城市最大沉降速率為 50 mm/y～60 mm/y，位于夏縣禹王鄉、尉郭鄉及裴介鎮中間地帶，從2019年～2021年的監測結果來看，禹王鄉、尉郭鄉及裴介鎮中間地帶的最大沉降速率為 42.25 mm/y，相較2012年速率減緩，這與運城市實施地下水壓采措施有關。但農業灌溉用水仍以大水漫灌為主，根據2020年的《運城市統計年鑒》可知，較上一年用水量增加7.3%，其中2019年底農業用水13億m3，部分農村地區仍采用地下水澆灌農田，對地面沉降起到了加速作用。可見超量開采是運城城區地面沉降形成的最主要原因。

圖14 地下水開采量與降落漏斗中心水位變化曲線

(3)地震與密集高程建筑的影響：運城地震活動比較頻繁，自公元1485年到現在，有記錄的Ms5級以上的破壞地震中，5.0級～5.9級10次、6.0級～6.9級2次，對運城地面沉降造成一定程度的影響。同時伴隨著人口的快速增長和經濟的發展，大規模城市工程建設的增多，也對地面沉降產生一定促進作用，截至2019年底，運城市的房屋建筑面積為 1 980.8萬m2，比上年增長18.1%，密集的高層建筑在大規模城市改造建設中地面沉降效應明顯，在研究區域中，圖15中姚孟街道周邊的沉降等就是新建樓房引起的，但其沉降量達到一定程度就會趨于穩定。

圖15 姚孟街道周邊(光學影像)

5 結 論

本文利用PS-InSAR和SBAS-InSAR技術處理了運城市30景哨兵數據，獲取了運城市2019-1～2021-2的地面形變信息，由于沒有獲取到監測時間段內的水準數據，通過對比分析與已有資料來驗證成果的可靠性。結果表明，運城市城區、夏縣、稷山、絳縣、新絳縣的部分地區沉降較為嚴重，PS的最大沉降速率達到了 28.50 mm/y，位于運城市市區，SBAS的最大沉降速率達到了 84.44 mm/y，其位于稷山縣。結合已有資料，地面沉降的主要原因在于地下水的超采，部分是城市工程建設增加等因素的影響，對運城市沉降防治提供了一定的參考。試驗區周圍有大量的農田，除冬季外，可能會產生嚴重的失相干，需要后續實驗進一步討論與分析。