基于SBAS-InSAR 技術的陽城地區采空區地表沉陷分析

2023-03-11 08:25:02孫志杰郝曉陽

地理空間信息 2023年2期

熊威，孫志杰，郝曉陽

（1.山西省交通科技研發有限公司黃土地區公路建設與養護技術交通行業重點實驗室，山西太原 030032；2.山西省煤炭地質物探測繪院資源環境與災害監測山西省重點實驗室，山西晉中 030053）

山西省作為傳統的煤炭資源大省，為國家的發展建設提供了源源不斷的能源供應，同時也形成了近5 000 km2的采空區，其中沉陷區約3 000 km2，給采空區周邊居民生活、公路運營、工程建設等造成了較大的安全隱患。為了對采空區的影響范圍和發育程度進行有效評估，有必要采取一種高效便捷的觀測手段對采空區進行識別和監測。傳統地表沉陷監測方法主要有GNSS、水準儀、全站儀等[1-3]，它們均是以點為單元進行測量，效率低下，無法覆蓋全區，需要耗費大量的人力財力。合成孔徑雷達干涉測量（interferometric synthetic aperture radar，InSAR）技術作為一項新型空間對地觀測技術，具有大面積、全天候、連續、高精度、高空間分辨率等優點，由此衍生而來的差分干涉測量（D-InSAR）技術在地震同震形變[4]、城市地面沉降[5]、礦區地表沉陷[6]等方面做了廣泛的研究，取得了一定的成效。由于D-InSAR 容易受到大氣效應、DEM誤差、基線誤差等因素的影響，其中，大氣效應對雷達信號傳播產生的延遲，在中緯度地區一年內的影響區間是20 cm，是最主要的誤差源[7]。為了消減各項誤差對監測結果的影響，在此基礎上發展了多時相InSAR 技術（MT-InSAR）來進行高精度形變監測，永久散射體（PS-InSAR）和小基線集（SBAS-In-SAR）是其中最具代表性且廣泛使用的兩種方法，二者均具有毫米級精度，且能回溯地表歷史形變過程，在采空區沉陷監測與成因分析[8-9]和形變預測[10-11]方面發揮了重要作用。PS-InSAR 方法基于一景主影像與其余輔影像配準后差分，獲取PS點的位移時間序列，適用于相干性較強的城鎮地區，且對數據量要求較高，SBAS-InSAR技術克服了PS-InSAR因選取一幅影像作為公共主影像而引起的部分干涉圖相干性較差的不足，同時降低了對SAR數據的需求量，運算效率較高，適用性更廣。

本文針對山西省陽城地區煤礦采空區現狀，利用SBAS-InSAR技術對2019年1月至2020年7月間共25景Sentinel-1A影像進行處理，根據時間序列形變結果分析黃土地區采空區沉陷的動態變化過程，以及引起沉陷變化的直接因素，為黃土地區采空區的特征識別及危險性分析提供技術支持。

1 SBAS-InSAR技術

當M≥N時，通過最小二乘法即可求解ν；當M＜N時，B為秩虧矩陣，可利用SVD方法對矩陣B進行分解，進而求得ν的最小范數二乘解。最后對ν在時間維度上積分即可得到相應時刻的累積形變量。

2 測區概況與數據處理

2.1 測區概況

測區位于山西陽城縣北部，晉運高速、濟陽高速在這里交匯（圖1），其中在晉運高速陽城段兩側分布有大大小小幾十個采煤區，開采時間從20 世紀80 年代至今，高速南側大部分開采時間早于2010年，北側大多為近年來形成的采空區，部分礦區甚至下穿高速。該地區主要為黃土丘陵地貌，地層從上往下主要為第四系全新統人工堆積層（Q4 ml）素填土，全新統沖洪積層（Q4 al+pl）粉質黏土、碎石土，第四系上更新統坡洪積層（Q3 dl+pl）黏質黃土和二疊系上統（P2）泥巖、砂巖。

圖1 測區地理位置圖

2.2 數據處理

本文利用2019-01-08—2020-07-01 的歐洲空間局Sentinel-1A 衛星數據對采空區的地表形變進行監測，數據詳情如表1 所示，并配合使用與影像數據相匹配的精密軌道數據，其位置精度優于5 cm。外部DEM 使用的是美國宇航局提供的SRTM-1，其分辨率為30 m，高程精度約為10 m。上述數據的使用能最大限度的減少DEM誤差和軌道誤差對干涉圖及解纏結果的影響。

表1 影像信息表

本文利用GMTSAR 軟件中SBAS 模塊對溝東隧道上覆邊坡的SAR數據進行處理，該軟件具有開源、自動處理、易于繪圖等優點[13]，在震后大范圍地表形變監測方面取得了較好的效果[14]。由于測區位于丘陵地帶，季節性的灌木叢生，為了減少時空失相干的影響，增加相干點數，設置時間基線閾值為180 d，空間基線閾值為最大基線長度的30%，組合形成的影像連接圖如圖2 所示，共形成163 對有效干涉像對，每幅影像至少與3 幅以上的影像產生干涉，網型健壯性較好。相干性閾值設置為0.1。在經過影像裁剪、配準、干涉圖生成、去地形相位、相位解纏、形變反演、地理編碼等一系列步驟之后最終得到累積形變結果如圖3 所示。其中距離向與方位向多視比為5∶1，成圖分辨率為20 m；解纏方法為經典的最小費用流法，解纏分解等級設置為1；最終選取的有效控制點數量為30個，均遠離采空區且均勻分布于全圖；初始反演模型選擇較為穩定的線性模型；軌道精煉選擇3 次多項式擬合方法。

圖2 影像連接圖

圖3 累積形變量圖

3 結果驗證與精度評價

為了驗證SBAS-InSAR技術監測結果的精度，利用前期布設在6#采空區周邊的水準線路進行對比分析。該條水準線路沿東西向的鄉村公路布設，包含1 個基準點KZ1 和12 個監測點（JC1～JC12），基準點位于曲堤村，點位穩定。線路總長3.3 km，平均點位間距0.3 km，按照二等精密水準的精度要求測設，監測頻率1次∕月，均在每月上旬完成。由于InSAR技術得到的是衛星視線向上的形變，水準監測得到的是垂直方向形變結果，兩者之間有如下轉換關系：

式中，θ為衛星視線向（LOS）與地表法向的夾角；dH為精密水準得到的地表垂向監測結果。將水準監測結果轉換到衛星視線向，并根據基準點沉降位移值做校正，再與InSAR結果做對比（圖4），并根據公式（7）計算兩者之間的相關性。

圖4 水準與時序InSAR數據之間相關性

通過式（7）計算結果可以看出兩類監測數據的線性相關程度較高，相關系數R2為0.982，表明二者具有較強的一致性。為進一步驗證時序InSAR數據的可靠性，將2 種手段得到的監測點位累積形變量進行對比分析，結果如圖5所示。

圖5 水準測量結果與時序InSAR形變結果對比

從圖5 中可以看出水準測量結果與時序InSAR 形變結果在不同點位的最大形變差值為+9.8 mm，最小形變差值為-0.5 mm，平均差值為3.0 mm，中誤差為6.7 mm。造成二者之間差異的主要原因有：①InSAR監測結果的像元分辨率為20 m×20 m，無法在空間位置上與水準點保持一致；②2 種監測方式獲取的監測結果時間上不能完全同步，有數天的差異；③時序In-SAR結果仍殘存有大氣、DEM、軌道誤差。從整體趨勢上來看，二者變化基本保持一致，時序InSAR技術監測精度小于1 cm，具有較高的可靠性。

4 結果分析

4.1 時序結果分析

從圖3可以看出晉陽高速、安陽高速兩旁分布有大大小小十幾個沉陷區，均位于采空區正上方，其中面積較大、形變量級＞50 mm、距離公路2.5 km以內的沉陷區有6個，1號沉陷區邊緣已經侵入公路界線，須謹防擴大。從2019年1月至2020年7月間最大累積形變量達到14.8 cm，位于6#采空區范圍內。選取10個典型沉陷區的最大形變點做進一步的時序分析（圖6），發現1、2、10號點形變有逐漸變緩的趨勢，其他7個沉陷區仍處于勻速或加速下沉的趨勢中。對沉陷點正上方的采空區開采時間進行統計分析，結果如表2所示。

圖6 沉降時序曲線

表2 采空區情況統計表

上述煤礦的開采方式均采用傾斜長壁采煤法，傾斜角小于10°，回采方式為后退式回采。從表2中可以看出1、2、10#煤礦開采時間較長，初采時間均在2010年前，2019年底局部區域開采完畢，地表沉陷在半年內趨于穩定。3#煤礦初采時間在2008年，采深較大，監測點位于2019年新掘采空區上方，形變趨勢較為平緩。4～8#采空區的初采時間為近4 a，形變速率和開采深度正相關，采深越大，形變速率越小，且目前仍在開采中，形變量持續增加，部分點還有加速下沉跡象。9#煤礦初采時間為2012 年，監測點位于2020年待采區域內，可以看出從2020年6月開始有加速下沉跡象。

4.2 沉陷區動態變化分析

為了深入分析采煤沉陷區的動態變化過程，以6#、8#沉陷區作為研究對象，繪制不同方向的多個剖面進行分析，其中剖面AA’與6#采空區多個采煤斷面掘進方向垂直，剖面BB’與6#采空區掘進方向平行，剖面DD’與8#采空區的掘進方向平行。

圖7 6#采空區典型剖面累積沉降變化

6#采空區包括從東至西多組不同時期的子采空區組成，每個子采空區的掘進方向均為由南至北，采用斜井開拓，走向長壁分層綜采，支護方式為支撐掩護式綜采支架，采深為320～470 m。為了細致地研究不同時期地面沉陷動態變化過程，判斷開采時間和掘進方向與地面沉陷之間的關系，選取3 條剖面上間隔一個季度的多期累積沉降量曲線進行分析。

從剖面AA’的沉降曲線上可以看出有3 個沉降中心，分別對應2019 年、2018 年和2017 年3 個子采空區橫斷面，開采時間越晚的采空區累積沉降量越大。2016 年采空斷面由于時間超過2 a，沉降趨勢基本穩定。剖面BB’從南向北累積沉降量越來越大，對應的開采時間由2018 年年初至年末，初采時間越晚，被記錄到的殘余形變越大。剖面CC’橫跨2 個子采空區，沉降曲線上相應有2 個沉降中心，且從2019 年4—7 月第二個沉降中心開始形成，與該區域的開采時間基本一致，在該時間段內第一個沉降中心也有加速趨勢，可能為第一個工作面上繼續進行分層開采造成的。從以上分析可以看出InSAR結果的沉降剖面上能較好地識別出不同時期采空區的沉降中心，且累積沉降量的變化趨勢能清晰反映開采方向、開采時間對地面沉降造成的影響。利用DD’剖面進一步分析開采方向上不同時間節點不同點位的沉降速率變化。8#采空區包括從南至北多組不同時期的子采空區，其中DD’剖面位于2019—2020 年采空區上，開采方式采用斜井開拓，長壁開采，采深350～500 m。沿掘進方向依次選取E、F、G 3 個監測點位進行沉降速率分析，E點從2019年2月開始沉降速率由小于10 mm∕a急速變大，F 點從2019 年4 月份開始沉降速率迅速增加，G 點的沉降速率峰值則出現在2019年7 月份，最大沉降速率分別達到68 mm∕a、79 mm∕a和78 mm∕a。沉降速率超過50 mm∕a 的時間周期為2～3個月，其他時間段內沉降速率會有小幅波動，可能為上覆巖層穩定性發生改變或是周邊煤礦采掘造成的影響。由于煤礦開采造成的地表沉陷響應較快，不同點位的沉降速率峰值基本能夠代表掘進位置變化。