InSAR技術(shù)在高寒礦區(qū)地災(zāi)治理監(jiān)測中的應(yīng)用

吳啟紅 謝飛鴻 楊何 董建輝 魏占璽 毋遠(yuǎn)召 黃程

摘 要：采用InSAR技術(shù)對大通煤礦沉陷區(qū)進(jìn)行監(jiān)測，得到了該礦區(qū)2014年10月至2021年11月的地面沉降變形情況和年均沉降速率等信息.結(jié)果表明，礦區(qū)主要沉降變形區(qū)集中在東北部和東南部2個片區(qū)，其沉降速率最大超過220 mm/yr；2017年后礦區(qū)沉降變形范圍和變形量值都呈顯著下降趨勢且沉降變形逐漸平穩(wěn)，前期治理工程對控制礦區(qū)沉降變形起到了重要作用，可為高寒礦區(qū)防災(zāi)減災(zāi)及治理工程評價提供重要參考.

關(guān)鍵詞：InSAR；礦區(qū)治理；地面沉降；沉降速率

中圖分類號：P694

文獻(xiàn)標(biāo)志碼：A

0 引 言

礦區(qū)地質(zhì)災(zāi)害主要包括崩塌、滑坡、泥石流、地面塌陷、地裂縫、礦井突水和矸石堆災(zāi)害等；礦區(qū)環(huán)境災(zāi)害主要包括水污染、土污染、水土流失、石漠沙化、地貌景觀破壞和植被破壞等.據(jù)不完全統(tǒng)計(jì)，在1970—2005年，各類礦山地質(zhì)災(zāi)害共計(jì)達(dá)到2 031次，而在2001—2010年僅10年時間就發(fā)生了礦山次生災(zāi)害12 367起，其中尤以地面塌陷、滑坡和崩塌居多［1］.礦區(qū)內(nèi)剝離大面積土壤植被及形成的煤矸石山，產(chǎn)生大量“黑色”污染和“灰色”污染，形成潛在不穩(wěn)定斜坡，對礦區(qū)的生產(chǎn)設(shè)施和施工人員造成嚴(yán)重威脅.而高寒礦區(qū)由于高海拔、低氧和低溫，其地質(zhì)環(huán)境尤為復(fù)雜及脆弱，如玉龍煤礦、木里煤礦和紅河溝煤礦等.因此，開展高寒礦區(qū)的變形監(jiān)測對保障礦區(qū)安全至關(guān)重要.

傳統(tǒng)的礦區(qū)沉降監(jiān)測方法主要采用精密水準(zhǔn)測量和全球定位系統(tǒng)（GPS）測量等.這些監(jiān)測方法存在監(jiān)測工作量大、費(fèi)時、費(fèi)財(cái)和測點(diǎn)難以保存等缺陷，同時其測量精度并不高，從而導(dǎo)致從其測量數(shù)據(jù)獲得的結(jié)果會出現(xiàn)一定的偏差，對后續(xù)工作造成影響［2］.合成孔徑雷達(dá)干涉測量（InSAR）因在地表變形監(jiān)測上具有全天時、全天候、高分辨率、高覆蓋率、直觀和形象等其他測量技術(shù)不具備的優(yōu)點(diǎn)，可以在部分領(lǐng)域代替?zhèn)鹘y(tǒng)的測量方式［3-4］.同時為去除大氣誤差的影響，加之合成孔徑雷達(dá)（SAR）影像數(shù)據(jù)的增加，在原有InSAR方法的基礎(chǔ)上，逐漸發(fā)展出Stacking-InSAR、D-InSAR、PS-InSAR和SBAS-InSAR等高級時序InSAR處理方法［5-6］.

InSAR技術(shù)最先被用于地球表面大范圍形變場的觀測，如地震位移測量和冰川漂移等.隨著理論發(fā)展與數(shù)據(jù)處理技術(shù)的進(jìn)步，InSAR技術(shù)逐漸在國內(nèi)外各行各業(yè)中得到廣泛應(yīng)用，如煤田、石油、天然氣開產(chǎn)礦區(qū)以及地下水開采區(qū)地面沉降監(jiān)測［7］.高海英等［8］通過小基線集InSAR技術(shù)識別黔西南州貞豐縣和安龍礦山地表形變；丁劉建等［9］基于SBAS-InSAR對滕州市附近礦區(qū)地面沉降進(jìn)行監(jiān)測；王志勇等［10］基于InSAR對濟(jì)寧礦區(qū)沉降進(jìn)行精細(xì)化監(jiān)測.這些研究主要針對礦區(qū)生產(chǎn)運(yùn)營階段進(jìn)行變形監(jiān)測，而對礦區(qū)治理過程及后期的變形監(jiān)測及應(yīng)用研究不足.InSAR技術(shù)無論是在礦區(qū)生產(chǎn)運(yùn)營階段、閉礦階段還是在治理過程中對災(zāi)害的預(yù)防及評價方面都具有重要意義，同時對于指導(dǎo)我國在類似地區(qū)的防災(zāi)減災(zāi)工作也具有重要參考價值.因此，為了研究InSAR技術(shù)在高寒礦區(qū)治理中應(yīng)用的有效性和實(shí)用性，本研究選取青海省大通煤礦區(qū)作為研究對象，采用InSAR技術(shù)對大通沉陷區(qū)進(jìn)行處理，得到了該礦區(qū)治理過程及治理后的地面沉降變形情況和年均沉降速率等信息，分析變形規(guī)律及評價治理工程的有效性.

1 研究區(qū)概況及數(shù)據(jù)

1.1 研究區(qū)位置

研究區(qū)位于青海省西寧市大通回族土族自治縣，地處青藏高原和黃土高原的過渡地帶，地理位置東經(jīng)101°37′～101°56′，北緯36°55′～36°58′之間，海拔2 450～2 850 m，距省會西寧市35 km.

1.2 治理工程情況

2012—2015年共分4期對沉陷區(qū)內(nèi)涉及橋頭鎮(zhèn)及良教鄉(xiāng)的7個村（橋頭鎮(zhèn)的元樹爾、小煤洞和大煤洞3個村，以及良教鄉(xiāng)的下甘溝、煤洞溝、白崖和上甘溝4個村）范圍內(nèi)因煤炭開采造成地表土地資源壓占、挖損、植被破壞、粉塵污染和地下含水層破壞等，以及后期剩余煤炭資源開采可能引發(fā)的采空塌陷范圍進(jìn)行治理，治理建設(shè)規(guī)模為1 385.02 hm2.雖然經(jīng)過4期地質(zhì)環(huán)境治理，但區(qū)內(nèi)依然存在地質(zhì)災(zāi)害和生態(tài)環(huán)境破壞等問題，因此，2020—2022年，又對歷史遺留和新產(chǎn)生的塌陷區(qū)進(jìn)行了綜合整治.

1.3 數(shù)據(jù)情況

1.3.1 Sentinel-1數(shù)據(jù)及獲取

Sentinel-1衛(wèi)星數(shù)據(jù)是歐洲航天局為接替Envisat-ASAR衛(wèi)星發(fā)射的科研SAR衛(wèi)星，數(shù)據(jù)免費(fèi)開放.該衛(wèi)星成像質(zhì)量高，軌道穩(wěn)定，有多極化能力，單顆衛(wèi)星重訪周期24 d，繼Sentinel-1A衛(wèi)星后，Sentinel-1B衛(wèi)星的發(fā)射使常用重訪周期變?yōu)?2 d，遇到突發(fā)情況時可以2顆衛(wèi)星合作將觀測重訪周期縮短至6 d.查詢可知，該地區(qū)Sentinel-1數(shù)據(jù)覆蓋良好，數(shù)據(jù)日期自2014年10月至今.Sentinel-1數(shù)據(jù)參數(shù)見表1.

Sentinel-1衛(wèi)星執(zhí)行全球?qū)Φ赜^測極化，以干涉寬幅成像為主，該模式有益于InSAR的干涉處理，且250 km×250 km的幅寬保證單景數(shù)據(jù)覆蓋觀測區(qū).該衛(wèi)星工作波段為C波段，雖然對植被的穿透能力較弱，但衛(wèi)星數(shù)據(jù)量多、成像質(zhì)量高，可以一定程度上彌補(bǔ)數(shù)據(jù)穿透能力弱的問題.本次計(jì)算使用2014年10月7日至2021年11月23日間共計(jì)173期降軌Sentinel-1數(shù)據(jù)，相對軌道號33，數(shù)據(jù)重訪周期為12/24 d.

1.3.2 DEM數(shù)據(jù)

本次計(jì)算使用ALOS World 3D DEM，其空間分辨率為30 m，該數(shù)據(jù)集是高精度全球數(shù)字地表模型數(shù)據(jù)，水平分辨率為 30 m，高程精度5 m，由高級陸地觀測衛(wèi)星ALOS上搭載的全色遙感立體測繪儀PRISM獲取.是目前世界上最精確的3D地圖之一，覆蓋全球所有的土地尺度.

2 InSAR數(shù)據(jù)處理

2.1 InSAR技術(shù)原理

InSAR技術(shù)主要基于雷達(dá)影像的相位信息獲取目標(biāo)點(diǎn)至雷達(dá)傳感器之間的距離，通過獲取2期甚至多期次影像并對其進(jìn)行干涉即可丈量目標(biāo)點(diǎn)在一定時間段內(nèi)的位移變化量，其屬于定量遙感技術(shù)范疇［11-12］.目前，常用星載SAR平臺獲取的SAR影像進(jìn)行差分干涉處理，其具有覆蓋面積大、飛行軌道穩(wěn)定和獲取數(shù)據(jù)成本較低等優(yōu)勢.星載SAR平臺的基本工作模式為重軌飛行，即1顆雷達(dá)衛(wèi)星在相鄰軌道上對同一地區(qū)2次SAR成像，每幅SAR圖像中的每個像素均記錄了強(qiáng)度值和絕對相位值.針對同一地區(qū)獲取的2幅SAR影像進(jìn)行差分干涉處理，可得到該地區(qū)地表形變信息.

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

2.2.1 SAR數(shù)據(jù)配準(zhǔn)

在進(jìn)行干涉計(jì)算前，首先要對生成的單視復(fù)數(shù)數(shù)據(jù)（SLC）在外部DEM數(shù)據(jù)的參與下進(jìn)行數(shù)據(jù)配準(zhǔn).Sentinel-1數(shù)據(jù)成像模式較為特殊，為多個子帶拼接構(gòu)成的，為保證子帶間干涉處理時不發(fā)生相位的躍變，SLC數(shù)據(jù)配準(zhǔn)精度要求較高，尤其是針對Sentinel-1數(shù)據(jù)而言，距離向（Azimuth）配準(zhǔn)精度需要達(dá)到1/1 000像元級別.

2.2.2 差分干涉處理

如圖1所示，選擇2018年1月19日SAR影像為主影像進(jìn)行圖像配準(zhǔn)計(jì)算.其中，Azimuth方向配準(zhǔn)精度小于1/1 000像元，滿足Sentinel-1數(shù)據(jù)配準(zhǔn)要求.干涉像對最小時間間隔為12 d，最大時間間隔為120 d，最小垂直基線長度為0.3 m，最大垂直基線長度為216.8 m，完全滿足干涉計(jì)算需求，共有干涉像對513個.

2.2.3 自適應(yīng)濾波處理

為平滑干涉圖噪聲，使用窗口為32的自適應(yīng)濾波進(jìn)行處理，增強(qiáng)干涉圖低質(zhì)量區(qū)的相干性，方便下一步相位解纏處理.

2.2.4 相位解纏

使用最小費(fèi)用流方法進(jìn)行相位解纏處理，相干性閾值設(shè)置為0.2，解纏相位均過度平滑，未出現(xiàn)相位跳躍的情況，但部分解纏相位圖受大氣效應(yīng)干擾嚴(yán)重.

2.2.5 大氣相位估計(jì)

估算高程相關(guān)的大氣延遲相位，并利用快速空間濾波的方式初步估算隨機(jī)大氣相位，而后去除大氣延遲相位，大氣相位干擾情況將被明顯改善.

2.2.6 Stacking-InSAR和SBAS-InSAR計(jì)算

利用去除大氣相位后的計(jì)算結(jié)果進(jìn)行2014—2021年度的Stacking-InSAR計(jì)算，估算年均變形速率，并分析工作區(qū)逐年沉降變形情況.同時通過SBAS-InSAR方法計(jì)算典型形變點(diǎn)的累計(jì)變形曲線.

2.3 誤差源分析

InSAR計(jì)算中，會受到多種因素的影響，而且在不同地區(qū)，使用不同類型數(shù)據(jù)，同樣的影響因素造成的誤差程度不同，在計(jì)算過程中要結(jié)合實(shí)際情況加以分析，針對不同地區(qū)和不同數(shù)據(jù)來調(diào)節(jié)參數(shù).誤差源主要有系統(tǒng)誤差、DEM誤差、失相干和大氣誤差.

3 InSAR結(jié)果分析

受制于Sentinel-1衛(wèi)星數(shù)據(jù)本身分辨率，工作區(qū)InSAR觀測結(jié)果空間分辨率為近15 m.工作區(qū)觀測到的變形主要是采礦后沉降變形，所以將InSAR計(jì)算得到的雷達(dá)視線向（LOS）變形轉(zhuǎn)為垂直向變形.

3.1 2014—2021年逐年變化分析

利用Stacking-InSAR方法計(jì)算工作區(qū)2014—2021年間逐年地面沉降變形情況，2014年的InSAR數(shù)據(jù)是從2014年10月開始，經(jīng)過運(yùn)算后，圖中代表的是整個自然年的變形分布，如圖2所示.

由圖2可知，大通礦區(qū)地面沉降情況自2014—2021年呈現(xiàn)在變形范圍和沉降量值上均表現(xiàn)為減弱的趨勢.礦區(qū)的沉降區(qū)主要集中在東北部和東南部，可以分為2個大的沉降區(qū)，其他零星分布.沉降時間主要發(fā)生在2016年前，在2014—2015年間年均沉降速率大，最大可達(dá)220 mm/yr以上，自2017年后大通礦區(qū)沉降變形呈現(xiàn)為明顯減弱趨勢，2018年后明顯沉降區(qū)主要集中在礦區(qū)東南角.

3.2 2014—2021年多年期時間序列分析

采用SBAS-InSAR方法針對大通礦區(qū)獲取其2014—2021年的形變情況，SBAS-InSAR年均形變速率結(jié)果如圖3所示.

由圖3可知，大通礦區(qū)主要沉降區(qū)可分為1#沉降區(qū)和2#沉降區(qū)，整個礦區(qū)面上為零星分布點(diǎn)狀變形，高變形區(qū)主要分布在礦區(qū)東南角，未出現(xiàn)大面積高速沉降區(qū)，礦區(qū)2014—2021年間年均變形速率約為58 mm/yr.選取礦區(qū)6處典型變形點(diǎn)做時間序列分析，其累計(jì)變形曲線如圖4所示.SBAS1和SBAS2位于礦區(qū)東北部的1#沉降區(qū).SBAS1位置累計(jì)變形量約340 mm，SBAS2位置累計(jì)變形量約320 mm，2條曲線都表明在2014—2017年初變形速率較大，后變形速率逐漸平穩(wěn).SBAS3和SBAS4位于礦區(qū)東南部的2#沉降區(qū).SBAS3位置累計(jì)變形量約270 mm，曲線表明2014—2017年初變形速率較大，后變形速率逐漸平穩(wěn)，到2020年變形速率又開始加大；SBAS4位置累計(jì)變形量已達(dá)400 mm，曲線表明該點(diǎn)2014—2017年初變形速率較大，后逐漸平穩(wěn)，至2021年4月后短暫出現(xiàn)加速變形情況，后逐漸平穩(wěn).SBAS5位置累積變形量約160 mm，曲線表明2014年6月至2017年6月間變形速率較大，后逐漸平穩(wěn)，并在2021年1月至2021年6月間出現(xiàn)加速變形情況.SBAS6位置累積變形量約40 mm，曲線表明該點(diǎn)2014—2019年變形較為平穩(wěn)，而2020年開始出現(xiàn)加速變形的情況，但形變量值較小.

在礦區(qū)內(nèi)做5條年均沉降速率剖面，剖面線位置如圖3所示，剖面結(jié)果如圖5所示.剖面A-A′表明沿該剖面線方向，地面沉降速率逐漸變大，在剖面線上1 200 m位置處出現(xiàn)轉(zhuǎn)折，變形速率開始減??；剖面B-B′表明沿該剖面方向變形速率較為平穩(wěn)，局部出現(xiàn)較大沉降變形區(qū)；剖面C-C′表明沿該剖面線方向約1 300 m前變形速率平穩(wěn)，1 300 m后出現(xiàn)較大速率變形區(qū)，在1 700m位置處變形速率最大；剖面D-D′表明沿該剖面線方向約1 000 m前變形速率平穩(wěn)，1 000 m后出現(xiàn)較大速率變形區(qū)；剖面E-E′表明沿該剖面線方向整體變形平穩(wěn)，局部單點(diǎn)存在較大變形，所以剖面呈離散狀.

3.3 2020年11月至2021年11月地面變形分析

2014—2021年間時間跨度較長，大通礦區(qū)地面變形在長時間跨度下可能會出現(xiàn)階段性的大幅度變化，所以使用2020年11月至2021年11月間33期影像對工作區(qū)做1年期的SBAS-InSAR沉降變形速率觀測，其結(jié)果如圖6所示.

由圖6可知，大通礦區(qū)2020年11月至2021年11月間最大年均變形速率不超過52 mm/yr，高變形區(qū)主要分布在礦區(qū)東南角，整個礦區(qū)面上為零星分布點(diǎn)狀變形，未出現(xiàn)大面積高速沉降區(qū).同樣選取礦區(qū)6處典型變形點(diǎn)做時間序列分析，其累計(jì)變形曲線如圖7所示.SBAS1與SBAS2位置累積變形量較小，數(shù)據(jù)較離散，變形較為平穩(wěn)，與圖4中趨勢一致；SBAS3和SBAS4位于礦區(qū)東南部，其變形曲線顯示該區(qū)域于2021年5月后開始加速變形，加速變形于2021年10月結(jié)束，后變形較為平穩(wěn)；SBAS5點(diǎn)位于礦區(qū)中部位置，其變形曲線表明該位置于2021年5月后開始加速變形，并于2021年7月后加速變形停止，后變形逐漸平穩(wěn)；SBAS6點(diǎn)位于礦區(qū)西部，其變形曲線表明該點(diǎn)一直處于勻速變形狀態(tài).

礦區(qū)內(nèi)沿坡向做5條年均變形速率剖面，剖面線位置如圖6所示，剖面線結(jié)果如圖8所示.與圖5年均速率變形分布圖一致，沿剖面線方向沉降變形表現(xiàn)為零星點(diǎn)狀分布，并未連接成為整體，所以剖面線方向上速率點(diǎn)較為離散.

3.4 對比分析

根據(jù)前期勘察設(shè)計(jì)資料和現(xiàn)場調(diào)查，2個集中沉降區(qū)主要位于煤礦采空區(qū)及影響區(qū)內(nèi).2012—2015年對礦山進(jìn)行4期次的地質(zhì)環(huán)境治理示范工程，主要治理措施包括煤矸石棄渣處置、廢棄礦井井口封堵、塌陷區(qū)土地整治、泵站及田間灌溉配套、溝道防洪、小煤洞渣堆景觀區(qū)、引水蓄水灌溉系統(tǒng)、排水系統(tǒng)、不穩(wěn)定斜坡治理、田間道路、林網(wǎng)及礦山公園等工程.經(jīng)過治理，沉降區(qū)范圍及沉降量逐漸減小，特別是2017年沉降基本趨于穩(wěn)定，可見治理工程對控制沉降起到了重要作用，治理效果明顯.2020年礦區(qū)內(nèi)部分沉降變形繼續(xù)增加，但變形速率較之前減小較多，減小約3/4.根據(jù)圖6和圖7變形特征圖，結(jié)合2021—2022年度綜合治理實(shí)施方案（見圖9）進(jìn)行對比，可見實(shí)施方案中12個不穩(wěn)定斜坡、3個滑坡和5個沉陷區(qū)在InSAR解譯圖中都有對應(yīng)的變形解譯特征點(diǎn).

4 結(jié) 論

1）利用SBAS-InSAR觀測大通礦區(qū)2014—2021年間長時間序列變化，表明礦區(qū)主要沉降變形區(qū)集中在東北部和東南部2個片區(qū)，大通礦區(qū)2014—2016年間沉降速率最大超過220 mm/yr，自2017年后礦區(qū)沉降變形范圍和變形量值都呈顯著下降趨勢且沉降變形逐漸平穩(wěn)，但是存在局部或某時間段內(nèi)發(fā)生突然變形加速的現(xiàn)象.

2）2014—2021年的5條沉降速率剖面表明主體沉降位于礦區(qū)中前部，后部（西部）沉降速率很小；2020年11月至2021年11月間的InSAR結(jié)果表明，高變形區(qū)主要分布在礦區(qū)東南角，整個礦區(qū)面上為零星分布點(diǎn)狀變形，未出現(xiàn)大面積高速沉降區(qū).

3）通過對比分析，表明前期治理工程對控制礦區(qū)沉降變形起到了重要作用，近期變形有所增加，但變形速率減小較大；可見InSAR技術(shù)對于隱患點(diǎn)的識別及變形監(jiān)測具有大范圍和高精度的特點(diǎn)，可以為防災(zāi)減災(zāi)及治理工程評價提供依據(jù).

參考文獻(xiàn)：

［1］陳家彪，沈冰，劉飛燕.西南地區(qū)礦山環(huán)境地質(zhì)問題研究［J］.礦產(chǎn)綜合利用，2007，146（4）：43-46.

［２］劉國林，張連蓬，成樞，等．合成孔徑雷達(dá)干涉測量與全球定位系統(tǒng)數(shù)據(jù)融合監(jiān)測礦區(qū)地表沉降的可行性分析［J］.測繪通報，2005，51（11）：13-16．

［3］蒲川豪，許強(qiáng)，趙寬耀，等.利用小基線集InSAR技術(shù)的延安新區(qū)地面抬升監(jiān)測與分析［J］.武漢大學(xué)學(xué)報（信息科學(xué)版），2021，46（7）：983-993.

［4］陸會燕，李為樂，許強(qiáng)，等.光學(xué)遙感與InSAR結(jié)合的金沙江白格滑坡上下游滑坡隱患早期識別［J］.武漢大學(xué)學(xué)報（信息科學(xué)版），2019，44（9）：1342-1354.

［5］李承航，張文春.基于InSAR技術(shù)的監(jiān)測應(yīng)用及研究進(jìn)展［J］.北方建筑，2020，5（4）：15-19.

［6］范洪冬，鄧喀中，承達(dá)瑜.InSAR拓展和融合技術(shù)在礦山開采監(jiān)測中的應(yīng)用［J］.金屬礦山，2008，382（4）：7-10.

［7］張訓(xùn)虎，章磊，郝樹賓，等.合成孔徑雷達(dá)干涉（InSAR）測量技術(shù)應(yīng)用及展望［J］.北京測繪，2014，115（2）：28-31.

［8］高海英，趙爭，章彭.時序InSAR的貴州地質(zhì)災(zāi)害監(jiān)測［J］.測繪科學(xué)，2020，45（7）：91-99.

［9］丁劉建，陶秋香，高騰飛，等.SBAS InSAR技術(shù)在礦區(qū)地面沉降監(jiān)測中的應(yīng)用［J］.中國科技論文，2019，14（3）：320-325.

［10］王志勇，張繼賢，黃國滿.基于InSAR的濟(jì)寧礦區(qū)沉降精細(xì)化監(jiān)測與分析［J］.中國礦業(yè)大學(xué)學(xué)報，2014，43（1）：169-174.

［11］Bamler R，Hartl P.Synthetic aperature radar interferometry［J］.Inverse Probl，1998，14（4）：1-54.

［12］Rosen P A，Hensley S，Joughin I R，et al.Synthetic aperture radar interferometry［J］.P IEEE，2000，88（3）：333-382.

（實(shí)習(xí)編輯：姚運(yùn)秀）

Abstract：

InSAR technology is used to monitor the subsidence area of Datong coal mine，and the information of land subsidence deformation and average annual subsidence rate in the mining area from October 2014 to November 2021 are obtained.The results show that the main subsidence deformation areas in the mining area are concentrated in the northeast and southeast，and the maximum subsidence rate is more than 220 mm / yr.After 2017，the settlement deformation range and deformation value of the mining area has shown a significant downward trend，and the settlement deformation has been gradually stable.The preliminary treatment work has played an important role in controlling the settlement and deformation of the mining area.The study can provide an important reference for disaster prevention and reduction and treatment engineering evaluation in alpine mining areas.

Key words：

InSAR;mining area treatment;land subsidence;sedimentation rate