基于卡爾曼濾波的D－InSAR和水準監測數據融合方法研究

李懷展 查劍鋒 米麗倩

1 中國礦業大學環境與測繪學院，徐州市大學路1號，221008

2 山東省國土測繪院，濟南市經十路2301號，250102

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

1 高等級公路多源數據對比分析

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

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

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

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

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

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

2 數據融合模型的建立

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

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

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

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

3 融合實驗與結果分析

3.1 實驗數據

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

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

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

圖3 監測點分布Fig.3 Location of bench mark

3.2 融合結果分析

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

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

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

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

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

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

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

表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監測結果與水準觀測值走勢十分接近，但兩種數據在數值上呈現出一定的差異，距離盆地中心越近，D－InSAR 與水準實測值的差別越大。基于集合卡爾曼濾波同化D－InSAR 監測結果與水準觀測值，實現高精度面狀監測的方法是可行的，同化結果相比于反演值和D－InSAR 監測結果有了很大改善。將該方法應用于濟寧鄒濟公路，總均方根誤差為17.7 mm，滿足高等級公路變形監測的精度要求。

［1］吳立新，高均海，葛大慶，等.工礦區地表沉陷D－InSAR監測試驗研究［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 的煤礦區開采沉陷遙感監測技術分析［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技術在大范圍滑坡監測中的應用［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］劉國祥，丁曉利.沿海地區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］費劍鋒，韓月琪.集合卡爾曼濾波數據同化在一維波動方程中的應用［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］張學峰，黃大吉，章本照，等.集合數據同化方法的發展與應用概述［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濾波在氣象數據同化中的發展 與 應 用［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 二軌法監測徐州大屯中心區 地表形 變［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 技術二軌法監測礦區地表沉陷的方法研究［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］王躍山.數據同化——它的緣起、含義和主要方法［J］.海洋預報，1999（1）：12－21（Wang Yueshan.Data Assimilation－Its Cause，Its Meaning and Main Procedures［J］.Marine Forecasts，1999（1）：12－21）