基于SBAS-InSAR技術的新疆某煤礦長時序地表形變監測與分析

摘" "要：新疆某煤礦長期開采導致礦區出現嚴重地表形變，嚴重威脅礦區的安全生產和區域可持續發展，有必要采用有效方法監測煤礦地表形變。為深入研究礦區地表形變規律，選取2017年1月23日至2022年12月23日共178景Sentinel-1A升軌雷達影像數據，利用SBAS-InSAR技術實測礦區2017—2022年地表形變信息，并引入GACOS數據去除大氣誤差影響，分析了礦區地表形變空間分布特征、時空演化規律及地表形變驅動因素。研究結果表明：①監測期間，研究礦區內地表形變以沉降為主，地表最大形變速率約為-113 mm/a，最大累計沉降量約為580 mm，主要地表沉降區域面積約為15.2 km2；②研究礦區地表沉降面積和累計沉降量不斷增大，地表沉降符合煤礦開采沉降曲線，地表沉降將會進一步增加；③研究礦區地表形變主要由采礦工程所致，降水也是影響該礦區地表變形的重要因素。本研究為煤礦區地表形變監測、地質災害預警與防控及礦區安全生產提供科學依據。

關鍵詞：SBAS-InSAR；GACOS；地表形變監測；地表沉降；時空演化規律；驅動因素

煤炭是我國能源安全的保障和壓艙石，煤炭開采在促進我國經濟發展的同時也引發了一系列地質災害問題[1]。隨著開采規模不斷擴大，礦區出現嚴重的地面塌陷問題，嚴重破壞了礦區建筑、交通道路、管線等基礎設施，加劇水土流失和土地沙漠化，對礦區安全生產和區域可持續發展造成嚴重威脅[2]。因此，采取高效準確的方法對礦區開展地表形變監測，研究地表形變特征和時空演化規律，并分析其影響因素，對礦區災害預警及治理具有重要意義。

地表形變監測的傳統方法主要有精密水準測量、全站儀測量及全球導航衛星系統測量等[3]。這些方法具有成熟可靠，監測精度高等優點，在監測礦區地表形變中發揮了重要作用[4]。但仍存在諸多不足，如監測周期長、成本高、監測點易丟失、時空分辨率低且無法進行大范圍沉降監測[5，6]。近年來，合成孔徑雷達干涉測量（InSAR）技術迅速發展并逐步應用于礦區地表形變的長時序監測，與傳統方法相比，該方法具有全天候、覆蓋廣、分辨率高、成本低等優點[7]。差分雷達干涉測量（D-InSAR）技術是利用同一地區的兩幅不同時間的SAR影像進行干涉處理，從而獲取高精度微小變形，監測精度理論上可達毫米級，但D-InSAR技術受時空失相干、大氣延遲及軌道誤差等干擾較大[8，9]。為解決D-InSAR技術缺陷，諸多學者提出了時序InSAR（TS-InSAR）技術，極大地提高了利用InSAR技術監測地表微小形變的能力，其主要包括永久散射體干涉測量（PS-InSAR）技術、小基線集干涉測量（SBAS-InSAR）技術和分布散射體干涉測量（DS-InSAR）技術等[10，11]。Ferretti等首次提出了PS-InSAR技術[12，13]，時間序列上的幅度和相位信息穩定性選取SAR圖像上的永久散射體目標，并利用形變模型反演地面形變速率。DS-InSAR技術是近年來研究的熱點，它利用分布散射體獲取地表形變信息，現在仍處于探索階段[14，15]。Berardino等首次提出了SBAS-InSAR技術[10，16]，該技術主要以多主影像的干涉對為基礎，利用高相干點獲取地表時間序列形變信息。

目前，國內外學者利用SBAS-InSAR技術對礦區地表形變監測與分析開展大量研究，并對算法進行相應的改進。杜建濤等利用SBAS-InSAR對蔚縣煤礦區進行地表形變監測，獲取礦區在2017—2018年間的地表形變空間分布特征，并對礦區地表沉降量級及面積進行統計分析[17]。刁鑫鵬等基于SBAS-InSAR技術對河北某礦區進行監測，反演了研究區2019年07月14日—2020年04月27日期間地表形變的空間分布形態和時間演變過程，并利用水準數據進行驗證，最后結合區域地質資料分析地表沉降原因[18]。Wang等以SAR圖像間的平均相干性作為基線約束指標[19]，對SBAS-InSAR干涉對進行優化選擇，并利用優化方法對江西萍鄉煤礦區的地表形變進行監測和分析。雖然SBAS-InSAR在礦區地表形變監測中得到了較好結果，但結果仍受大氣延遲影響，針對這一問題，Yu等提出使用一種通用大氣校正模型對InSAR干涉圖進行校正[20]。Liu等利用SBAS-InSAR技術對三山島金礦開采沉陷進行分析和研究[21]，并利用通用大氣校正在線服務（GACOS）對大氣相位進行校正，降低數據誤差對研究結果的影響。

目前，已有部分學者利用SBAS-InSAR技術對礦區地表形變開展研究，但針對新疆煤礦區的研究相對較少，且部分研究結果受大氣延遲的干擾，精度和可靠度較差。因此，本文以新疆某煤礦為研究對象，選取178景Sentinel-1A升軌雷達影像數據，利用SBAS-InSAR技術獲取礦區2017—2022年的地表形變信息，并結合GACOS數據對干涉圖進行大氣校正，進而分析地表形變的空間分布特征及時空演化規律，并進一步研究地表形變的相關驅動因素，科學指導礦區未來的安全生產與發展。

1" 研究區概況與數據源

研究區為新疆某煤礦，位于阿克蘇地區庫車縣北部，占地面積約為33 km2。研究區所處地勢較平坦，平均海拔1 914 m，地表植被覆蓋率較低，受時間失相干影響較小。

研究數據源主要有Sentinel-1A影像數據、AW3D30數字高程（DEM）數據、GACOS數據和降水數據等。其中，選取了2017年1月23日—2022年12月23日共178景Sentinel-1A升軌影像數據用于提取研究礦區地表形變信息。利用空間分辨率為30 m的AW3D30 DEM數據模擬地形相位，去除地形相位和平地效應的影響。此外，采用GACOS數據對所有差分干涉圖進行大氣相位校正，削弱大氣延遲對監測結果的影響。

2" 研究方法

2.1" "SBAS-InSAR技術基本原理

2002年Berardino等首次提出SBAS-InSAR技術，其主要原理是：N+1幅覆蓋同一區域的SAR影像組成若干子集合，設置適當的時空極限閾值形成不同的干涉對，從而生成差分干涉圖[16， 22]。假設不考慮高程、大氣及噪聲等誤差因素，干涉圖中像元干涉相位為：

[δφk（i，j）=φ（tA，i，j）-φ（tB，i，j）≈4πλ[d（tA，i，j）-d（tB，i，j）]]" " （1）

式（1）中：[δφk（k=1，…，M）]為干涉相位值；[φ（tA，i，j）]和[φ（tB，i，j）]分別表示像元[tA]和[tB]時刻的相位；[d（tA，i，j）]和[d（tB，i，j）]分別表示像元在[tA]和[tB]時刻沿雷達視線向的形變。

計算得到的差分干涉圖，利用奇異值分解法將若干子集合聯合求解獲取研究區的形變速率，進而求解出研究區的形變量。

2.2" 基于GACOS的大氣校正技術基本原理

本文利用GACOS生成的大氣校正圖對SBAS-InSAR中的差分干涉圖進行大氣校正，提高監測的精確度和可靠度。GACOS的基本原理是利用迭代對流層分解（ITD）模型生成90 m分辨率的大氣校正圖，利用ITD模型可從對流層天頂總延遲（ZTD）中分離出分層延遲和湍流延遲，ZTD被定義為[23，24]：

[ZTDk=T（xk）+L0e-βhk+εk] （2）

式（2）中：[ZTDk]為[k]處的天頂總延遲；[xk]為站點坐標向量；[T（xk）]為[xk]處的湍流分量；[L0]為海平面的分層分量延遲；[e-β]為分層分量； [hk]為高度標尺；[εk]為剩余未建模殘余誤差。

2.3" 數據處理流程

本文利用SBAS-InSAR技術對所獲取的178景Sentinel-1A影像數據進行時間序列處理，求解出形變速率，進而求解出研究區整個觀測時間段的形變量，SBAS-InSAR的技術主要包括以下5個步驟（圖1）：①基線估算和生成連接圖。綜合考慮研究區的地形、氣候、時間基線和空間基線等因素，確定時間基線閾值為72天，空間基線閾值為臨界基線的2%，形成708個干涉對，其中2020年2月1日獲取的影像為超級主影像，時空基線見圖2；②數據干涉處理。首先對主輔影像進行多級精確配準，其次利用復數像對共軛相乘得到干涉圖，并利用DEM數據去除地形相位。采用Goldstein濾波方法進行濾波處理，去除干涉結果中的噪聲。采用Delaunay最小費用流算法（Delaunay MCF）對濾波后的干涉圖解纏，獲得解纏相位圖，解纏相干系數閾值設置為0.2；③軌道精煉與重去平。解纏后仍存在殘余恒定相位和相位坡道，在相位相對穩定且無明顯形變的區域選擇適量的地面控制點（GCP），然后利用多項式擬合優化方法去除殘余相位的影響；④SBAS反演。第一次反演：利用線性模型估算形變速率和殘余地形相位；第二次反演：在第一次反演基礎上，利用大氣濾波，去除大氣相位與殘余相位；⑤地理編碼。反演得到的形變信息由SAR坐標系轉為地理坐標系，并投影到雷達視線向和垂直向。

3" 結果與分析

3.1" 地表沉降空間分布特征分析

利用SBAS-InSAR技術對研究礦區2017年1月23日—2022年12月23日的Sentinel-1A升軌影像數據進行處理，得到研究礦區長時序地表形變速率圖（圖3），其中，負值表示地表沉降，正值表示地表抬升，藍色越深表示地表沉降越大。InSAR監測結果表明，2017年1月—2022年12月研究礦區范圍內多處發生地表沉降，地表沉降分布范圍較廣，主要集中在A1～A2和B1～B4重點監測點位置附近，沉降區域呈條帶狀分布。研究礦區地表最大形變速率約為-113 mm/a，最大累計沉降量約為580 mm，主要地表沉降區域面積約為15.2 km2。從表1可看出，地表形變速率最大區域位于A1處，該區域由于地表形變速率過大，出現較為嚴重的失相干問題。

3.2" 地表沉降時空演化特征分析

3.2.1" 地表沉降空間演化特征分析

為研究礦區地表沉降空間演化特征，以2017年6月4日為起始時間，每6個月時間為一個間隔，共選取12幅研究礦區的累計形變量圖（圖4）。從圖4可看出，2017年6月—2018年6月礦區地表沉降面積增加明顯，B3重點監測點附近區域沉降明顯，最大累計沉降量約為165 mm，主要沉降面積約7 km2；2018年12月—2019年12月礦區地表沉降面積持續增加，A1重點監測點附近沉降面積增加明顯，B1和B3重點監測點附近局部區域（藍色區域）沉降量變化迅速，累計沉降量增加明顯，最大累計沉降量約338 mm，主要沉降面積約為11.01 km2；2020年6月—2021年6月礦區整體地表沉降面積增加放緩，A1、B1和B3重點監測點附近局部區域沉降繼續增大（藍色加深），B4區域沉降面積增加明顯，最大累計沉降約為463 mm，主要沉降面積約13.7 km2；2021年12月—2022年12月，礦區沉降面積總體保持穩定，A1、B1和B3重點監測點附近沉降繼續增加，累計沉降量約580 mm，主要沉降面積約15.2 km2。在整個監測期間，研究礦區逐漸形成2個沿煤礦工作面走向的主要沉降漏斗，隨時間推移，礦區范圍內沉降漏斗未來可能逐步連片形成一個較大范圍的沉降區。

3.2.2" 地表沉降時間演化特征分析

為研究礦區地表沉降隨時間的演化，選取礦區6個重點監測點（圖3），提取監測點2017—2022年長時序形變信息，獲取重點監測點累計形變量隨時間變化折線圖（圖5）。從圖5可看出，各個重點監測點均呈現出不同速率的非線性地表沉降，各重點監測點累計沉降量曲線基本符合煤礦開采沉陷曲線，總體來看，地表累計沉降量呈逐年增加的趨勢。其中，A1點在監測時間段內地表累計沉降量最大，約為568 mm，該區域在整個監測期間開采時間較長，開采強度較高，地表形變最為嚴重；B4點在監測時間段內地表形變累計沉降量最小，約253 mm，該區域開采時間較短，開采強度較低，因而地表形變較為緩慢。

3.2.3" 地表沉降未來趨勢

運用線性擬合方法對重點監測點累計形變量進行擬合（圖6）。由圖6分析可知，2017—2022年地表累計形變與時間基本呈線性關系，隨時間增加，地表累計沉降逐漸增大。據此推測礦區地表沉降會進一步增加，且在未來很長時間內不會停止，需引起相關部門的高度重視，做好防治工作，防止嚴重地質災害的發生。

3.3" SBAS-InSAR結果精度驗證

由于缺乏水準測量數據和GNSS等精密外部測量數據，本文利用內符合精度和實地調查來對比驗證本次監測結果精度和可靠性[25，26]。結果表明，不同時段SAR影像干涉測量獲得的形變速率均方差多在3 mm/a內，內符合精度滿足《地質災害InSAR監測技術指南》（T/CAGHP 013—2018）要求。同時，據現場實地調查發現，礦區東北部出現因地表沉降造成的地表拉張裂縫和地面塌陷（圖7），證實了通過SBAS-InSAR技術監測地表沉降的可行性。

3.4" 地表形變驅動因素分析

3.4.1" 采礦對地表形變的影響

煤礦區地表形變多因地下煤礦開采導致，從圖4可看出，隨采礦的進行，地表沉降面積不斷擴大，地表累計沉降量不斷增大，地表沉降區域基本成條帶狀分布，符合現代煤礦以綜采工作面為主的開采方式。

為進一步分析采礦對地表形變的影響，選取B研究區為研究對象，分別沿1-1′和2-2′繪制地表形變剖面，分析其在時間序列上的沉降變化情況，提取地表形變剖線上各時期累計沉降量繪制曲線圖（圖8）。從沿1-1′和2-2′剖線各期地表形變圖可看出，2017年6月—2018年6月煤礦開采工作面由中北部向西不斷推進，B研究區西部地表沉降逐漸增大；2018年12月礦區通過產能核增，礦區整體開采規模和開采強度增大，因此，2018年12月—2021年6月研究礦區整體地表沉降迅速增大；2021年6月—2021年12月研究礦區以西部工作面開采為主，地表累計沉降繼續增加；2021年12月—2022年12月研究礦區中西部開采基本結束，考慮對采空區進行放頂處理，因此中西部地表沉降保持穩定。綜上，2017—2022年采礦工程是引起并加快研究礦區大面積地表沉降的主要因素。

3.4.2" 降水對地表形變的影響

大氣降水同樣是影響地表形變的重要因素。為研究地表形變與降水量之間的關系，本文選取了礦區所在地區2017—2020年月降水量數據與A1和B1監測點地表沉降速率進行對比分析（圖9）。據圖9分析發現，年降水量較多月份重點監測速率點地表沉降速率呈上升趨勢，附近區域地表沉降速率加快。但圖9中A1監測點在2018年5月地表沉降速率降低，考慮該區為煤矸石堆積，地表累計沉降量減小，導致地表沉降速率降低。綜上可知，降水也是影響地表形變的重要因素之一。

4" 結論

（1） 整個監測期間，研究礦區內地表形變以沉降為主，地表沉降分布較廣，局部沉降明顯，地表最大形變速率約為-113 mm/a，最大累計沉降量約為580 mm，主要地表沉降區域面積約為15.2 km2。

（2） 隨煤礦不斷開采，研究礦區地表沉降的面積和程度不斷增大，礦區沿煤礦開采工作面形成2個沉降漏斗，并在未來可能連片形成一個整體沉降區。

（3） 研究礦區各重點監測點呈現不同程度的非線性變形，累計沉降量曲線基本符合煤礦開采沉陷曲線，推測礦區地表沉降在未來會進一步增加，相關部門需高度重視，做好地表沉降的防治工作，防止嚴重地質災害的發生。

（4） 據研究礦區地表沉降驅動因素發現，采礦工程是導致地表沉降的主要因素，降水是影響礦區地表形變的重要因素。

參考文獻

[1] 袁亮.我國煤炭主體能源安全高質量發展的理論技術思考[J].中國科學院院刊，2023，38（1）：11-22.

[2] 仝云霄，楊俊泉，王雪，等.基于時序InSAR的山西大同煤田地表沉降監測及時空演化分析[J].中國地質，2024.51（1）：170-183.

[3] 劉輝，陳斯滌，朱曉峻，等.基于D-InSAR技術的煤礦工業廣場動態沉降特征研究[J].煤田地質與勘探，2023，51（5）：99-112.

[4] 周大偉，安士凱，吳侃，等.礦山開采損害InSAR/UAV融合監測關鍵技術及應用[J].煤炭科學技術，2022，50（10）：121-134.

[5] 張豐，刁鑫鵬，譚秀全，等.基于SBAS技術的煤礦采空區探查與穩定性評價[J].煤炭科學技術，2022，50（3）：208-214.

[6] 柴華彬，胡吉彪，耿思佳.融合實測數據的地表沉降SBAS-InSAR監測方法[J].煤炭學報，2021，46（S1）：17-24.

[7] 許強，董秀軍，李為樂.基于天-空-地一體化的重大地質災害隱患早期識別與監測預警[J].武漢大學學報（信息科學版），2019，44（7）：957-966.

[8] DU Y，ZHANG L，FENG G，et al.On the Accuracy of Topographic Residuals Retrieved by MTInSAR[J].IEEE transactions on geoscience and remote sensing，2017，55（2）：1053-1065.

[9] 李達，鄧喀中，高曉雄，等.基于SBAS-InSAR的礦區地表沉降監測與分析[J].武漢大學學報（信息科學版），2018，43（10）：1531-1537.

[10] 朱建軍，李志偉，胡俊.InSAR變形監測方法與研究進展[J].測繪學報，2017，46（10）：1717-1733.

[11] 明小勇，田義超，張強，等.基于Sentinel-1A欽防地區地面沉降監測與分析[J].自然資源遙感，2024（1）：35-48.

[12] FERRETTI A，PRATI C，ROCCA F.Nonlinear subsidence rate estimation using permanent scatterers in differential SAR interferometry[J].IEEE transactions on geoscience and remote sensing，2000，38（5）：2202-2212.

[13] A.F，C.P，F.R.Permanent scatterers in SAR interferometry[J].IEEE Transactions on Geoscience and Remote Sensing，2001，39（1）：8-20.

[14] FERRETTI A，FUMAGALLI A，NOVALI F，et al.A New Algorithm for Processing Interferometric Data-Stacks：SqueeSAR[J].IEEE transactions on geoscience and remote sensing，2011，49（9）：3460-3470.

[15] 范銳彥，焦健，高勝，等.InSAR時序分析高相干目標選取方法比較研究[J].地球信息科學學報，2016，18（6）：805-814.

[16] BERARDINO P，FORNARO G，LANARI R，et al.A new algorithm for surface deformation monitoring based on small baseline differential SAR interferograms[J].IEEE transactions on geoscience and remote sensing，2002，40（11）：2375-2383.

[17] 杜建濤，閆麗，趙超英.蔚縣礦區地面沉陷InSAR多維形變監測[J].煤田地質與勘探，2020，48（1）：168-173.

[18] 刁鑫鵬，孫全帥，白志輝，等.采動區地表異常形變SBAS識別與損害原因揭示[J].煤炭科學技術，2022，50（11）：168-175.

[19] WANG S，ZHANG G，CHEN Z，et al.Surface deformation extraction from small baseline subset synthetic aperture radar interferometry （SBAS-InSAR） using coherence-optimized baseline combinations[J].GIScience and remote sensing，2022，59（1）：295-309.

[20] YU C，LI Z，PENNA N T，et al.Generic Atmospheric Correction Model for Interferometric Synthetic Aperture Radar Observations[J].Journal of Geophysical Research：Solid Earth，2018，123（10）：9202-9222.

[21] LIU J，MA F，LI G，et al.Evolution Assessment of Mining Subsidence Characteristics Using SBAS and PS Interferometry in Sanshandao Gold Mine，China[J].Remote sensing （Basel，Switzerland），2022，14（2）：290.

[22] 張香凝，賀黎明，劉翠芝，等.基于SBAS-InSAR技術的煤礦開采沉陷監測與分析[J].遙感技術與應用，2022，37（4）：1021-1028.

[23] YU C，LI Z，PENNA N T.Interferometric synthetic aperture radar atmospheric correction using a GPS-based iterative tropospheric decomposition model[J].Remote Sensing of Environment，2018，204：109-121.

[24] YU C，PENNA N T，LI Z.Generation of real-time mode high-resolution water vapor fields from GPS observations[J].Journal of Geophysical Research：Atmospheres，2017，122（3）：2008-2025.

[25] 劉琦，岳國森，丁孝兵，等.佛山地鐵沿線時序InSAR形變時空特征分析[J].武漢大學學報（信息科學版），2019，44（7）：1099-1106.

[26] 姚鑫，張路青，李凌婧，等.工程場址區域地殼穩定性In SAR評價研究[J].工程地質學報，2021，29（1）：104-115.

Long-Term Surface Deformation Monitoring and Analysis of A Coal

Mine in Xinjiang Based on SBAS-InSAR Technology

Zhang Xuehui1， Cui Zhendong2，3，4， Zhang Zhongjian1， Zhao Leilei5， Wei Tao2，3，4， Liu Dongxu5， Wang Longcan5

（1.School of Engineering and Technology，China University of Geosciences， Beijing，100083，China;2.Key Laboratory of Shale Gas and Geoengineering，Institute of Geology and Geophysics， Chinese Academy of Sciences， Beijing，100029，China;3.Institutions of Earth Science，Chinese Academy of Sciences， Beijing，100049，China;4.College of Earth

and Planetary Sciences， University of Chinese Academy of Sciences，Beijing，100049，China;5.School of Mechanics

and Civil Engineering，China University of Mining and Technology-Beijing， Beijing，100083，China）

Abstract： The long-term mining of a coal mine in Xinjiang causes serious surface deformation in the mining area， which poses a serious threat to the mining area's ability to produce safely and to the area's ability to develop sustainably.Effective procedures must be adopted to keep track of the mining area's surface deformation.This work selected 178 Sentinel-1A ascending radar image data from January 23， 2017， to December 23， 2022， to further examine the surface deformation law of the mining area.The mining area's surface deformation data from 2017 to 2022 was measured using SBAS-InSAR technology.To eliminate the impact of atmospheric inaccuracy， GACOS data was added.Analysis was done on the surface deformation in the mining region， including its driving elements， driving features， and spatial-temporal evolution law. According to the findings： （1） Surface subsidence is the predominant cause of surface deformation in the studied area during the monitoring period.The greatest rate of surface deformation is around -113 mm/a， the maximum cumulative subsidence is approximately 580mm， and the primary surface subsidence area is approximately 15.2 km2; （2） The mining area's cumulative sinking and surface subsidence are both expanding.The surface subsidence follows the coal mining subsidence curve， and it will continue to rise in the future; （3） Mining engineering is primarily to blame for the surface deformation of the mining region， although precipitation is also a significant influence in this deformation.This study provides a scientific foundation for surface deformation monitoring， geological disaster warning， prevention， and control， and mine safety production in coal mining areas.

Key words： SBAS-InSAR; GACOS; Surface deformation monitoring; Surface subsidence; Spatial-temporal evolution law; Driving factors