基于D-InSAR的礦區(qū)地面沉降監(jiān)測

劉付剛， 張洪全， ?？■危?江曉林， 董 軍

(黑龍江科技大學(xué) 電子與信息工程學(xué)院，哈爾濱 150022)

基于D-InSAR的礦區(qū)地面沉降監(jiān)測

劉付剛， 張洪全， ?？■?， 江曉林， 董 軍

(黑龍江科技大學(xué) 電子與信息工程學(xué)院，哈爾濱 150022)

為獲得礦區(qū)地表的具體形變量和形變趨勢，根據(jù)礦區(qū)地表沉陷形變特征，分析了合成孔徑雷達(dá)差分干涉測量(D-InSAR)技術(shù)中SAR影像配準(zhǔn)、去平地效應(yīng)、干涉圖濾波以及相位解纏等關(guān)鍵問題。結(jié)合實(shí)際礦區(qū)SAR影像數(shù)據(jù)，利用相應(yīng)技術(shù)的經(jīng)典處理算法，獲取了煤炭開采引起的地面沉降數(shù)據(jù)。實(shí)驗結(jié)果與實(shí)地勘察結(jié)論基本一致，證明D-InSAR技術(shù)在地面沉降監(jiān)測方面具有強(qiáng)大優(yōu)勢。該研究可為礦區(qū)生產(chǎn)、整體規(guī)劃與數(shù)字化礦山建設(shè)提供技術(shù)支持。

D-InSAR； 去平地效應(yīng)； 干涉濾波； 地面沉降監(jiān)測

0 引 言

煤炭的不合理開采破壞了礦區(qū)原有的地質(zhì)結(jié)構(gòu)，極易引起地面形變和沉陷，導(dǎo)致人類生存環(huán)境質(zhì)量的嚴(yán)重破環(huán)。常規(guī)的礦區(qū)地面沉降監(jiān)測主要有GPS定位和水準(zhǔn)測量，由于礦區(qū)監(jiān)測地域廣、地表形變漸變性和長期性等特點(diǎn)，導(dǎo)致地表形變監(jiān)測和預(yù)測效果不佳，甚至失效。D-InSAR技術(shù)是20世紀(jì)90年代發(fā)展起來的一種可以高精度監(jiān)測大面積微小地面形變的新技術(shù)。在礦區(qū)地面監(jiān)測應(yīng)用方面，各發(fā)達(dá)國家科研人員利用D-InSAR技術(shù)取得了一系列的突出成果。Marco等于1995年利用ERS數(shù)據(jù)監(jiān)測美國Belridge油田的地面形變，并將監(jiān)測分析結(jié)果同傳統(tǒng)水準(zhǔn)法對比，在年均30～40 cm的下沉量中，最大高程誤差不足5 mm；此后，世界各國學(xué)者針對不同的煤礦開采區(qū)，如德國圖林根州鈾礦開采區(qū)、波蘭Upper Silesian煤田開采區(qū)、澳大利亞某煤田開采區(qū)等礦區(qū)的地面沉降進(jìn)行監(jiān)測分析研究[1]，得到了毫米級的垂直方向形變數(shù)據(jù)，進(jìn)一步證明了D-InSAR技術(shù)在礦區(qū)地面形變沉降監(jiān)測中的作用和良好的監(jiān)測效果。近年來，我國科研工作者利用D-InSAR技術(shù)開展礦區(qū)地面沉降研究也取得了可喜的成果[2-5]。2017年3月，李振洪團(tuán)隊利用我國“高分三號”衛(wèi)星獲取了高質(zhì)量重軌干涉合成孔徑雷達(dá)實(shí)驗數(shù)據(jù)，生成了我國第一幅衛(wèi)星干涉SAR影像，并從影像中提取到了亞厘米級地面沉降沉息，開啟了我國自主雷達(dá)衛(wèi)星InSAR形變監(jiān)測的新進(jìn)程[6-7]。筆者針對合成孔徑雷達(dá)差分干涉測量中關(guān)鍵問題進(jìn)行分析和討論，并結(jié)合實(shí)際礦區(qū)SAR影像數(shù)據(jù)進(jìn)行實(shí)驗，以驗證D-InSAR在礦區(qū)地表形變監(jiān)測方面的應(yīng)用效果。

1 D-InSAR成像基本原理

D-InSAR技術(shù)有別于傳統(tǒng)的基于點(diǎn)觀測的地面測量方法，該技術(shù)融合了合成孔徑雷達(dá)成像技術(shù)和干涉測量技術(shù)，利用單軌道雙天線模式(兩副天線同時觀測)或單天線重復(fù)軌道模式(兩次近平行觀測)獲得同一地區(qū)的兩景數(shù)據(jù)，通過獲取同一目標(biāo)對應(yīng)的兩個回波信號之間的相位差，結(jié)合軌道數(shù)據(jù)來獲取高精度、高分辨率的地面高程信息，精確測量地表某一點(diǎn)的三維空間位置及微小變化。星載重復(fù)軌道成像示意如圖1所示[8]。

圖1中b為基線距，基線距與水平方向的傾角為α,S1和S2分別表示主輔圖像傳感器，主圖像入射角為θ，H為主傳感器相對地面高度，R和R+ΔR分別為主輔圖像斜距，P為地面目標(biāo)點(diǎn)，其高程為h。設(shè)地面目標(biāo)點(diǎn)P兩次成像時的圖像分別為c1和c2,則可表示為

c1=A1ejφ(R),

(1)

c2=A2ejφ(R+ΔR),

(2)

式中:A1、A2——c1、c2的振幅；

φ(R)、φ(R+ΔR)——c1、c2的相位。

圖1 InSAR成像示意

忽略大氣效應(yīng)等因素，可以認(rèn)為相位信息由兩部分組成，一是與信號往返路徑相對應(yīng)的相位，二是散射特性引起的相位。

假設(shè)主輔圖像傳感器成像散射特性不變，對兩圖像信號進(jìn)行復(fù)共扼相乘可得到復(fù)干涉紋圖，即

(3)

式中：*——共軛； λ——雷達(dá)波的波長。

由式(1)～(3)可以得到干涉圖的相位為

(4)

根據(jù)圖1所示，利用三角函數(shù)可得

(5)

z=H-Rcosθ,

(6)

展開式(5)，由于ΔR遠(yuǎn)遠(yuǎn)小于R，故可得

ΔR≈bsin(θ-α),

(7)

將式(7)帶入(4)，可得

(8)

對干涉相位進(jìn)行解纏可以得到式(8)中的N?，F(xiàn)只考慮干涉圖相位的主值有

(9)

將式(9)帶入到(6)得

至此，可以求出目標(biāo)點(diǎn)的高度值h。由此可見，InSAR技術(shù)屬于主動微波遙感系統(tǒng)，工作機(jī)理和方式與傳統(tǒng)光學(xué)不同，它不依賴于可見光, 能全天時全天候工作,精度相對較高。

2 研究區(qū)概況和數(shù)據(jù)來源

濟(jì)寧市北部礦井，煤田南北長4 km，東西寬6 km。礦區(qū)內(nèi)地勢東高西低，地貌較為復(fù)雜，地面標(biāo)高均值是+38 m，氣候溫和，歷年最大積雪0.15 m，最大凍土0.31 m。兗石鐵路和新兗鐵路貫穿東西，交通便利。文中所用數(shù)據(jù)為歐空局ENVISAT衛(wèi)星提供的2008年12月至2009年7月的6景C波段ASAR數(shù)據(jù)，影像覆蓋濟(jì)寧市及北部煤田區(qū)域。影像數(shù)據(jù)所采用的雷達(dá)為ASAR型，軌道號為35571和35573～35577，軌跡為447，幀號為2889。

3 D-InSAR關(guān)鍵技術(shù)分析

D-InSAR的關(guān)鍵技術(shù)包括圖像配準(zhǔn)、產(chǎn)生干涉圖、濾波、平地效應(yīng)消除、相位解纏和差分干涉處理等。D-InSAR的經(jīng)典方法為三軌法，采用三景雷達(dá)圖像，以其中之一作為主圖像，另外兩景為從圖像，產(chǎn)生干涉圖。圖2給出了三軌法差分干涉測量數(shù)據(jù)流程。

圖2 三軌法差分干涉測量數(shù)據(jù)流程

3.1 SAR圖像配準(zhǔn)

兩幅圖像在空間和灰度上的映射即為圖像配準(zhǔn)。在選取公共主影像之后，為得到干涉條紋圖，必須使得主從兩幅SAR圖像上的點(diǎn)對應(yīng)地面上的同一點(diǎn)，即將從影像配準(zhǔn)到主影像的影像坐標(biāo)系內(nèi)。由于地形對相位的影響非常敏感，加之時間基線和空間基線的影響，很容易導(dǎo)致一些影像由于相干性過低而無法得到較高的圖像配準(zhǔn)質(zhì)量。

圖像對的配準(zhǔn)包括粗配準(zhǔn)和精配準(zhǔn)兩個步驟。粗配準(zhǔn)就是利用衛(wèi)星軌道參數(shù)或人工選取少量的控制點(diǎn)計算主、輔圖像之間在方位向和距離向上的坐標(biāo)映射關(guān)系，其配準(zhǔn)精度約為30個像元。粗配準(zhǔn)常采用基于衛(wèi)星DEOS精密軌道參數(shù)的配準(zhǔn)方法[9]。精配準(zhǔn)又分為像元級配準(zhǔn)和亞像元級配準(zhǔn)。兩者的基本步驟和方法大體相同，一般來說，亞像元級圖像配準(zhǔn)精度要優(yōu)于1/8個像元才能符合SAR干涉處理的精度要求，此時所造成的失相干很小(約為4%)。

3.2 去平地效應(yīng)

影像配準(zhǔn)后對應(yīng)的點(diǎn)共軛相乘，可以獲得干涉相位圖。干涉相位圖并不能直觀地體現(xiàn)地形高低的變化，這是因為干涉相位圖的相位信息里包含著平地(高度不變)所對應(yīng)的相位量。這種由高度不變的平地引起干涉相位在距離向和方位向呈周期性變化的現(xiàn)象稱為平地效應(yīng)，它在時域內(nèi)顯示沿方位向平行的明暗相間的密集條紋。平地效應(yīng)所引起的條紋和地形高低變化所引起的條紋相互疊加，掩蓋了地形高低變化對應(yīng)的相位信息，增加了相位解纏難度，必須進(jìn)行去除平地效應(yīng)處理。平地相位去除后，干涉圖的相位信息就能反映真實(shí)的地表高程相位和形變相位。針對ASAR數(shù)據(jù)，衛(wèi)星搭載的平臺飛行軌道穩(wěn)定性較ERS-I衛(wèi)星要差，飛行軌道之間存在著一定的夾角，其方位向也存在著平地效應(yīng)。采用精密軌道參數(shù)和成像區(qū)域中心點(diǎn)位置去除平地效應(yīng)干涉圖如圖3所示。

a 去除前

b 去除后

由圖3可見，去平地效應(yīng)之后的干涉紋圖條紋密度較去除前稀疏了很多，并且條紋的清晰度有明顯改善，去平地效應(yīng)后，有利于噪聲濾波和相位解纏的進(jìn)行。

3.3 干涉圖濾波

考慮到SAR圖像噪聲(通常定義為滑動窗口內(nèi)的異常最大值或最小值)對圖像重構(gòu)的影響，還需要對成像后干涉圖進(jìn)行濾波處理。目前常用的濾波方法分為空間域濾波和頻率域濾波[10]?？臻g域濾波方法是利用圖像像素的空間相關(guān)性對相干斑進(jìn)行濾波，如多視平均處理方法、均值濾波處理方法、中值濾波處理方法等。頻率域的方法主要是利用小波變換，比較著名的有基于小波變換和多尺度分析的濾波方法。理想的SAR圖像濾波算法應(yīng)該能夠在抑制斑點(diǎn)噪聲的同時保持圖像的邊緣和紋理信息。采用小波變換的濾波方法得到的干涉圖如圖4所示。

a 濾波前

b 濾波后

對比圖4a、b可以看出，采用小波變換的濾波方法濾波具有一定的效果，一些噪聲明顯濾除，濾波后的干涉圖效果明顯改善。利用工具Cursor Location/Value，查看圖4的像元值，可以看到相干性系數(shù)分布在0～1，圖像中越白的點(diǎn)，相干性系數(shù)值越大，說明該區(qū)域的相干性越高。圖5給出了截取濾波前后兩幅圖中第250行數(shù)據(jù)的剖面。由圖5可見，經(jīng)過濾波去噪處理后圖像曲線變化平穩(wěn)，數(shù)據(jù)處理過程中濾掉了大量的噪聲數(shù)據(jù)。

a 濾波前

b 濾波后

3.4 相位解纏

去平地效應(yīng)和濾波處理之后的干涉圖中，相位值并非真實(shí)的絕對相位值，而是以非線性的方式混合纏繞了2π的整數(shù)倍，如式(4)所示，假設(shè)二維矩陣φi,j表示纏繞相位，為了解纏該二維矩陣，對矩陣中的每個點(diǎn)(i,j)加減2π的整數(shù)倍，可以得到一個連續(xù)的解纏相位函數(shù)φi,j,即

φi,j=φi,j+2πki,j，i[0,M-1]，j[0,N-1],

其中,Ki,j為整數(shù)，-π≤φi，j<π。

相位解纏的準(zhǔn)確性直接關(guān)系到地表形變量的測量精度。目前，相位解纏方法主要分為三類：路徑跟蹤法、最小范數(shù)法基于網(wǎng)絡(luò)規(guī)劃的相位解纏算法。利用加權(quán)最小二乘法解纏干涉圖如圖6。

a 解纏前

b 解纏后

圖6中解纏后相位圖中干涉條紋明顯比原始數(shù)據(jù)少，其相位值是一個比原始纏繞相位值更加平滑的估計值，解纏后的相位值更接近真實(shí)的相位值。

3.5 實(shí)驗結(jié)果分析

綜合上述實(shí)驗數(shù)據(jù)及實(shí)驗過程，根據(jù)圖4 的干涉圖可以得出礦區(qū)局部地區(qū)的地表形變結(jié)果。在2008年12月至2009年7月之間，煤田北部東西長0.8 km，南北長1.5 km塌陷較為突出，最大塌陷處達(dá)3.0 cm，其他地區(qū)的沉降量處于[0.7,2.2]區(qū)間內(nèi)，沉降趨勢較緩慢。最大沉降地東北方向1.4 km處有兩個村莊，地下煤層的開采引起的地質(zhì)形變將會影像到村莊的建筑物，需要加強(qiáng)注漿充填開采的力度。經(jīng)過實(shí)地勘察和調(diào)查資料對干涉測量結(jié)果進(jìn)行驗證，發(fā)現(xiàn)分析結(jié)果與實(shí)際情況基本一致。由此可見，D-InSAR技術(shù)監(jiān)測礦區(qū)地表沉降效果較好，與傳統(tǒng)的地表形變監(jiān)測方法相比較，該技術(shù)具有覆蓋范圍廣、監(jiān)測成本低、監(jiān)測效率高的明顯優(yōu)勢。

4 結(jié)束語

文中采用歐空局ENVISAT衛(wèi)星提供的C波段ASAR數(shù)據(jù)，利用D-InSAR技術(shù)中的典型算法對濟(jì)寧北部煤田區(qū)域進(jìn)行了地表形變監(jiān)測，并將實(shí)驗結(jié)果通過實(shí)地勘察和資料查證進(jìn)行驗證，得出了與實(shí)際結(jié)果一致的結(jié)論。通過分析可知，D-InSAR技術(shù)在礦區(qū)地表沉降形變監(jiān)測應(yīng)用中完全可行，并且具有明顯的技術(shù)優(yōu)勢，利用該技術(shù)能夠為礦產(chǎn)資源的開采提供科學(xué)合理的決策參考，同時可以科學(xué)預(yù)測可能出現(xiàn)的地質(zhì)災(zāi)害，保證和諧礦區(qū)、綠色礦區(qū)的可持續(xù)發(fā)展?？梢灶A(yù)見，未來隨著SAR衛(wèi)星的不斷專業(yè)化、InSAR處理軟件和算法的日趨成熟，D-InSAR技術(shù)的實(shí)際應(yīng)用水平將會越來越高。

[1] Carnec C，Massonnet D，King C. Two examples of the use of SAR interferometry on displacement of small spatial extent[J].Geophysical Research Letters, 1996，23(4): 3579-3582.

[2] 楊澤發(fā), 易輝偉, 朱建軍, 等. 基于InSAR時序形變的礦區(qū)全盆地沉降時空演化規(guī)律分析[J]. 中國有色金屬學(xué)報, 2016, 26(7): 1515-1522.

[3] 陸燕燕, 柯長青, 陳德良, 等. PS-InSAR在沛縣礦區(qū)地表沉降監(jiān)測中的應(yīng)用[J]. 地理空間信息, 2016, 14(5): 96-99.

[4] Saghafi A, Pinetown K, Javanmard H. Developing a new method to identify the source of gas emissions into longwall and goaf from sSurrounding strata[C]// Coal Operators’ Conference 2015. Wollongong, Australia: The University of Wollongong, The Australasian Institute of Mining and Metallurgy and Mine Managers Association of Australia, 2015.

[5] 汪寶存, 郭凌飛, 王軍見, 等. 礦區(qū)地表形變InSAR監(jiān)測——以永城市為例[J]. 測繪與空間地理信息, 2016, 39(6): 24-27.

[6] 丁赤飚, 劉佳音, 雷 斌, 等. 高分三號SAR衛(wèi)星系統(tǒng)級幾何定位精度初探[J]. 雷達(dá)學(xué)報, 2017, 6(1): 11-16.

[7] 劉淮宇. 高分專項應(yīng)用推廣全面啟動[J]. 太空探索, 2016(2): 7-7.

[9] Rudzinski L, Mirek K, Mirek J. Source mechanism analysis of strong mining induced seismic event and its influence on ground deformation observed by InSAR technique[C]//EGU General Assembly 2016, Vienna, Austria： European Geophysical Society and the European Union of Geoscienes, 2016，18: 6759.

[10] Chaussard E, Kerosky S. Characterization of black sand mining activities and their environmental impacts in the philippines using remote sensing[J]. Remote Sensing, 2016, 8(2): 100-116.

(編校 王 冬)

D-InSAR-based monitoring of land subsidence in mining area

LiuFugang,ZhangHongquan,ChangJunxin,JiangXiaolin,DongJun

(School of Electronics & Information Engineering, Heilongjiang University of Science & Technology, Harbin 150022, China)

This paper provides an insight into the specific mining surface deformation and deformation trend. The research building on the characteristics behind subsidence deformation in mines, works towards the following objectives：analyzing the key issues, such as SAR image registration, flat-earth removal, interferogram filtering, and phase unwrapping, all of which are involved in synthetic aperture radar differential interferometry(D-InSAR)；combined with SAR image data derived from actual mining area, employing the classical processing algorithm of corresponding technology to carry out tests and thereby derive the data on land subsidence caused by coal mining; and providing the results basically consistent with the field investigation results to verify the strong advantage D-InSAR technology has in land subsidence monitoring. The research could provide technical support for mining area production, overall planning and digital mine construction.

D-InSAR; flat-earth removal; interference filter; ground subsidence monitoring

2017-04-03

黑龍江省教育廳科學(xué)技術(shù)研究項目(12541716)

劉付剛(1981-)，男，山東省德州人，副教授，博士，研究方向：信號檢測處理、InSAR遙感遙測，E-mail：liufugang_36@163.com。

10.3969/j.issn.2095-7262.2017.03.013

TD173.4

2095-7262(2017)03-0265-05

A