基于Sentinel-1影像的2019～2020年河北隆堯地裂縫及周邊形變監測與分析

馮曉松 楊成生 趙超英 牛玉芬 高 涵

1 長安大學地質工程與測繪學院，西安市雁塔路126號，710054 2 河北工程大學礦業與測繪工程學院，河北省邯鄲市太極路19號，056038 3 邯鄲市自然資源空間信息重點實驗室，河北省邯鄲市太極路19號，056038 4 云南省地震局信息中心，昆明市茨壩藍桉路2號，650201

地裂縫是一種由內部和外部地質動力過程所引起的淺表破裂[1]，通常表現為土壤地表層在各種人為因素(采礦、澆灌、抽水等)[2-4]及其他自然因素(構造斷裂[4]、地震[5])影響下，土層的連續性被力的相互作用打斷，并以一種重大地裂縫和地表裂縫破壞現象呈現出來[6]。地裂縫直接或間接破壞了人類和野生動物的生存環境[7]，并且地裂縫周邊地區的地表往往會出現沉降現象[8]。

隆堯地裂縫位于華北平原隆堯縣內，形成于1966年邢臺地震[9]，震后裂縫進入非活躍期，2003年裂縫再次暴露[10]，2006年以來裂縫規模及活躍度迅速增加[11]?，F今，隆堯地裂縫是華北平原規模最大、災害程度最為嚴重的構造地裂縫[12]，分布在近EW向的20多個村莊里[13]。地表裂縫最大寬度為20 cm，可見深度達數m[14]。Yang等[15]研究了2007～2011年間隆堯地裂縫的形變特征，并進行建模處理，但隨著《全國地面沉降防治規劃(2011～2020年)》、《京津冀平原地面沉降綜合防治總體規劃》等文件的發布及地下水控采方案和南水北調工程的實施，隆堯地裂縫的特征可能發生了改變，但目前缺乏對其現階段活動性的研究。

有學者使用全球導航衛星系統GNSS、全站儀等手段對地裂縫進行了大量的形變監測研究，但傳統測量技術成本昂貴且數據獲取不連續，嚴重阻礙了對地裂縫的進一步研究[16]。干涉合成孔徑雷達InSAR技術是一項基于地球物理學、信號學、計算機等學科發展而來的新型技術，被廣泛應用于火山監測、地震監測、滑坡監測、地裂縫監測等地質災害監測領域。本文采用時序InSAR技術探測隆堯地裂縫及其周邊地區的形變，并進一步分析現階段隆堯地裂縫的活動特征。

1 研究區域及實驗數據

1.1 研究區概況

中國大陸位于印度、亞歐和太平洋等板塊大量匯聚形成的沖擊區和擠壓區，華北平原長期遭受拉張應力[17]，因此1 000多條裂縫遍布華北平原[18]。隆堯地裂縫位于華北平原的邢臺市，處于華北地塊的核心位置，擁有3個地質塊：臨清地塹、邢-衡隆起、寧晉-衡水斷層，形成了一個“三叉戟”的形狀。該區域分布有NE向和EW向斷層，其中主要裂隙是太行山山前斷裂和巨鹿EW向斷裂。

隆堯地裂縫是華北平原最長、最活躍的地裂縫，位于隆堯縣城與魏家鎮之間(圖1 紅色實線)，全長約36 km，地表裂縫最大寬度為20 cm，可見深度達數m[19]，沿途的道路、河流、農田、房屋等都被其切斷，造成了巨大的經濟損失。在隆堯地裂縫周邊分布有太行山山前斷裂、臨城斷裂、巨鹿西部斷裂、巨鹿東部斷裂、雞澤斷裂和明華斷裂，這些斷裂與隆堯地裂縫相互作用，使得隆堯地裂縫周邊地區的地表沉降尤為復雜。由于研究區域較為廣闊，本文采用InSAR技術對河北隆堯地裂縫及其周邊的地表形變進行監測。

圖1 隆堯地裂縫及其周圍地區地形Fig.1 Topographic of Longyao ground fissure and its surrounding area

1.2 研究數據

哨兵一號Sentinel-1采用的是C波段，其中心頻率為5.405 GHz，C波段在X波段和S波段之間，兼具兩者的特性。本文采用31景VV極化的Sentinel-1數據，所選時段為2019-01-03～2020-12-23，如表1所示。數據覆蓋范圍如圖1所示，圖中藍色實線代表裁剪后的雷達影像范圍，即本文的研究區域。本文使用的外部DEM數據是由美國地質勘探局USGS提供的SRTM DEM，分辨率為30 m。

2 研究方法

2.1 SBAS-InSAR

由于研究區域植被覆蓋度較高，失相干現象較為嚴重，為了通過不同時間基線干涉圖提供的多余觀測來降低失相干的影響，并去除干涉圖中的各種誤差(如大氣誤差、DEM殘差等)，本文首先采用SBAS-InSAR技術獲得隆堯地裂縫及周邊地區的時間形變序列。SBAS方法是在利用小基線干涉圖集進行形變測量的前提下，合成符合條件的小基線干涉圖對，然后基于形變速率的最小范數準則，通過SVD方法獲得系數矩陣的偽逆，最終獲得相干目標的形變速率及時間序列[20]。計算公式如下：

(1)

表1 實驗數據

2.2 Stacking InSAR

為獲取研究區域含有較小誤差的形變速率，使用Stacking InSAR技術對解纏后的差分干涉圖進行處理。Stacking InSAR技術假設大氣在時間上的擾動是隨機的，且地表形變為線性形變，因此可以將所有符合運算條件、經過解纏的干涉圖進行線性疊加，從而減少信號穿過對流層、電離層所產生的大氣誤差對監測結果的影響，提高監測精度[21]。公式如下：

(2)

式中，Vi為第i幅干涉圖的相位變化速率，φi為第i幅干涉圖的解纏相位，ΔTi為第i幅干涉圖的成像時間間隔。

3 數據處理及結果分析

3.1 數據處理

數據處理流程如下：1)首先對原始數據進行預處理，轉化為SLC格式，然后進行基線估計，空間基線設置為200 m，時間基線設置為100 d，時空基線如圖2所示；2)將得到的干涉對進行共軛相乘，得到差分干涉圖，并對干涉圖進行濾波以提高影像的信噪比；3)通過最小費用流法對濾波后的干涉圖進行相位解纏，得到地表像元的真實相位；4)分別通過SBAS-InSAR技術和Stacking InSAR技術對研究區域的數據進行時序處理。經過上述數據處理得到研究區域的形變時間序列和形變速率，其中形變速率如圖3所示，2019-01-03～2020-12-23的累積形變量如圖4所示。

圖2 時空基線Fig.2 Spatial-temporal baseline

圖3 2019～2020年隆堯地裂縫周邊地區形變速率Fig.3 Deformation rate of Longyao ground fissure from 2019 to 2020

圖4 2019～2020年隆堯地裂縫周邊地區累積沉降量Fig.4 Cumulative subsidence around Longyao groundfissure from 2019 to 2020

3.2 隆堯地裂縫及其周邊地區形變特征分析

3.2.1 隆堯地裂縫

由圖3和4可見，隆堯地裂縫的南北兩側有一定的形變速率差異，為了更加直觀形象地觀察隆堯地裂縫南北兩側的速率差，對隆堯地裂縫南北兩側進行剖線分析。為與Yang等[15]2007～2011年的ALOS PALSAR結果進行比較，本文選擇相近位置的剖線(圖3黑色虛線EF)。隆堯地區形變速率剖線所得結果如圖5所示，圖中淺灰色實線為隆堯地裂縫，實線左側為地裂縫北側、右側為地裂縫南側。由圖可見，相對于隆堯地裂縫南側，隆堯地裂縫北側呈抬升趨勢，二者的形變速率差達4 cm/a。

圖5 沿剖線EF隆堯地裂縫南北兩側形變速率差異Fig.5 Deformation rate of both sides of Longyao ground fissure along profile EF

將本文采用Sentinel-1衛星影像獲取的2019～2020年隆堯地裂縫南北側形變速率剖面，與采用日本地球觀測衛星合成孔徑雷達ALOS PALSAR-1獲取的2007～2011年形變速率剖面[14]進行對比分析發現，前者的隆堯地裂縫南側沉降速率比后者的相對沉降速率約低1 cm/a。搜集資料后發現，國務院在2012年發布了《全國地面沉降防治規劃(2011～2020年)》，隨后自然資源部又組織編制了《京津冀平原地面沉降綜合防治總體規劃》。隨著地下水控采措施和南水北調工程的實施，京津冀地區的沉降狀況在2016年以后得到了大幅度的緩解，沉降范圍和沉降中心形變量級的降低可能影響了隆堯地裂縫的活動性。此外，由剖線圖可以明顯看出，隆堯斷層兩側形變特征存在明顯差異，這一現象與Yang等[15]的研究中L波段獲取的形變特征類似。

3.2.2 地裂縫周邊地區

由圖3和4可見，除巨鹿縣縣城所在位置有輕微的抬升外，境內其他地區均處于沉降狀態。其中，固城店鎮、魏家鎮、官莊鎮是明顯的沉降中心，3個區域年沉降速率均達6～8 cm/a。對固城店鎮、魏家莊鎮、官莊鎮的形變時間序列進行統計分析(圖6)可知，固城店鎮、魏家莊鎮及官莊鎮處于周期沉降狀態，即當年2～9月處于快速沉降狀態，當年9月至次年2月又存在一定程度的抬升，但其總體形變趨勢仍是沉降。

圖6 固城店鎮、魏家莊鎮及官莊鎮單點形變時序Fig.6 Time series subsidence of Guchengdian, Weijiazhuang and Guanzhuang

調查發現，固城店鎮、魏家莊鎮、官莊鎮都種植了大面積的農作物，包括冬小麥、玉米、棉花、蔬菜等，每年需要抽取大量的地下水進行灌溉。由于河北獨特的氣候條件，當地的主要農作物(小麥、玉米)為一年兩收模式，冬小麥通常在當年10月播種并進行初次灌溉，次年開春時進行第2次灌溉。由此推測，當地對地下水的大量抽取及周邊斷層的影響造成了2～10月的沉降狀態，之后隨著灌溉抽水用量的急劇減少及冬季降雪和周圍地下水回流，區域表現為緩慢的地表抬升[22]。此外，隆堯地裂縫東南端有一個較大的形變區域(圖3和4中區域D)，查閱資料發現，當地于2003-11成立了隆堯工業園，主要生產方便面和飲料，其工業過程消耗大量的地下水，可能是造成當地地面沉降的主要原因。

4 隆堯地裂縫模型構建

研究表明，彈性均勻空間中的蠕變位錯[6,23-24]可用來解釋活動斷層的大地觀測。一個活動斷層的地表變形可以被建模為反正切，為斷層滑動速率和鎖定深度提供了估計方法。垂直走滑斷層運動的一個簡單模型是均質彈性半空間的，可以描述為:

(3)

式中，V(x)為沿斷層垂直剖面估計的各點平行斷層速度，x為離斷層的距離，S為斷層滑動速率，D為位錯深度，α為靜態偏移量。

為更好地了解隆堯地裂縫的地面變形，使用上述模型對InSAR形變剖面進行建模。通過計算得知，地裂縫深度為5 m，幾乎破裂至地表，滑移速率為27 mm/a，相比于2007～2011年的48 mm/a，滑移速率減緩約21 mm/a，說明現階段隆堯地裂縫雖仍處于活躍狀態，但隨著對京津冀地區地下用水的宏觀調控，地裂縫活躍度降低。

5 結 語

本文采用SBAS-InSAR技術，基于31景Sentinel-1衛星影像對邢臺市隆堯地裂縫及其周邊區域進行時間序列分析和物理模型構建。結果表明，隨著對京津冀地區地下水使用的宏觀調控，隆堯地裂縫及其周邊地區形變趨勢有所減緩，活躍度降低。具體結論如下：

1)現階段隆堯地裂縫兩側相對形變有所降低，相較于2007～2011年降低約1 cm/a；

2)隆堯地裂縫周邊地區如固城店鎮、魏家莊鎮、官莊鎮是明顯的沉降中心，3個區域年沉降速率均達6～8 cm/a，且呈現周期性的形變特征；

3)物理模型構建結果表明，隆堯地裂縫幾乎破裂至地表，滑移速率約為27 mm/a，相較于2007～2011年減緩約21 mm/a。