祝傳廣,鄧喀中,張繼賢,張永紅,范洪冬,張立亞
(1.中國礦業(yè)大學國土環(huán)境與災害監(jiān)測國家測繪局重點實驗室,江蘇徐州 221116;2.中國測繪科學研究院,北京 100036;3.煤炭資源清潔利用與環(huán)境保護湖南省重點實驗室,湖南湘潭 411100)
基于多源SAR影像礦區(qū)三維形變場的監(jiān)測
祝傳廣1,2,鄧喀中1,張繼賢2,張永紅2,范洪冬1,張立亞3
(1.中國礦業(yè)大學國土環(huán)境與災害監(jiān)測國家測繪局重點實驗室,江蘇徐州 221116;2.中國測繪科學研究院,北京 100036;3.煤炭資源清潔利用與環(huán)境保護湖南省重點實驗室,湖南湘潭 411100)
針對DInSAR(differential interferometric synthetic aperture radar)技術僅能獲取雷達視線向(line of sight,LoS)形變的不足,研究了融合多衛(wèi)星平臺求解三維形變場的模型與算法。該算法基于多衛(wèi)星軌道模式下具有不同成像幾何的多源SAR影像聯(lián)合求解礦區(qū)地表形變場。研究結果表明:采用該算法反演的下沉值與水準測量結果相互吻合,均方根誤差為±4 mm,吻合程度優(yōu)于單一影像源反演結果;垂直向位移場與等值線均表明下沉盆地向老采空區(qū)偏移,說明老采空區(qū)可能活化;東西向水平位移場與等值線符合開采沉陷地表移動規(guī)律,而且對于不同的成像模式,東西向水平移動的影響亦不同;由于衛(wèi)星航向角的正弦值近乎為0,使得三維算法對南北向位移不敏感。
合成孔徑雷達;多源影像;雷達視線向;水平移動;開采沉陷
Key words:synthetic aperture radar;multi-source imagery;line of sight;horizontal displacement;mining subsidence
開采地下煤炭資源會破壞巖層內部應力平衡力,使得地表出現(xiàn)下沉、傾斜等變形。為保護地面建筑物并指導地下開采工作,進行地表變形觀測是必要的。目前,觀測方法主要為水準測量與GPS測量,雖然精度較高(靜態(tài)GPS測量可達mm級精度,水準測量可優(yōu)于mm級),但是觀測點的空間密度低、費用高。在Gabriel成功利用合成孔徑雷達差分干涉(differential SAR interferometrical,DInSAR)技術獲取地表形變[1],證明DInSAR技術進行變形監(jiān)測的可行性之后,該技術憑借廣空間覆蓋、高空間采樣率以及快速重復觀測等優(yōu)點,廣泛應用于地震、火山、滑坡等自然災害[2-6]的監(jiān)測中。隨后,為研究長時間跨度的形變信息,基于時間序列的DInSAR技術被陸續(xù)提出,如PSInSAR技術[7-8]、SBAS-InSAR技術[9-10]等,并且在監(jiān)測地下水資源[11-12]開采等引起的大范圍地表位移方面取得了成功。
然而,DInSAR技術僅能獲取目標在雷達視線向(line of sight,LoS)的一維形變。在研究目標的垂直向形變時,通常忽略水平移動,然后根據(jù)雷達波入射角將LoS向位移投影到垂直向。這對于大范圍沉降是可行的,如抽取地下水導致的地表沉降范圍通常比較大,使得小范圍內的水平移動可以忽略。但是在小范圍內發(fā)生大的水平移動是礦區(qū)地表形變的特點,因此需要采用新的技術手段研究礦區(qū)的變形規(guī)律。
DInSAR技術獲取的LoS向位移可以分解成垂直向、東西向與南北向3個位移分量,當利用至少3種具有不同成像幾何的影像源獲取目標在不同LoS向的位移后,就可以聯(lián)合解算目標的三維形變場?;谶@種思想,根據(jù)兩種衛(wèi)星傳感器ASAR與PalSAR,選取了3種具有不同成像幾何的SAR影像源,分別獨立解算出LoS向位移,然后聯(lián)合求解研究區(qū)域的三維形變場,并研究分析了東西向水平移動對單一影像源反演的LoS向位移的影響。
對于重復軌道干涉測量模式,如果地面目標k在兩次成像期間發(fā)生了位移,則其干涉相位φ可記為
其中,ω表示纏繞算子;φlos為LoS向位移相位;φatm, φflat,φtopo與φn分別為大氣變化、平地效應、地形與噪聲相位。結合各相位[13-14]自身的特點,可以利用DInSAR技術[15]從式(1)中分離出φlos,也就獲取了目標在LoS向的位移Dlos。在三維空間中,Dlos可以分解成垂直向分量Dv、東向分量De及北向分量Dn[16-17],如圖1所示。
圖1 LoS向位移三維分解模型Fig.1 Decomposition of the displacement along the LoS
圖1中,航向角α指衛(wèi)星飛行方向與北方向間的夾角,由北方向開始指向飛行方向。θ是研究區(qū)域中心處目標k處的雷達波入射角,下標sl表示目標k處的LoS在地面上的投影,Dn,sl表示北向位移分量Dn在sl上的投影,De,sl表示東向位移分量De在sl上的投影。Dsl則是Dn,sl與De,sl的矢量和,其在LoS向的投影為Dsl,los。Dv,los表示垂直向位移分量Dv在LoS向的投影,Dsl,los與Dv,los的矢量和即為Dlos。因此,根據(jù)圖1所示的三維分解原理,目標k的LoS向位移Dlos可記為
式(2)中有3個未知量即Dv,De與Dn,并且這3個未知量是目標k在三維空間中的真實位移,與雷達波的入射角、波長等無關。因此,當獲取同一目標k在n≥3個不同LoS向的位移時,式(2)可記為
2.1 區(qū)域概況
研究區(qū)域在3種SAR影像中的位置以及開采工作面、水準點的分布如圖2所示,工作面采用走向長壁全部垮落法開采。其中,工作面wf1的開采時間為2009-06-10—2009-12-31,走向長度510 m,傾角5°,煤層厚2.6 m,平均采深1 163 m;工作面wf2的開采時間為2009-08-15—2010-03-20,走向長度660 m,傾角3°,煤層厚2.5 m,平均采深1 140 m,其余工作面的開采時間在2009年之前。研究區(qū)域內有水準點43個,平均間距35 m。
圖2 研究區(qū)域位置以及工作面、水準點的分布Fig.2 Location of the study area and the layout of working face and leveling points
2.2 實驗數(shù)據(jù)
基于2種衛(wèi)星傳感器ASAR與PalSAR,筆者收集了3種軌道模式下的SAR影像,其成像時間與部分幾何參數(shù)見表1。由表1可以看出,ASAR數(shù)據(jù)組成的2個干涉對,其垂直基線與時間基線均較短,可以有效提高相干性以及測量精度;對于PalSAR數(shù)據(jù)組成的干涉對,其時間基線較短,雖然空間基線達580 m,但是相對于6 km的臨界基線仍可認為比較短。
表1 3種SAR影像的成像參數(shù)Table 1 Imaging paremeters of three types of SAR data
首先,表1中的3個干涉對分別利用DInSAR技術提取研究區(qū)域在3個LoS向的位移。因為3個干涉對的時間跨度并不完全重疊,需要通過插值獲得相同時間段(2010-01-13—2010-02-28)內的位移值。然后,根據(jù)式(4)計算研究區(qū)域內所有目標在上述時間段內的三維形變場,如圖3所示。在三維空間中,以向上、向東、向北為正方向。
圖3 垂直向、東西向與南北向位移場與等值線Fig.3 The displacement and contour in vertical,eastwest and north-south direction
圖3(a)為垂直向位移場,圖3(d)為與之對應的下沉等值線圖,可以看出,在2010-01-13—2010-02-28,研究區(qū)域下沉量值在0~60 mm,下沉盆地中心位于工作面wf1的上方,但是下沉盆地并沒有表現(xiàn)為類似工作面wf1的橢圓形狀,而是向兩側偏移,分析其原因可能是受到工作面wf1的影響導致其附近的老采空區(qū)活化,從而使得下沉盆地范圍向老采空區(qū)一側偏移,其中向北方向的偏移還有可能受到工作面wf2的開采影響。
由東西向位移場圖3(b)以及東西向位移等值線圖3(e)可以看出:研究區(qū)域的水平移動在-15~20 mm。在下沉盆地中心,水平移動量值比較小。自盆地中心向西與向東,水平移動量值均表現(xiàn)為先增大后減小的趨勢,并且在下沉盆地西側表現(xiàn)為向東移動,在東側則表現(xiàn)為向西移動,符合開采沉陷地表移動規(guī)律。
由圖3(c)與3(f)代表的南北向位移場與等值線圖可以看出,只在工作面wf1的西側與南側有5~10 mm的北向位移,其余大部分區(qū)域的南北向位移不明顯且較為雜亂,究其原因,由于現(xiàn)在的SAR衛(wèi)星均采用極軌飛行方式,飛行方向與北方向的夾角很小,式(2)中航向角α的正弦值近乎為0,因而使得本文算法對南北向位移不敏感,反演結果受噪聲影響嚴重。
3.1 與實測下沉值的對比分析
利用研究區(qū)域內43個水準點的觀測資料對本文采用的三維算法獲取的垂直向位移進行驗證分析,如圖4(a)所示??梢钥闯?兩者反映的下沉趨勢相互一致,并且在量值上也互相吻合,最大相差10 mm,絕對值相差在5 mm以內的點占77%,基于43個水準點統(tǒng)計兩者間的均方根誤差RMSE為±4 mm。由此可以說明本文所采用的模型與算法是有效的,同時也表明本文利用DInSAR技術提取的3個LoS向位移結果也是可靠的。為進一步對比分析,根據(jù)傳統(tǒng)算法忽略水平移動,采用式(5)計算垂直向位移。
式中,入射角θ與衛(wèi)星軌道模式以及目標所在位置有關,但由于本文研究區(qū)域較小使得因位置差異引起的θ的變化可以忽略不計,因此對同一軌道模式下的SAR影像數(shù)據(jù)采用相同的θ值。
根據(jù)式(5)將3個LoS向位移投影到垂直向,結果如圖4(b)所示。可以看出:傳統(tǒng)計算結果與實際下沉趨勢是一致的,但是根據(jù)PalSAR數(shù)據(jù)計算的下沉盆地向西發(fā)生了偏移,而由2個ASAR數(shù)據(jù)計算的下沉盆地形狀以及下沉量值均相似,并且與實際下沉的吻合程度優(yōu)于PalSAR數(shù)據(jù)結果。為定量比較2種算法,本文基于最大相差、絕對值相差在5 mm以內點的比例以及均方根誤差3項標準對2種方法進行了統(tǒng)計分析,結果見表2??梢钥闯?①定量分析的結果亦表明三維算法優(yōu)于傳統(tǒng)計算方法;②由式(5)可知,隨著雷達波入射角θ的增大,3組SAR數(shù)據(jù)反演結果與實際下沉的偏離應該越大。統(tǒng)計結果則證實了該結論的正確性。
圖4 垂直向位移與水準資料對比Fig.4 Comparison between the vertical displacements and leveling data
表2 不同算法得到的垂直向位移間的統(tǒng)計Table 2 Statistics of vertical displacements calculated using different methods
3.2 不同LoS向位移間的對比分析
圖5顯示的是43個水準點在3個LoS向的位移以及利用三維算法得到的東西向水平移動。從圖中可以看出:①利用2組ASAR數(shù)據(jù)提取的LoS向位移趨勢互相一致,量值也近乎相同,這是由該2種SAR數(shù)據(jù)的成像幾何參數(shù)相差較小所致。②ASAR與PalSAR兩種數(shù)據(jù)結果相差則較大,這是因為兩者的雷達波入射角與衛(wèi)星航向角相差均比較大所致。并且結合實際下沉曲線可以發(fā)現(xiàn),雷達視角越大,LoS向位移與實際值相差也越大。③由東西向水平移動曲線可以看出,在盆地中心,水平移動近乎為0,并且3個LoS向位移值均小于實測值。這點可以由式(2)證明,即在盆地中心,東西向與南北向位移均較小,使得LoS向位移僅是垂直向位移的一個分量。④在下沉盆地中心西側地表向東移動,并表現(xiàn)為先增大后減小的趨勢,在盆地中心東側地表向西移動,也表現(xiàn)為先增大后減小的趨勢,符合開采沉陷地表移動規(guī)律。并且由于本實驗選用的傳感器均采用右側視成像,使得在下沉盆地中心西側,東向的水平移動使得目標在ASAR影像上表現(xiàn)為靠近衛(wèi)星,在PalSAR影像上則表現(xiàn)為遠離衛(wèi)星;而在下沉盆地東側,西向的水平移動使得目標在ASAR影像上表現(xiàn)為遠離衛(wèi)星,在PALSAR影像上則表現(xiàn)為靠近衛(wèi)星。同時由圖5可以看出,在下沉盆地西側利用ASAR數(shù)據(jù)獲取的LoS向位移量值小于實測下沉值,利用PALSAR數(shù)據(jù)獲取的LoS向位移量值則大于實測下沉值,而在下沉盆地中心東側情況則相反。究其原因,當?shù)乇沓尸F(xiàn)下沉即Dv朝下時,由圖1中Dlos的矢量合成原理可知,靠近衛(wèi)星的東西向水平移動會減小Dlos的量值,而遠離衛(wèi)星的東西向水平移動則會使之增大。
圖5 三維算法得到的東西向位移以及不同LoS向的位移Fig.5 The displacement in east-west direction calculated using 3-D method and along different LoS
(1)提出了一種基于多源SAR影像聯(lián)合求解礦區(qū)地表三維變形的模型與算法。通過與實測水準資料的對比證明,三維算法獲取的垂直向位移曲線與水準測量結果的吻合程度優(yōu)于傳統(tǒng)算法,為利用多源SAR影像監(jiān)測礦區(qū)變形提供了理論基礎。
(2)由三維算法獲取的垂直向位移場與等值線可以發(fā)現(xiàn),下沉盆地向老采空區(qū)偏移,說明受工作面wf1的開采影響,附近的老采空區(qū)可能發(fā)生活化。
(3)利用三維算法獲取的東西向位移場與等值線符合開采沉陷地表移動規(guī)律,并且結合東西向位移對不同LoS向的位移進行對比分析,結果表明,當?shù)乇沓尸F(xiàn)為下沉時,朝向衛(wèi)星的東西向水平移動會減小LoS向的位移量值,而遠離衛(wèi)星的東西向水平移動則會使之增大。并且,雷達波入射角越大,對東西向水平移動的影響越大。
(4)因為衛(wèi)星采用極軌飛行方式,使得航向角α的正弦值近乎為0,本文提出的三維算法對南北向位移不敏感,反演精度低。歐洲空間局(ESA)與日本宇宙航空研究開發(fā)機構(JAXA)為本文提供了SAR數(shù)據(jù)(C1P.8371, C1P.10185),國家科學數(shù)據(jù)服務平臺為本文提供了SRTM-3數(shù)據(jù),在此一并致謝。
[1] Gabriel A K,Goldstein R M,Zebker H A.Mapping small elevation changes over large areas:differentialdar interferometry[J].Journal of Geophysical Research,1989,94(B7):9183-9191.
[2] Currenti Gilda,NapoliRosalba,Del Nero Ciro.Toward a realistic deformation model of the 2008 magmatic intrusion at Etna from combined DInSAR and GPS observations[J].Earth and Planetary Science Letters,2011,312(1-2):22-27.
[3] Cascini Leonardo,Fornaro Gianfranco,Peduto Dario.Advanced lowand full-resolution DInSAR map generation for slow-moving landslide analysis at different scale[J].Engineering Geology,2010,112(1-4):29-42.
[4] Shirzaei M,Buergmann R,Oncken O,et al.Response of forearc crustal to the megathrust earthquake cycle:InSAR evidence from Mejillones Peninsula,Northern Chile[J].Earth and Planetary Science letters,2012,333:157-164.
[5] Parks M M,Biggs J,Mather T A,et al.Co-eruptive subsidence at GalerasidentifiedduringanInSARsurveyofColombian Volcanoes(2006-2009)[J].Journal of Volcanology and Geothermal Research,2011,202(3-4):228-240.
[6] Hsieh C S,Shih T Y,Hu J C,et al.Using differential SAR interferometry to map land subsidence:A case study in the Pingtung Plain of SW Taiwan[J].Naural Hazards,2011,58(3):1311-1322.
[7] Ferretti A,Prati C,Rocca F.Non-linear subsidence rate estimation using permanent scatterers in differential SAR interferometry[J].IEEE Transactions on Geoscience and Remote Sensing,2000,38: 2202-2212.
[8] Ferretti A,Prati C,Rocca F.Permanent scatterers in SAR interferometry[J].IEEE Transactions on Geoscience and Remote Sensing, 2001,39:8-30.
[9] Mora O,Mallorqui J J,Duro J,et al.Longterm subsidence monitoring of urban areas using differential interferometric SAR technique[A].Preoceedings of IEEE Geoscience and Remote Sensing Symposium 2001[C].Sydney:IEEE,2001:1104-1106.
[10] Berardino P,Fornaro G,Lanari R,et al.A new algorithm for surface deformation monitoring based on small baseline differential interferograms[J].IEEE Transactions on Geoscience and Remote Sensing,2002,40:2375-2383.
[11] Samsonov S,Beavan J,Gonzalez,P J,et al.Ground deformation in the Taupo Volcanic Zone,New Zealand,observed by ALOS PALSAR interferometry[J].Geophysical Journal International,2011,187(1):147-160.
[12] 張永紅,吳宏安,孫光通.時間序列InSAR技術中的形變模型研究[J].測繪學報,2012,41(6):864-869.
Zhang Yonghong,Wu Hongan,Sun Guangtong.Deformation model of time series interferometric SAR techniques[J].Acta Geodaetica et Cartographica Sinica,2012,41(6):864-869.
[13] Via G D,Crosetto M,Crippa B.Resolving vertical and east-west horizontal motion from differential interferometric synthetic aperture radar:the L’Aquila earthquake[J].Journal of Geophysical Research,2012,117(B02310):1-14.
[14] Hanssen R F.Radar interferometry data interpreation and error analysis[M].New York:Springer,2001.
[15] 范洪冬,鄧喀中,祝傳廣,等.基于時序SAR技術的采空區(qū)上方高速公路變形監(jiān)測及預測方法[J].煤炭學報,2012,37(11): 1841-1846.
Fan Hongdong,Deng Kazhong,Zhu Chuanguang,et al.Study on deformation monitoring and predicting methods for expressway above goaf based on time series SAR technique[J].Journal of China Coal Society,2012,37(11):1841-1846.
[16] Fialko Y,Simons M.The complete(3-D)surface displacement field in the epicentral area of the 1999Mw 7.1 Hector Mine earthquake, California,from space geodetic observations[J].Geophysical Research Letters,2001,28(16):3063-3066.
[17] Aly M H,Cochran E S.Spatio-temporal evolution of yellowstone deformation between 1992 and 2009 from InSAR and GPS observations[J].Bulletin of Volcanology,2009,73(9):1407-1419.
Three-dimensional deformation field detection based on multi-source SAR imagery in mining area
ZHU Chuan-guang1,2,DENG Ka-zhong1,ZHANG Ji-xian2, ZHANG Yong-hong2,FAN Hong-dong1,ZHANG Li-ya3
(1.Key Laboratory for Land Environment and Disaster Monitoring of SBSM,China University of Mining and Technology,Xuzhou 221116,China;2.Chinese Academy of Surveying and Mapping,Beijing 100036,China;3.Hunan Province Key Laboratory of Coal Resources Clean-utilization and Mine Environment, Xiangtan 411100,China)
Only one dimensional displacement in the line of sight(LoS)can be detected using DInSAR technique.In order to overcome this limitation,the three-dimensional(3-D)model and algorithm which fusing multi-source SAR images acquired with different geometric parameters in different satellite platform was proposed and validated to determine the 3-D displacement due to underground mining.Comparing with leveling observation indicates that the displacement derives from the 3-D algorithm coincide with leveling data more closely,and the root mean square error (RMSE)is±4 mm.May be due to the activation of old goafs,both the displacement and contour in vertical direction indicates that the subsidence basin extends to the old goaf.The displacement and contour in east-west direction accord with the law of surface movement caused by mining,and the effects on the displacement along LoS will be significantly different when imaging mode different.Due to the small sine value of heading angle,the 3-D algorithm is not sensitive to the displacement in the north-south direction.
P237;TD173
A
0253-9993(2014)04-0673-06
祝傳廣,鄧喀中,張繼賢,等.基于多源SAR影像礦區(qū)三維形變場的監(jiān)測[J].煤炭學報,2014,39(4):673-678.
10.13225/j.cnki.jccs.2013.0424
Zhu Chuanguang,Deng Kazhong,Zhang Jixian,et al.Three-dimensional deformation field detection based on multi-source SAR imagery in mining area[J].Journal of China Coal Society,2014,39(4):673-678.doi:10.13225/j.cnki.jccs.2013.0424
2013-04-07 責任編輯:王婉潔
江蘇高校優(yōu)勢學科建設工程資助項目(SZBF2011-6-B35);國家自然科學基金資助項目(41272389,41071273)
祝傳廣(1984—),男,山東成武人,博士研究生。E-mail:zhuchuanguang@163.com。通訊作者:鄧喀中(1957—),男,四川資中人,教授,博士生導師。E-mail:kzdeng@cumt.edu.cn