基于時序InSAR技術(shù)的礦區(qū)形變監(jiān)測與分析

陳志達，龐校光

（1.紹興市上虞區(qū)自然資源監(jiān)測中心，浙江 紹興 312399；2.中國科學院精密測量科學與技術(shù)創(chuàng)新研究院 大地測量與地球動力學國家重點實驗室，湖北 武漢 430077；3.中國科學院大學，北京 100049）

InSAR技術(shù)是20世紀后期發(fā)展起來的衛(wèi)星遙感新方法，可以獲取高精度地形和地表形變信息[1-2]。與常規(guī)地面水準測量、全球定位系統(tǒng)（GPS）、光學遙感等測繪手段相比，InSAR技術(shù)進行地面形變監(jiān)測時優(yōu)勢明顯[3-5]，且已在地表形變[6-8]、礦區(qū)[9-10]、滑坡[11]及地震形變[12]等研究與監(jiān)測中得到了廣泛的應用。

礦區(qū)開采引起的地面沉降是非常復雜的過程，其危害十分嚴重，研究礦區(qū)塌陷具有重要的意義[13-14]。2011年，Jiang[15]等利用時序InSAR對礦區(qū)煤火引起的地表形變進行了監(jiān)測，結(jié)果顯示時序InSAR技術(shù)對礦區(qū)煤火引起的地表形變具有可靠的應用潛力。2014年，Dong[16]等利用差分干涉測量短基線集時序分析技術(shù)（SBAS-InSAR）對淮南地下煤礦開采引起的地面沉降進行時空分析，證明了L波段InSAR可用于監(jiān)測礦區(qū)地下開采活動引起的地面沉降。同年，Zhao[17]等通過InSAR技術(shù)對大同不同采礦沉降階段及時間演變進行分析，發(fā)現(xiàn)1992—2011年間約80%煤礦發(fā)生形變。2015年，LI[3]等對徐州礦區(qū)進行實驗，再次驗證了In-SAR技術(shù)用于礦區(qū)形變監(jiān)測結(jié)果的可靠性。2019年，朱建軍[18]利用InSAR技術(shù)求取了煤礦開采沉陷的三維形變場，并進行了實驗模擬。

山西省大同市渾源縣煤炭資源豐富，從20世紀70年代以來，煤礦無序開采導致周邊地區(qū)大范圍地面沉降，嚴重破壞了周圍的生態(tài)環(huán)境，直接威脅著當?shù)鼐用竦纳柏敭a(chǎn)安全。為此，本文利用覆蓋大同市渾源縣46景SAR影像，采用時序InSAR技術(shù)進行分析處理，獲取了2007—2013年間渾源縣平均年形變速率結(jié)果和時序結(jié)果，并根據(jù)實地調(diào)查資料分析了地面沉降時空變化特征及其與采空區(qū)分布的關(guān)系。

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

1.1 研究區(qū)概況

渾源縣（39.67°～39.83°N、113.72°～113.92°E）地處黃土高原的邊緣地帶，位于山西省東北部大同盆地東南邊緣，東臨廣靈縣，南鄰繁峙縣，西毗應縣，北接云州區(qū)。渾源縣地處桑干河支流渾河中上游，全縣南北長57 km，東西寬約52 km，總面積1 968 km2。渾源縣地形復雜，東、南、北三面環(huán)山，境內(nèi)起伏較大，海拔最高的臥羊場山主峰海拔2 333 m，最低的渾河西部河床海拔1 026 m，相對高差1 307 m。境內(nèi)山脈較多，海拔均在2 000 m以上，山勢陡峭，溝壑縱橫[19]。

渾源縣礦產(chǎn)資源豐富，儲量大，包含煤炭、花崗巖、鐵礦等，已探明的煤礦儲量約17億t。自上世紀70年代以來，粗放式煤礦開采破壞了當?shù)氐纳鷳B(tài)環(huán)境，威脅著人民的生命財產(chǎn)安全[17,20]。本文選取整個渾源縣為研究區(qū)域，圖1為渾源縣概況及數(shù)據(jù)覆蓋情況。

圖1 研究區(qū)域和數(shù)據(jù)覆蓋情況

1.2 實驗數(shù)據(jù)

ALOS-1于2006-01-24發(fā)射升空，是日本新一代的對地觀測衛(wèi)星。其攜帶的L波段SAR適用于地表破壞嚴重的煤礦開采區(qū)域，可有效提取渾源縣礦區(qū)地表形變范圍和形變量信息。RADARSAT-2是一顆搭載C波段傳感器的高分辨率商用雷達衛(wèi)星，于2007-12-14發(fā)射升空，是RADARSAT-1的后續(xù)衛(wèi)星，在新的圖像獲取能力及性能方面，有長足進展。該衛(wèi)星具有最高1 m高分辨率成像能力，多種極化方式。

本文 獲取2007年6月—2011年1月36景ALOS PALSAR影 像 和2012年1月—2013年4月10景RADARSAT-2影像數(shù)據(jù)，以獲取渾源縣長期地表高程變化速率及時間序列。ALOS/PALSAR數(shù)據(jù)分辨率為10 m，入射角為34.3°。RADARSAT-2數(shù)據(jù)分辨率為30 m，入射角為21.25°。DEM數(shù)據(jù)采用NASA提供的SRTM，其空間分辨率為90 m，用以去除干涉對中的地形相位。此外，還搜集5大煤礦區(qū)域、99座小煤礦坐標以及煤礦的地質(zhì)報告等煤礦開采資料，用于InSAR監(jiān)測結(jié)果的驗證與分析。通過獲取渾源縣地面形變范圍和形變量信息，結(jié)合煤礦開采存檔、野外實地調(diào)查等輔助數(shù)據(jù)，進行相互驗證。

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

2.1 SBAS-InSAR基本原理

SBAS-InSAR方法原理是將D-InSAR形變結(jié)果作為觀測值，通過SVD方法將多個小基線集聯(lián)合進行最小二乘解求解，獲得的形變時間序列可有效反演出地表形變范圍和形變量信息。對某一區(qū)域不同時間獲取N+1幅SAR影像，選擇主影像后將其他影像配準至主影像，為避免時空失相干造成觀測誤差及影響監(jiān)測結(jié)果可靠性，設置時空基線閾值以保證干涉對的有效性，后利用外部DEM去除參考橢球相位與地形相位，對形變相位進行解纏及地理編碼得到地理坐標系下研究區(qū)域形變結(jié)果[1-2,21]。兩景SAR影像干涉相位可用φint表示為：

式中，φflat為參考橢球相位；φtopo為地形相位；φatm為大氣延遲相位；φorbit為軌道誤差；φnoise為噪聲相位；φdef為地表形變相位。對去除各項誤差后的差分干涉對地表形變相位進行解纏后，可獲取視線向的位移量研究區(qū)域各像元形變量。

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

本文對ALOS/PALSAR與RADARSAT-2影像數(shù)據(jù)分別處理，干涉對主輔SAR數(shù)據(jù)精配準精度均控制在0.1個像元內(nèi)，空間基線均為最大基線的45%，最大基線分別為1 984.26 m和-271.49 m，時間基線分別為365 d和180 d，最終分別生成70和24個干涉對。對各干涉對分別進行精密軌道數(shù)據(jù)優(yōu)化空間基線，消除軌道誤差，去除地平效應，利用外部SRTM-DEM去除地形相位，干涉圖濾波，相位解纏，消除大氣相位延遲等處理得到地表形變相位，對獲取的形變干涉圖采用多基線方法獲取時間序列[22-23]，流程圖如圖2所示。

圖2 SBAS-InSAR方法流程

3 結(jié)果與討論

3.1 地表形變結(jié)果

利用覆蓋渾源縣的ALOS/PALSAR和RADARSAT-2衛(wèi)星雷達遙感數(shù)據(jù)，采用SBAS-InSAR技術(shù)獲取了渾源縣地表形變結(jié)果，包括渾源縣年均形變速率與重點區(qū)域的累積形變量。根據(jù)InSAR地表形變結(jié)果，結(jié)合99座小煤礦范圍、5大煤礦位置和野外地面調(diào)查資料，圈定了49個采空區(qū)分布范圍，如圖3所示。針對49個采空區(qū)獲取了采空區(qū)歷年累計區(qū)域平均形變量統(tǒng)計結(jié)果，其時序形變過程如圖4所示。

圖4 渾源縣49個形變區(qū)域中心2007-2013年間時序結(jié)果

從時序圖上可知各采空區(qū)7 a間地表高程持續(xù)下降，最大累計形變達-211.21 cm。其中15、17、18、19、20等區(qū)域平均形變量在2010年前后發(fā)生了劇烈變化，從圖3可知，5個區(qū)域位于黃花灘鄉(xiāng)-大磁窯鎮(zhèn)-青瓷窯鄉(xiāng)交界處，屬渾源縣東部沉降群。上述5個區(qū)域自2010年后煤礦開采規(guī)模與速度均有所加劇，其中20號區(qū)域平均形變量為當年渾源縣最為嚴重區(qū)域。2012年地面沉降進一步加劇，區(qū)域平均沉降量超過-50 cm的有3個區(qū)域，分別為20號的-62.01 cm、36號的-60.34 cm、41號的-56.32 cm。2013年平均沉降量超過50 cm的區(qū)域為13號和20號，其值分別為-53.33 cm、-52.45 cm。

圖3 渾源縣年均沉陷速率及沉陷區(qū)分布范圍

為驗證InSAR監(jiān)測結(jié)果的可靠性，結(jié)合實地調(diào)查結(jié)果與煤礦開采記錄資料等輔助數(shù)據(jù)對結(jié)果進行驗證。實地調(diào)查了沙圪坨鎮(zhèn)西坪村、大仁莊鄉(xiāng)吳家圪坨村、蔡村鎮(zhèn)白道村、青磁窯鄉(xiāng)青磁窯村，如圖5所示，結(jié)合形變時間序列結(jié)果可知7 a間累計形變量分別達-90.0 cm、-141.22 cm、-125.51 cm、-75.77 cm。受煤礦開采影響，煤礦區(qū)域周圍居民遭遇不同程度的地面沉降危害，村莊房屋建筑出現(xiàn)不同程度受損情況。

圖5 部分村莊房屋損毀情況

3.2 采空區(qū)形變分布特征

從圖3可知，渾源縣7 a間在東部、北部、西北部分別形成了三大沉降漏斗群，其中東部范圍最廣，包括大磁窯鎮(zhèn)-黃花灘鄉(xiāng)-青瓷窯鄉(xiāng)北部-千佛嶺鄉(xiāng)北部-大仁莊鄉(xiāng)中部等廣闊區(qū)域，該區(qū)域分布著果子園、古瓷窯、白莊等眾多煤礦。對各區(qū)域時序形變結(jié)果中歷年不同形變量比例進行統(tǒng)計得到圖6，從圖中可知，隨著時間變化，渾源縣49個區(qū)域中歷年形變量大的比重不斷增加。2007—2013年間，沉降量＞-10 cm的區(qū)域從75.5%下降至2%，年沉降量超過-50 cm區(qū)域從無到有，其中2012年共3個，占6.1%；2013年2個，占4.1%。但從2007年看形變量較2008年大，分析其原因可能是當年數(shù)據(jù)從2007年6月起，僅半年數(shù)據(jù)，獲取的形變量受季節(jié)性變化影響較大，難以準確反映出全年的形變規(guī)律；2013年最大沉降值較2012年更小，可能有兩個原因，一是2013年采用的數(shù)據(jù)不足半年，由于數(shù)據(jù)量不足及季節(jié)性影響可能無法完整反映出渾源縣2013年形變結(jié)果；另一個原因是2012年政府對煤礦開采加以控制，地面形變情況在2013年得到改善，針對原因不確定性，后續(xù)采用完整數(shù)據(jù)驗證具體影響因素。

圖6 渾源縣歷年不同形變量所占比例

根據(jù)時序InSAR結(jié)果可知渾源縣地表形變在空間上呈散點式漏斗形狀，面積小而多，均呈持續(xù)下沉狀態(tài)。結(jié)合時序InSAR結(jié)果獲取了渾源縣采空區(qū)中心隨時間變化的空間演變結(jié)果（圖7）。從圖上可知渾源縣北部與西部在2007—2013年間中心形變發(fā)生明顯形變，其中部分采空區(qū)中心形變具有明顯階段性特征，以北部41號采空區(qū)為例，該區(qū)形變以2011年為分界線；2007—2011年間呈穩(wěn)步線性趨勢變化，發(fā)生突然加速。渾源縣東部在2007—2010年間采空區(qū)中心形變較為平穩(wěn)，無明顯突變，但2011—2013年該區(qū)域地表形變異常，通過與時序結(jié)果對比發(fā)現(xiàn)在2010年后地表發(fā)生明顯變化的15、17、18、19、20號采空區(qū)均位于東部漏斗群，7 a間最大累計形變超過-2.3 m。

圖7 渾源縣采空區(qū)中心歷年累計變化量

2007—2013年間，渾源縣三大沉降漏斗群均出現(xiàn)不同程度形變。在時間緯度，49個采空區(qū)均呈持續(xù)的、大量級的形變，最大累計沉降超過-2.3 m；在空間緯度，采空區(qū)分布于11個鄉(xiāng)鎮(zhèn)，其中以東部的黃花灘鄉(xiāng)、大磁窯鎮(zhèn)、青瓷窯鄉(xiāng)、千佛嶺鄉(xiāng)、大仁莊鄉(xiāng)最為嚴重。7 a間渾源縣地面沉降無論從沉降量級還是沉降范圍均發(fā)生了巨大的變化。

3.3 討論

為分析時序InSAR監(jiān)測結(jié)果可靠性，本文搜集了渾源縣采煤沉陷區(qū)村莊治理、大同市原煤產(chǎn)量與渾源縣煤炭銷量[24]等資料。首先對涉及采煤沉陷區(qū)治理的村莊分析，發(fā)現(xiàn)村莊主要分布在大仁莊鄉(xiāng)、大磁窯鎮(zhèn)、東坊城鄉(xiāng)、黃花灘鄉(xiāng)、青瓷窯鄉(xiāng)、千佛嶺鄉(xiāng)、蔡村鎮(zhèn)、沙圪坨鎮(zhèn)等鄉(xiāng)鎮(zhèn)，分布于渾源縣東部、北部、西北部，其空間分布上與根據(jù)時序InSAR結(jié)果圈定的49個采空區(qū)空間位置分布相一致。說明煤礦開采是渾源縣出現(xiàn)形變漏斗的主要原因。

此外，通過對大同市原煤產(chǎn)量與渾源縣煤炭銷量進行分析發(fā)現(xiàn)，自2008年起，當?shù)孛旱V開采呈現(xiàn)爆炸式增長。大同市原煤銷量從2007年的7 394萬t增長至2010年的9 403萬t，2015年達122 59萬t；渾源縣煤炭銷量從2008年的15.46萬t，增長至2011年1 060萬t，年均增長幅度達309%。當?shù)孛旱V的產(chǎn)量與銷量變化情況與部分采空區(qū)形變量在2010年前后發(fā)生突變趨勢相吻合（圖8）。經(jīng)調(diào)查發(fā)現(xiàn)，自2008年起，煤炭價格飛速上漲，煤礦開采日益嚴重，尤其開采難度較低的露天煤礦已初具規(guī)模，使得多數(shù)區(qū)域年形變量開始呈指數(shù)加速增長狀態(tài)。

圖8 大同市原煤產(chǎn)量與渾源縣煤炭銷量及沉降中心年均形變時序圖

4 結(jié)論

本文分析了2007年6月—2011年1月36景ALOS PALSAR影像和2012年1月—2013年4月10景RADARSAT-2影像數(shù)據(jù)，獲取了渾源地區(qū)地表形變速率結(jié)果和時序形變結(jié)果。通過對結(jié)果分析可知：

1）根據(jù)時序InSAR地表形變結(jié)果，結(jié)合99座小煤礦范圍、5大煤礦位置和野外地面調(diào)查資料，圈定49個采空區(qū)范圍及各采空區(qū)歷年累計區(qū)域平均形變量統(tǒng)計結(jié)果。7 a間渾源縣地表呈持續(xù)虧損狀態(tài)，累計形變量不斷加劇，且在部分區(qū)域2010年后形變呈突然加劇狀態(tài)。

2）2007—2013年間，渾源縣出現(xiàn)嚴重的地表形變，在東部、北部、西北部分別形成了三大沉降漏斗群，其中范圍最廣的東部區(qū)域分布著果子園、古瓷窯、白莊等眾多煤礦。隨著時間推移，渾源縣49個區(qū)域中歷年形變量大的比重不斷增加，沉降量級與沉降范圍均呈惡化趨勢。

3）根據(jù)時序InSAR結(jié)果圈定的形變區(qū)空間分布與采煤沉陷區(qū)開采量一致性較高，表明煤礦開采是造成地表沉陷加劇的主要原因。