基于SBAS-InSAR技術(shù)的礦區(qū)沉降監(jiān)測研究

摘要：本文介紹了SBAS-InSAR 技術(shù)的數(shù)據(jù)處理流程，結(jié)合ALOS PALSAR-2影像數(shù)據(jù)，分析了礦區(qū)沉降形變特征。研究結(jié)果表明：礦區(qū)沉降區(qū)域空間模式主要表現(xiàn)為橢圓形，且隨著時(shí)間推移沉降范圍在逐漸擴(kuò)大，沉降量級也逐漸增加，研究區(qū)域最大垂直向形變量達(dá)到了-0.20m，垂直向的最大年平均形變速率可達(dá)-0.15m/y。從空間屬性上來看，最大形變峰值處位于沉降盆地幾何中心，符合礦山開采沉陷理論，也較好地反映了監(jiān)測結(jié)果的可靠性。

關(guān)鍵詞：SBAS-InSAR；礦區(qū)；沉降監(jiān)測

一、引言

傳統(tǒng)的礦區(qū)沉降監(jiān)測方法有GPS與水準(zhǔn)測量，這兩種傳統(tǒng)的測量方法具有毫米級的監(jiān)測精度，但整體監(jiān)測成本高，維護(hù)工作量大，且監(jiān)測目標(biāo)為單點(diǎn)，只能通過單點(diǎn)預(yù)估出大概的沉降 [1]。合成孔徑雷達(dá)干涉測量技術(shù)（Interferometric Synthetic Aperture Radar，InSAR）是用于獲取地表三維信息的主動(dòng)式微波遙感新興技術(shù)[2]，它是通過獲取同一地物目標(biāo)的兩幅合成孔徑雷達(dá)（Synthetic Aperture Radar，SAR）復(fù)數(shù)圖像，并將兩幅 SAR 影像共軛相乘來得到包含有地表高程信息的干涉相位圖。然后通過數(shù)據(jù)處理，去除干涉相位圖中的高程信息來獲得干涉相位圖中的形變信息[4-6]。該方法能提供區(qū)域的大范圍連續(xù)性形變監(jiān)測，具有非接觸、全天時(shí)、全天候、空間連續(xù)覆蓋等優(yōu)點(diǎn)。

本文將InSAR技術(shù)應(yīng)用于礦區(qū)地表形變監(jiān)測，結(jié)合礦區(qū)特點(diǎn)采取SBAS處理方法，并且結(jié)合實(shí)地觀測資料對其檢驗(yàn)分析，證明SBAS-InSAR技術(shù)在礦區(qū)地表沉陷監(jiān)測的可靠性。

二、SBAS-InSAR的基本原理

SBAS-InSAR的數(shù)據(jù)處理流程如圖1所示。①將首先構(gòu)建小基線干涉對；②根據(jù)二軌差分干涉原理，對所有干涉對進(jìn)行干涉、去平地效應(yīng)、去地形、濾波、相位解纏等處理；③在高相干點(diǎn)上建立線性形變速率，高程誤差的線性模型，估算線性形變速率和高程殘差；④采用合適濾波算法對大氣相位和地形殘余相位進(jìn)行去除；⑤將形變反演結(jié)果轉(zhuǎn)換到合適地理坐標(biāo)系下；⑥通過LOS向和垂直向的轉(zhuǎn)換因子，獲取垂直向的形變信息。

三、研究區(qū)域概況

從地形地貌上來看，礦區(qū)地處山脈北段，山脊線呈長條狀，該山脊的延伸方向與背斜軸向基本一致。區(qū)內(nèi)山脊較為平穩(wěn)，兩側(cè)為迭瓦狀的順向坡，地形坡度為9°-25°，一般表現(xiàn)為西坡陡，東坡緩。此外，礦區(qū)內(nèi)無大的地表積水體，多為流量較小的溪溝。

另外，根據(jù)多年氣象觀測資料，該區(qū)屬于亞熱帶濕潤型氣候，氣候溫和潮濕，四季分明，雨量適中。區(qū)域年平均氣溫為17.3℃，最高氣溫達(dá)42.3℃，最低氣溫為-4.7℃。年平均總降水量為1208.3mm，年均水面蒸發(fā)量1085.0mm。

四、數(shù)據(jù)處理及結(jié)果分析

（一）數(shù)據(jù)來源

本實(shí)驗(yàn)數(shù)據(jù)為2018年10月18日至2020年5月28日的5景ALOS PALSAR-2影像（其中2019年9月5日的影像為SM1成像模式，其余四景為SM3成像模式）。

（二）數(shù)據(jù)處理

1.根據(jù)小基線原則生成SAR數(shù)據(jù)干涉對和時(shí)空基線連接圖

本實(shí)驗(yàn)選取2018年10月18日的影像為主影像，通過設(shè)置時(shí)空基線閾值，共生成6個(gè)小基線干涉對。干涉對的具體信息見表2。

2.差分干涉處理

干涉對組合確定后，根據(jù)二軌差分干涉原理，對所有干涉對進(jìn)行干涉、去平地效應(yīng)、去地形、濾波、相位解纏等處理。

由于本實(shí)驗(yàn)所選用的ALOS PALSAR-2掃描模式（SM1）幅寬為50KM，且研究區(qū)域地形多起伏，故在干涉時(shí)采用10：10多視，以犧牲分辨率為代價(jià)提高信噪比及處理效率?？紤]到時(shí)效性，選用TanDEM-X 90m DEM數(shù)據(jù)模擬地形相位。根據(jù)計(jì)算，影像覆蓋區(qū)域最大高差約 480米。

在獲得地形模擬相位后，即可從干涉圖中扣除地形相位以獲取地表形變相位（實(shí)際上還混疊有其他誤差）。然而由于各種失相干作用導(dǎo)致差分干涉圖中存有不同程度的噪聲，影響之后相位解纏質(zhì)量和沉降監(jiān)測精度，此處采用Goldstein 進(jìn)行濾波。濾波后的影像能夠有效削弱干涉相位中的絕大部分噪聲，從而為后續(xù)相位解纏提供更高質(zhì)量的基礎(chǔ)數(shù)據(jù)。在解纏處理中，本實(shí)驗(yàn)選用最小費(fèi)用流解纏分別對各差分干涉圖進(jìn)行解纏。

3.線性形變建模與估計(jì)

對所有構(gòu)建的干涉對實(shí)現(xiàn)差分處理后，依據(jù)相干系數(shù)選擇時(shí)序穩(wěn)定的高相干目標(biāo)。其后在高相干點(diǎn)上建立線性形變速率，高程誤差的線性模型，然后在線性模型的基礎(chǔ)上構(gòu)建方程組，利用矩陣的奇異值分解（SVD）法，求解方程組并估算線性形變速率和高程殘差。

4.殘余相位估算和誤差去除

通過形變速率去除后，殘余相位中主要包含了噪聲和大氣延遲湍流部分。為此，本實(shí)驗(yàn)采用高通濾波和低通濾波削弱大氣延遲，從而提高時(shí)序形變精度。其中大氣低通濾波算法和大氣的空間特性相關(guān)，由于大氣在空間上具有高相干性，采用空間窗口進(jìn)行濾波時(shí)，窗口越小，其濾波效果越強(qiáng)；大氣高通濾波算法和大氣的時(shí)間特性有關(guān)，由于大氣在時(shí)間上具有低相干性，采用時(shí)間窗口進(jìn)行濾波時(shí)，窗口越大，其濾波效果越強(qiáng)。因此，本實(shí)驗(yàn)選擇空間域的低通濾波和時(shí)間域的高通濾波算法對大氣相位和地形殘余相位進(jìn)行去除。

5.地理編碼

上述過程均是在SAR坐標(biāo)系下進(jìn)行的，地理編碼可以將形變反演結(jié)果轉(zhuǎn)換到地理坐標(biāo)系下，以便和實(shí)測結(jié)果進(jìn)行對比和分析。圖2、圖3為2018年10月18日至2020年5月28日礦區(qū)的LOS向的累計(jì)形變圖以及形變速率圖。

6.垂直向形變分解

InSAR技術(shù)獲取的是雷達(dá)視線方向的形變，該方向的形變往往不能準(zhǔn)確反映礦區(qū)的形變特征，因此，可以通過LOS向和垂直向的轉(zhuǎn)換因子（入射角的余弦值 cosθ），獲取垂直向的形變信息。為顧及InSAR成像幾何結(jié)構(gòu)對沉降監(jiān)測的影響，此處忽略水平移動(dòng)將LOS向形變直接轉(zhuǎn)化至垂直向，形變結(jié)果如圖4、圖5所示。

（三）結(jié)果分析

圖2中白色實(shí)線標(biāo)注的區(qū)域?yàn)榈V區(qū)所在區(qū)域。結(jié)果清晰地顯示了20余個(gè)地表沉降區(qū)域（紅色區(qū)域），其空間模式主要表現(xiàn)為橢圓形，且隨著時(shí)間推移沉降范圍在逐漸擴(kuò)大，沉降量級也逐漸增加，最大LOS向形變量級達(dá)到了-0.16m。可以推測，這些地表沉降很可能為地下開采所致，因?yàn)槠鋾r(shí)空演化模式與地下開采形變演化模式非常相似。

圖3為礦區(qū)及周邊地區(qū)2018年10月18日至2020年5月28日期間的LOS向年平均形變速率圖。從中可以看出，形變速率的空間分布與形變累積分布情況基本相同，LOS向最大年平均形變速率可達(dá)-0.12m/y。

圖4為垂直向形變圖。由圖可知，主要形變區(qū)域的空間分布與LOS 向形變區(qū)域分布基本相同，垂直向的最大累積形變量可達(dá)-0.20m。

圖5為礦區(qū)及周邊地區(qū)2018年10月18日至2020年5月28日期間的垂直向年平均形變速率圖。由圖可知，垂直向的最大年平均形變速率可達(dá)-0.15m/y，數(shù)值高于LOS向的形變速率。

為直觀分析礦區(qū)開采引起地表形變時(shí)序特征，本實(shí)驗(yàn)在礦區(qū)開采形成的地表沉陷盆地上方繪制剖面線D-D（圖6所示），隨后提取各期沉陷量并繪制剖面圖，具體見圖7。

如圖7所示，（a）和（b）分別表示時(shí)序LOS向累積形變及時(shí)序垂直向累積形變，從時(shí)間屬性上來看，自2018年10月18日開始，剖面累積形變隨時(shí)間推移逐漸加大，在LOS 向上形變峰值可達(dá)-0.10m，轉(zhuǎn)化至垂直向后最大峰值可達(dá)-0.12m，并且對照（a）和（b）可知垂直向累積形變趨勢與LOS向累積形變趨勢一致。此外，2019年9月5日至2020年5月28日期間，剖面累積形變增加顯著（LOS向最大形變由0.06m增至0.10m）。

從空間屬性上來看，結(jié)合累積形變與距離的關(guān)系可知最大形變峰值處大致位于沉降盆地幾何中心，較為符合礦山開采沉陷理論。此外，對照時(shí)序InSAR獲取的全盆地沉降演化圖可知剖面線橫跨兩個(gè)地表移動(dòng)盆地，在圖8中體現(xiàn)出了兩個(gè)沉降峰值，也較好地反映了監(jiān)測結(jié)果的可靠性。

五、結(jié)束語

本文基于ALOS PALSAR-2數(shù)據(jù)，利用SBAS-InSAR技術(shù)，得到某礦區(qū)2018-2020年各時(shí)間段的地表形變情況。結(jié)果顯示，研究區(qū)域存在20余個(gè)地表沉降區(qū)域，其空間模式主要表現(xiàn)為橢圓形，且隨著時(shí)間推移沉降范圍在逐漸擴(kuò)大，沉降量級也逐漸增加，垂直向的最大累積形變量可達(dá)-0.20m，垂直向的最大年平均形變速率可達(dá)-0.15m/y。監(jiān)測結(jié)果符合礦山開采沉陷理論并與實(shí)際采空區(qū)位置基本吻合，證明監(jiān)測結(jié)果具有一定的有效性。

作者單位：況長虹? ? 四川 成都

參? 考? 文? 獻(xiàn)

[1] 劉冬，劉星年，陳超.基于CR-InSAR的煤礦區(qū)地表沉陷監(jiān)測研究[J].現(xiàn)代測繪，2015，38（4）：26-30.

[2] 楊成生.基于D-InSAR技術(shù)的煤礦沉陷監(jiān)測[D].長安大學(xué)， 2008.

[3] 劉永坦.雷達(dá)成像技術(shù)[M]. 哈爾濱工業(yè)大學(xué)出版社， 2014.

[4] 王遠(yuǎn)堅(jiān)，齊麟，姜岳.基于D-InSAR技術(shù)的金礦老采空區(qū)下沉監(jiān)測與位置反演[J]. 金屬礦山， 2019（10）：68-73.

[5] 范軍林，劉飛鵬，涂梨平.基于D-InSAR技術(shù)的南昌市地表形變監(jiān)測研究[J].礦產(chǎn)與地質(zhì)，2019，33（3）：546-550.

[6] 陳優(yōu)良，王之玉，栗正剛，等.基于D-InSAR的礦山地表沉陷監(jiān)測研究[J].化工礦物與加工，2019，48（8）：10-14.