李懷展 查劍鋒 米麗倩

1 中國礦業(yè)大學環(huán)境與測繪學院，徐州市大學路1號，221008

2 山東省國土測繪院，濟南市經(jīng)十路2301號，250102

D－InSAR 技 術是 沉 降 監(jiān) 測 的 重 要 手 段［1－4］，其監(jiān)測結果受多種誤差源影響，沉降最大值附近精度較低；而水準測量雖然精度較高，但只能得到有限個監(jiān)測點的變形值。兩種數(shù)據(jù)具有很好的互補性，可嘗試應用數(shù)據(jù)同化的方法將D－InSAR 和水準監(jiān)測數(shù)據(jù)融合。數(shù)據(jù)同化是指在考慮數(shù)據(jù)時空分布以及觀測場和背景場誤差的基礎上，在數(shù)值模型的動態(tài)運行過程中融合新的觀測方法，對陸地、大氣、海洋的動態(tài)模型軌跡進行一定的約束，從而提高變量的估計精度。作為數(shù)據(jù)同化的一種新方法，集合卡爾曼濾波近10a來取得很大進展［4－8］。本文以高等級公路高精度沉降變形監(jiān)測為例，提出一種基于集合卡爾曼濾波同化的DInSAR－水準數(shù)據(jù)融合的思路。

1 高等級公路多源數(shù)據(jù)對比分析

鄒濟公路在南屯煤礦井田內(nèi)長約8 400 m，煤礦開采對其有顯著影響。為確保公路的安全運行，必須開展采動區(qū)高等級公路變形監(jiān)測及預報。對鄒濟公路進行了水準測量和InSAR 監(jiān)測，干涉雷達影像對及參數(shù)如表1。采用“二軌法”處理雷達影像［9－10］，并與實測水準數(shù)據(jù)進行對比。各期D－InSAR與水準累加結果對比如圖1。

從圖1可以看出，疊加后的D－InSAR 監(jiān)測結果與水準觀測值走勢上十分接近，能夠從整體上反映公路沉降變化的趨勢。但兩種數(shù)據(jù)在數(shù)值上呈現(xiàn)出一定的差異，D－InSAR 監(jiān)測值在整體上要小于水準實測值，且兩種數(shù)據(jù)呈現(xiàn)出“中間差別大、兩邊差別小”的現(xiàn)象，即距離盆地中心越近，DInSAR 與水準實測值的差別越大。這主要同雷達影像的空間分辨率和波長有關。D－InSAR 最大形變監(jiān)測梯度公式為：

表1 干涉雷達影像對及參數(shù)Tab.1 Interferometric SAR images and their parameters

其中，λ為波長，ps為像元大小。對于波長為3 cm 的TerraSAR－X，經(jīng)多視后像元大小為9m，干涉像元可探測到的最大形變梯度為1.6 mm／m。兩次觀測間隔內(nèi)地表的下沉差越大，形變的梯度則越大，當形變梯度超過最大可探測梯度時，這些地表沉降值差別過大的點便無法準確地在地表形變圖上體現(xiàn)出來，使得沉降越大的區(qū)域D－InSAR的誤差越大。

圖1 D－InSAR與水準測量數(shù)據(jù)對比Fig.1 Comparison of D－InSAR and leveling data

綜合以上分析可知，D－InSAR 監(jiān)測能周期性地獲取監(jiān)測區(qū)域的面元變形信息，可達到cm級精度，但與常規(guī)水準監(jiān)測相比精度不高；而水準測量數(shù)據(jù)雖然精度較高，能達到mm級精度，但只能得到點元變形信息。所以，在實際應用中，應探索把這兩種數(shù)據(jù)進行有效融合，達到數(shù)據(jù)的高空間分辨率與高精度的統(tǒng)一。

2 數(shù)據(jù)融合模型的建立

同化的本意是把不同的事物變得相近或相同［11］。設A為模型模擬結果，B為真實觀測值，則同化可有以下5層含義：第一，A不變，B變得和A相近或相同；第二，A基本不變，B向A同化；第三，A和B以相互靠攏的方式變得相近或相同；第四，B基本不變，A向B同化；第五，B不變，A變得和B相近或相同。如果以多源數(shù)據(jù)的融合集成利用為同化目標，則應該采用第三種同化方案。

考慮到基于點觀測的水準實測數(shù)據(jù)比雷達影像的處理結果更為可靠、數(shù)值也更為準確，利用集合卡爾曼濾波進行衛(wèi)星影像監(jiān)測與水準實測數(shù)據(jù)的融合時，采用以水準實測數(shù)據(jù)（觀測）為主的同化方式。

由雷達影像處理得到的D－InSAR 監(jiān)測值反映的是區(qū)域面元信息，而路面水準實測值獲取的是有限個監(jiān)測點元的變形，要將兩種數(shù)據(jù)進行融合，必須對數(shù)據(jù)進行有效的插值，實現(xiàn)二者的空間匹配。先對水準數(shù)據(jù)的監(jiān)測點結果進行插值，得到沉降實測估計值在整個區(qū)域的變形，再與區(qū)域的D－InSAR數(shù)據(jù)進行融合計算。

同化系統(tǒng)如圖2所示，數(shù)據(jù)包括雷達影像數(shù)據(jù)、水準實測數(shù)據(jù)，同化方法為Burgers等提出的擾動觀測集合卡爾曼濾波算法［12］，數(shù)據(jù)和同化方法組成了整個同化系統(tǒng)。

3 融合實驗與結果分析

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

本文選取9308 工作面為研究對象，選取2012－01－05～04－02的D－InSAR 累積沉降結果進行實驗。由于9308工作面水準監(jiān)測從2012－03－02開始，無法給出公路監(jiān)測點自01－05的水準累積沉降結果，因此將2012－01－05～02－29的D－In－SAR 監(jiān)測數(shù)據(jù)作為初始累積沉降，與后期的水準測量結果進行累加，得到03－11～04－02的水準累加值（水準值）。

首先根據(jù)水準值，利用遺傳算法反演出概率積分法預計參數(shù)，并在此參數(shù)基礎上，采用概率積分法預計出工作面的區(qū)域沉降值（反演值），作為觀測場引入到同化系統(tǒng)中，通過集合卡爾曼濾波同化算法融合D－InSAR值與反演值，并將同化結果與真實值相比較，從而判斷數(shù)據(jù)同化方法融合數(shù)據(jù)的準確性和可靠性，將D－InSAR值與反演值經(jīng)集合卡爾曼濾波計算得出融合結果。

受9308工作面采動影響的鄒濟公路監(jiān)測點共24個，隨機選取其中11個作為驗證點，其余為實驗點，如圖3所示，其中“Δ”表示實驗點，“○”表示驗證點。

圖3 監(jiān)測點分布Fig.3 Location of bench mark

3.2 融合結果分析

將數(shù)據(jù)導入同化系統(tǒng)，利用集合卡爾曼濾波進行整個面元的數(shù)據(jù)融合計算。將融合結果內(nèi)插后得出公路上各觀測站點的同化值，并與水準值、反演值以及D－InSAR值進行對比，如圖4所示。

圖4 驗證點數(shù)據(jù)結果對比Fig.4 The comparison diagram of verification point’s data

從圖4可以看出，4種數(shù)據(jù)的變化趨勢總體上一致，但是數(shù)值大小呈現(xiàn)出一定的差異?？傮w而言，反演值與水準值的中心區(qū)域差別較小，兩邊差別較大；D－InSAR 值與水準值的中心區(qū)域差別較大，兩邊差別較?。煌祫t處于兩者之間，這與理論分析相吻合。實際上，圖4的結果驗證了同化的含義，即作為分析場的同化值是將背景場與觀測場向中間靠近，其數(shù)值應該介于作為背景場的D－InSAR值與作為觀測場的反演值之間。地表沉降的同化結果雖然與水準值仍存在一定差距，但相比于反演值和D－InSAR結果有了很大改善。

為了更準確地衡量系統(tǒng)的同化效果，在RMSE、r以及MAPE（平均絕對百分比誤差）3個指標下比較了反演值、D－InSAR 值和同化值。其中，

式中，Yi表示第i個測點的水準值，Xi為第i個測點的其他測量數(shù)據(jù)，N為測點數(shù)。

由表2 可得，在RMSE、r、MAPE 三個指標下，大部分同化值結果要比D－InSAR 與反演值好，同化值的3個指標各自有0、1、1個觀測時刻的數(shù)據(jù)效果位于D－InSAR 與反演值之間。各個時刻同化值的均方根誤差均小于22mm，滿足高等級公路變形監(jiān)測的精度要求。針對表2的驗證數(shù)據(jù)，圖5進一步從散點圖和誤差分布直方圖對同化結果與反演值以及D－InSAR值進行比較。

綜合表2與圖5可得，與水準值相比，同化值的效果最好，反演值次之，D－InSAR 值最差。目前高等級公路沉降觀測沒有相應規(guī)范。現(xiàn)場調(diào)研情況顯示，公路下沉200mm 以內(nèi)，采動變形對公路影響輕微，公路上產(chǎn)生不大于2mm 的裂縫，路面基本無隆起現(xiàn)象，行車無不適感，不用維修。根據(jù)國際測量工作者聯(lián)合會第13屆會議提出的精度要求，變形觀測中誤差均方根誤差應小于允許變形值的1／10～1／20，高等級公路采動變形破壞測量精度應達20 mm。同化值的總均方根誤差為17.7 mm，滿足高等級公路變形監(jiān)測的精度要求。

表2 公路驗證點精度分析Tab.2 Precision analysis table of verification poing’s data

圖5 同化值與反演值及D－InSAR值的比較Fig.5 Comparison of the assimilation results and inversion，D－InSAR

4 結 語

D－InSAR監(jiān)測結果與水準觀測值走勢十分接近，但兩種數(shù)據(jù)在數(shù)值上呈現(xiàn)出一定的差異，距離盆地中心越近，D－InSAR 與水準實測值的差別越大?；诩峡柭鼮V波同化D－InSAR 監(jiān)測結果與水準觀測值，實現(xiàn)高精度面狀監(jiān)測的方法是可行的，同化結果相比于反演值和D－InSAR 監(jiān)測結果有了很大改善。將該方法應用于濟寧鄒濟公路，總均方根誤差為17.7 mm，滿足高等級公路變形監(jiān)測的精度要求。

［1］吳立新，高均海，葛大慶，等.工礦區(qū)地表沉陷D－InSAR監(jiān)測試驗研究［J］.東北大學學報，2005（8）：778－782（Wu Lixin，Gao Junhai，Ge Daqing，et al.Experimental Study on Surfaee Subsidenee Monitoring with D－InSAR in Mining Area［J］.Journal of Northestern University，2005（8）：778－782）

［2］吳立新，高均海，葛大慶，等.基于D－InSAR 的煤礦區(qū)開采沉陷遙感監(jiān)測技術分析［J］.地理與地理信息科學，2004（2）：22－25（Wu Lixin，Gao Junhai，Ge Daqing，et al.Technical Analysis of the Remote Sensing Monitoring for Coal－Mining Subsidence Based on D－InSAR［J］.Geography and Geo－Information Science，2004（2）：22－25）

［3］王桂杰，謝謨文，邱騁，等.D－InSAR技術在大范圍滑坡監(jiān)測中的應用［J］.巖土力學，2010（4）：1 337－1 344（Wang Guijie，Xie Mowen，Qiu Cheng，et al.Application of D－InSAR Technique to Landslide Monitoring［J］.Rock and Soil Mechanics，2010（4）：1 337－1 344）

［4］劉國祥，丁曉利.沿海地區(qū)D－InSAR 形變探測：精度與可應用性分 析［J］.測 繪 通 報，2006（9）：9－13（Liu Guoxiang，Ding Xiaoli.Detecting of Surface Deformation in Coastlands by D－InSAR：Accuracy and Applicability Analysis［J］.Bulletin of Surveying and Mapping，2006（9）：9－13）

［5］費劍鋒，韓月琪.集合卡爾曼濾波數(shù)據(jù)同化在一維波動方程中的應用［J］.氣象科技，2005（2）：109－114（Fei Jianfeng，Han Yueqi.Application of Ensemble Kalman Filter Assimilation in One－Dimensional Wave Equation［J］.Meteorological Science and Technology，2005（2）：109－114）

［6］張學峰，黃大吉，章本照，等.集合數(shù)據(jù)同化方法的發(fā)展與應用概述［J］.海洋學研究，2007（1）：88－94（Zhang Xuefeng，Huang Daji，Zhang Benzhao，et al.The Developments and Applications of Ensemble－Based Data Assimilation Methods［J］.Journal of Marine Science，2007（1）：88－94）

［7］陳春，吳振森，孫樹計，等.集合卡爾曼濾波在電離層短期預報中的應用［J］.空間科學學報，2010（2）：148－153（Chen Chun，Wu Zhensen，Sun Shuji.Application of the Ensemble Kalman Filter in Short－Term Ionospheric Forecast［J］.Chinese Journal of Space Science，2010（02）：148－153）

［8］高山紅，吳增茂，謝紅琴.Kalman濾波在氣象數(shù)據(jù)同化中的發(fā)展 與 應 用［J］.地 球 科 學 進 展，2000（5）：571－575（Gao Shanhong，Wu Zengmao，Xie Hongqin.The Development and Applications of Kalman Filters in Meteorological Data Assimilation［J］.Advance in Earth Sciences，2000（5）：571－575）

［9］朱江，汪萍.集合卡爾曼平滑和集合卡爾曼濾波在污染源反演中的應用［J］.大氣科學，2006（5）：871－882（Zhu Jiang，Wang Ping.Ensemble Kalman Smoother and Ensemble Kalman Filter Approaches to the Joint Air Quality State and Emission Estimation Problem［J］.Chines Journal of Atmospheric Sciences，2006（5）：871－882）

［10］何敏，陸曉燕，何秀鳳.利用D－InSAR 二軌法監(jiān)測徐州大屯中心區(qū) 地表形 變［J］.地 理 空 間 信 息，2011（5）：3－5（He Min，Lu Xiaoyan，He Xiufeng.Monitoring Ground Deformation in Datun Central District of Xuzhou Using Two Pass D－InSAR［J］.Geospatial Information，2011（5）：3－5）

［11］魏長婧，汪云甲，閆建偉.D－InSAR 技術二軌法監(jiān)測礦區(qū)地表沉陷的方法研究［J］.煤炭技術，2012（7）：121－122（Wei Changjing，Wang Yunjia，Yan Jianwei.Study of Application of D－InSAR Technology in Monitoring Mining Area Surface Settle by Two－Track Method［J］.Coal Technology，2012（7）：121－122）

［12］王躍山.數(shù)據(jù)同化——它的緣起、含義和主要方法［J］.海洋預報，1999（1）：12－21（Wang Yueshan.Data Assimilation－Its Cause，Its Meaning and Main Procedures［J］.Marine Forecasts，1999（1）：12－21）