一種基于灰色預(yù)測理論和抗差自適應(yīng)Kalman濾波的滑坡監(jiān)測算法

楊 旭，楊 旭，李 佳，王建國

(1.北京航空航天大學(xué)電子信息工程學(xué)院，北京 100089；2.北京航空航天大學(xué)合肥創(chuàng)新研究院，合肥 230000)

0 引言

山體滑坡災(zāi)害具有隱蔽性強、突發(fā)性強、破壞力大的特點。傳統(tǒng)的山體滑坡監(jiān)測方法主要有宏觀地質(zhì)觀測法、簡易觀測法、設(shè)站觀測法、儀表觀測法等，它們普遍存在的問題是需要人工定期到現(xiàn)場進行數(shù)據(jù)的采集，缺少實時性，監(jiān)測精度較低，并且對滑坡的預(yù)測預(yù)警能力差[1-3]。

近年來，在地質(zhì)形變監(jiān)測技術(shù)不斷發(fā)展和完善的過程中，大地精密測量法、全球衛(wèi)星導(dǎo)航系統(tǒng)(Global Navigation Satellite System，GNSS)高精度地質(zhì)形變監(jiān)測法、近景攝影測量法、時間域反射測試技術(shù)、遙感監(jiān)測技術(shù)等方法常被用于地質(zhì)形變監(jiān)測中。目前主流的山體滑坡監(jiān)測方法是通過合成孔徑雷達(dá)(Synthetic Aperture Radar，SAR)技術(shù)進行監(jiān)測，但是該方法受天氣狀況和植被覆蓋的影響大，實時性差，難以實現(xiàn)長期監(jiān)測[4-6]。

實時動態(tài)(Real-Time Kinematic，RTK)技術(shù)可實現(xiàn)實時定位，并且在靜態(tài)情況下的定位精度可達(dá)mm級。與傳統(tǒng)GNSS形變監(jiān)測技術(shù)相比，RTK定位技術(shù)具有十分顯著的優(yōu)勢。RTK技術(shù)被廣泛應(yīng)用于安全監(jiān)測、土地勘探、城市布局、水利等領(lǐng)域[7-11]。在我國的金川二礦區(qū)地裂縫[9]、黃土坡臨江I號崩滑體[10]等地質(zhì)災(zāi)害的變形監(jiān)測中，均采用了高精度的全球定位系統(tǒng)(Global Positioning System，GPS)靜態(tài)相對定位技術(shù)，監(jiān)測精度可以達(dá)到mm級；并逐漸開始嘗試采用GPS定位技術(shù)監(jiān)測滑坡體的垂直形變，發(fā)現(xiàn)在基線較短、觀測時間較長和觀測條件較好的情況下，其精度也可以達(dá)到mm級[11]。

作為形變監(jiān)測中的重要數(shù)據(jù)處理算法，Kalman濾波算法廣泛地應(yīng)用于建筑、橋梁、地質(zhì)災(zāi)害的形變監(jiān)測中，并取得了非常不錯的監(jiān)測效果[12-15]。C.D.Ince等通過GPS定位技術(shù)和Kalman濾波方法對土耳其的伊斯坦布爾和馬爾馬拉進行形變監(jiān)測，并取得預(yù)期的監(jiān)測效果[12]；朱建軍等使用變形監(jiān)測數(shù)據(jù)通過Kalman濾波模型動態(tài)分析鏈子崖危巖體位移變形規(guī)律[13]；高雅萍、張勤使用自適應(yīng)Kalman濾波在地理信息系統(tǒng)(Geographic Information System，GIS)滑坡監(jiān)測中也取得了預(yù)期監(jiān)測效果[14]。

影響山體滑坡的因素復(fù)雜多變，有已知因素和未知因素，具有典型的灰色特征。灰色預(yù)測理論在解決信息不完全的系統(tǒng)問題時具有十分顯著的優(yōu)勢。對于具有這樣特征的系統(tǒng)，可使用灰色預(yù)測理論對其進行深入的分析[16-19]。李曉鴿采用灰色預(yù)測理論對魯?shù)乩娬緦ν饨煌ü坊卤O(jiān)測數(shù)據(jù)進行了擬合和預(yù)測研究，得到了很好的預(yù)測效果[16]；王旭昭等使用改進的灰色預(yù)測理論對新灘滑坡和黃茨滑坡進行預(yù)測分析，取得了良好預(yù)測效果[17]；此外，K. Yin[18]、王東岳[19]等均使用灰色預(yù)測理論進行了山體滑坡預(yù)測，并取得良好預(yù)測效果。

針對當(dāng)前的山體滑坡監(jiān)測技術(shù)精度低、實時性差、缺少完善的預(yù)測預(yù)警機制的問題，本文提出了一種基于抗差自適應(yīng)Kalman濾波和灰色預(yù)測理論的山體滑坡監(jiān)測技術(shù)。該技術(shù)通過抗差自適應(yīng)Kalman濾波進行包括RTK定位技術(shù)、土工帶傳感器技術(shù)、無人機近景攝影測量技術(shù)在內(nèi)的多源數(shù)據(jù)融合，然后通過灰色預(yù)測理論與蠕變切線角判據(jù)實現(xiàn)對山體滑坡的預(yù)測預(yù)警功能。

1 抗差自適應(yīng)Kalman濾波

1.1 抗差自適應(yīng)Kalman濾波[12-15]

在目標(biāo)系統(tǒng)中，設(shè)Xi為由一組隨機變量組成的向量，稱為該系統(tǒng)在i時刻的狀態(tài)向量；設(shè)Zi為由一組觀測變量組成的向量，稱為該系統(tǒng)在i時刻的觀測向量。當(dāng)目標(biāo)系統(tǒng)為離散時間系統(tǒng)時，系統(tǒng)的狀態(tài)向量隨時間的變化規(guī)律可以通過狀態(tài)方程描述，目標(biāo)系統(tǒng)的狀態(tài)向量和觀測向量的關(guān)系可以通過觀測方程描述

(1)

其中，Xi-1、Xi分別是目標(biāo)系統(tǒng)在i-1和i時刻的狀態(tài)向量；Zi是目標(biāo)系統(tǒng)在i時刻的觀測向量；U稱為控制向量；W稱為系統(tǒng)動態(tài)噪聲向量；V稱為系統(tǒng)觀測噪聲向量；A稱為系統(tǒng)從i-1時刻到i時刻的狀態(tài)轉(zhuǎn)移矩陣；B稱為控制向量系統(tǒng)矩陣；H稱為觀測值系統(tǒng)矩陣。

抗差自適應(yīng)Kalman濾波技術(shù)將抗差估計和自適應(yīng)估計相結(jié)合，以減小觀測誤差和系統(tǒng)建模誤差對濾波結(jié)果的影響，其技術(shù)原理可以通過如下方程加以描述：

時間更新方程

(2)

Pi/i-1=APi-1/i-1AT+Q

(3)

狀態(tài)更新方程

(4)

(5)

(6)

1.2 觀測噪聲等價協(xié)方差矩陣的構(gòu)造

典型的觀測噪聲等價協(xié)方差矩陣D的構(gòu)造方法[14]有：Huber構(gòu)造法、丹麥法、IGG法、抗差貝葉斯估計法等。本文選用IGG1構(gòu)造法，具體內(nèi)容如下

(7)

其中，一般取k0=2.5，k1=3.5。Vi為觀測殘差，表達(dá)式如下

(8)

1.3 自適應(yīng)因子的構(gòu)造

本文采用兩段函數(shù)模型[14]進行自適應(yīng)因子的構(gòu)造

(9)

式中，c為常量，一般取c=1～2.5，0

2 灰色預(yù)測理論

2.1 灰色預(yù)測模型(GM模型)

2.1.1 基本GM(1,1)模型

設(shè)X(0)為含有n個元素的序列，其中x(0)(k)≥0,k=1,2,…,n。上角標(biāo)處括號內(nèi)的數(shù)字為序列模型的階數(shù)，第二個括號內(nèi)的數(shù)字為序列中元素的序號。對序列X(0)進行一次累加生成得到它的1-AGO序列X(1)。對序列X(1)進行一次緊鄰均值生成得到序列Z(1)。稱方程x(0)(k)+ax(1)(k)=b為GM(1,1)模型。其中a為發(fā)展系數(shù)，b為灰色作用量。令

(10)

Y=[x(0)(2),x(0)(3),…,x(0)(n)]T

(11)

(12)

則有

(13)

方程x(0)(k)+ax(1)(k)=b的解

(14)

取x(1)(0)=x(0)(1)，則

(15)

2.1.2 新陳代謝GM(1,1)模型

2.2 精度檢驗

(16)

(17)

在實際應(yīng)用中，根據(jù)實際情況設(shè)定閾值A(chǔ)0，當(dāng)A

(18)

在實際應(yīng)用中，根據(jù)實際情況設(shè)定閾值P0，當(dāng)P>P0時，則認(rèn)為該模型為小概率誤差合格模型。關(guān)于均值方差比閾值A(chǔ)0和小概率誤差閾值P0的設(shè)定準(zhǔn)則與模型精度等級的關(guān)系如表1所示。

表1 模型精度等級參考[19]

3 蠕變切線角判據(jù)

通過大量的山體滑坡實例可知，山體滑坡的形變過程可分為3個時期：初始變化時期、穩(wěn)定變化時期和劇烈變化時期。劇烈變化時期進一步細(xì)分為：初加速時期、中加速時期和臨破壞時期[19,21-22]。如圖1所示。

圖1 山體滑坡演化過程中的形變曲線

本文基于上述分析建立蠕變切線角預(yù)警判據(jù)[20]。通過一定時間的觀測得到穩(wěn)定變化時期監(jiān)測點的形變速度v0，將穩(wěn)定變化時期的山體形變運動看作是勻速運動，可知v0是定值，定義

(19)

其中，S(k)為從起始監(jiān)測時刻到當(dāng)前監(jiān)測時刻監(jiān)測點處的累計位移量；T(k)為具有時間量綱的縱坐標(biāo)值。由T(k)和t(k)繪制T-t曲線，如圖2所示，通過T-t曲線可以得到切線角的表達(dá)式

圖2 經(jīng)過處理后得到的T-t曲線

(20)

其中，αk為切線角；t(k)為第k個監(jiān)測時刻；Δt為采樣時間間隔。

根據(jù)行業(yè)內(nèi)的長期觀測研究得到山體滑坡過程的預(yù)警判據(jù)表，如表2所示。該預(yù)警判據(jù)表綜合考慮了山體滑坡監(jiān)測點的形變速率V和切線角α包含的信息內(nèi)容，對山體滑坡預(yù)警等級進行了5層劃分，分別代表不同的山體滑坡可能性和危險性，從而實現(xiàn)對山體滑坡的等級預(yù)警功能。

課文大意：《綜合英語》第四冊中第七單元the Monster是一篇典型的描寫文，栩栩如生地刻畫了十九世紀(jì)德國劇作家瓦格納的人物形象。作為劇作家，瓦格納是一個天才，畢生創(chuàng)作十三部歌劇和音樂劇，其中十一部仍繼續(xù)上演，八部被列入世界音樂名劇之列。但現(xiàn)實中，瓦格納可以稱得上是一個怪物，他病態(tài)，神經(jīng)質(zhì)，喜歡妄自尊大，絲毫不具責(zé)任感，喜新厭舊，甚至偷人妻子，賴賬不還。(張伊娜,216)這就是本文的主人公，一個集所有缺點于一身的怪物，同時他又是音樂界的天才，具備超凡的稟賦。作者這種先抑后揚的寫法極好地刻畫了這個杰出的人物。

表2 山體滑坡預(yù)警判據(jù)表[20]

4 實驗測試

4.1 滑坡監(jiān)測系統(tǒng)

本文構(gòu)建以RTK為核心的空天地一體化立體監(jiān)測網(wǎng)絡(luò)，實現(xiàn)對滑坡災(zāi)害的監(jiān)測評估。

滑坡監(jiān)測系統(tǒng)如圖3所示，使用RTK進行地質(zhì)沉降監(jiān)測，在監(jiān)測區(qū)域布置土工帶傳感器，采集監(jiān)測區(qū)域的位移信息并上傳到云端以待處理，利用無人機攝影測量技術(shù)提供空域增強信息，輔助RTK進行滑坡形變監(jiān)測。北斗短報文通信技術(shù)可以實現(xiàn)全天候、全覆蓋域、高可靠性的通信，其通信時延約為0.5s，點對點通信時延約為1～5s，同時該技術(shù)的設(shè)備成本低，適用于廣泛的環(huán)境監(jiān)測活動，所以RTK形變監(jiān)測信息通過北斗短報文進行傳遞。

圖3 山體滑坡監(jiān)測系統(tǒng)

設(shè)計后處理軌跡平滑算法剔除異常測量信息，使用數(shù)據(jù)融合[23]技術(shù)對RTK定位數(shù)據(jù)、土工帶傳感器數(shù)據(jù)、無人機三維成像數(shù)據(jù)進行綜合分析后，在滑坡監(jiān)測點處計算得到位移、速度、加速度測量信息。本文選用抗差自適應(yīng)Kalman濾波技術(shù)實現(xiàn)數(shù)據(jù)的融合分析，獲得精確的形變監(jiān)測結(jié)果。在獲得經(jīng)過融合處理的形變監(jiān)測數(shù)據(jù)后，通過灰色預(yù)測理論對未來時刻的形變數(shù)據(jù)進行合理預(yù)測，得到未來時刻監(jiān)測點處的形變信息，最后結(jié)合蠕變切線角判據(jù)實現(xiàn)對山體滑坡災(zāi)害的預(yù)測預(yù)警功能。

設(shè)備布設(shè)的具體情況如圖4所示。在重點監(jiān)測區(qū)域進行監(jiān)測設(shè)備布設(shè)時，為合理兼顧設(shè)備成本和監(jiān)測效果，主要將監(jiān)測設(shè)備布設(shè)在可能滑坡體的頂部和底部，尤其是一些裂縫處。考慮到皖南山區(qū)植被覆蓋程度大、連續(xù)性強，植被覆蓋斷裂痕跡處也需重點布設(shè)監(jiān)測設(shè)備。RTK接收機安置在網(wǎng)絡(luò)的中心位置。該土工帶傳感器與RTK接收機的排布網(wǎng)絡(luò)與無人機攝影測量得到的形變監(jiān)測數(shù)據(jù)相結(jié)合，三種形變監(jiān)測數(shù)據(jù)在空間上具備了數(shù)據(jù)融合的條件。

圖4 監(jiān)測區(qū)域數(shù)據(jù)采集設(shè)備布設(shè)

4.2 滑坡監(jiān)測模型

對山體滑坡進行監(jiān)測分析時，監(jiān)測點的形變可以使用運動模型進行分析。運動模型的描述方程如下

(21)

(22)

位移x(i)通過土工帶傳感器直接測得

x(i)=dLtu(i)

(23)

(24)

加速度a(i)通過定期無人機攝影測量形變量和飛行間隔時間計算得到

(25)

為保證監(jiān)測系統(tǒng)的實時性，在平穩(wěn)形變階段，監(jiān)測區(qū)域加速度信息基本穩(wěn)定，無人機攝影測量周期TUAV可設(shè)置較大，TUAV≈30d≈15T。在加速形變易發(fā)生滑坡階段，應(yīng)按照實際情況進行調(diào)整，對于監(jiān)測價值極高的區(qū)域，無人機攝影測量周期TUAV需大幅縮短；對于監(jiān)測價值一般的區(qū)域，可屏蔽無人機攝影測量數(shù)據(jù)，僅僅使用土工帶傳感器測量數(shù)據(jù)和RTK接收機測量數(shù)據(jù)進行重點監(jiān)測，此時將加速度信息作為一種隨機噪聲進行處理。在損失一定監(jiān)測精度的情況下降低監(jiān)測成本。

4.3 數(shù)據(jù)采集

4.3.1 RTK定位數(shù)據(jù)

RTK定位數(shù)據(jù)包含RTK接收機的經(jīng)度、緯度、天線高度、接收機的工作模式和觀測衛(wèi)星數(shù)量等信息。通過RTK定位得到監(jiān)測點的經(jīng)度、緯度和天線高度，則監(jiān)測點位移量dLRTK與豎直位移量dL1和水平位移量dL2滿足如下關(guān)系

(26)

dL1=H1-H2

(27)

(28)

其中θ1、φ1分別為在起始時刻監(jiān)測點的經(jīng)緯度；θ2、φ2為在觀測時刻監(jiān)測點的經(jīng)緯度；R為地球半徑，取平均值6371km。

已知目前RTK接收機動態(tài)定位精度可達(dá)cm級，平面精度最高可達(dá)1cm，高程精度可達(dá)2cm。觀測時間足夠長的情況下，RTK接收機的靜態(tài)定位精度可達(dá)mm級，平均誤差在1～2mm。考慮到山體滑坡在平穩(wěn)形變時期形變速度極為緩慢，因此對于RTK接收機而言可以視為靜態(tài)定位，形變監(jiān)測精度可達(dá)mm級，使mm級的數(shù)據(jù)融合具有可行性。

4.3.2 土工帶傳感器數(shù)據(jù)

基于土工帶電阻值計算其形變量dL的計算公式如下

(29)

(30)

式中，R(mΩ)是土工帶在觀測時刻的電阻測量值；R0(mΩ)是土工帶在布置完成時的初始電阻值；L0(mm)是土工帶在布置完成時的初始長度；k是比例系數(shù)；T(℃)為土壤溫度。

因為土工帶排布網(wǎng)絡(luò)的相對尺寸較小，所以忽略山體凹凸性，在對一套監(jiān)測設(shè)備進行分析時，可以看作將土工帶排布網(wǎng)絡(luò)安置在監(jiān)測點處的切平面上。土工帶傳感器排布方式如圖5所示。10條土工帶分成縱橫兩組，各有5條土工帶。每一條土工帶上有10個土工帶傳感器。10條土工帶按照網(wǎng)格狀排列，空心點是土工帶傳感器，實心點是2個土工帶傳感器重合點，RTK接收機安置在網(wǎng)絡(luò)的中心位置。

圖5 土工帶傳感器與RTK接收機的排布網(wǎng)絡(luò)

在山體上布設(shè)土工帶排布網(wǎng)絡(luò)時，使方形網(wǎng)絡(luò)的一邊沿山體等高線，該方向定義為橫向；2個臨邊垂直于山體等高線，該方向定義為縱向。對土工帶數(shù)據(jù)進行中值濾波和平滑處理后得到監(jiān)測點的橫向位移x和縱向位移y，則由土工帶傳感器測得的監(jiān)測點的形變量dLtu為

(31)

土工帶傳感器測試設(shè)備如圖6所示，測試結(jié)果如表3所示。測試結(jié)果顯示，土工帶傳感器的形變監(jiān)測數(shù)據(jù)可以達(dá)到mm級的監(jiān)測精度，為高精度地質(zhì)形變監(jiān)測提供了可能，使mm級的數(shù)據(jù)融合具有可行性。

圖6 土工帶傳感器地質(zhì)采集設(shè)備

表3 土工帶傳感器設(shè)備測試

4.3.3 無人機攝影測量數(shù)據(jù)

通過無人機近景攝影測量技術(shù)[24-25]獲得監(jiān)測區(qū)域的三維點云數(shù)據(jù)和傾斜攝影處理結(jié)果。使用ArcMap對多期測量數(shù)據(jù)中的影像控制點進行匹配后，在某一期的三維模型中加入另一期的相對三維模型，計算出前后兩期攝影測量中滑坡監(jiān)測點處的形變量dLUAV，并獲得形變區(qū)域的位置和大小。目標(biāo)區(qū)域的數(shù)字高程地圖(Digital Elevation Map，DEM)和數(shù)字水平地圖分別如圖7和圖8所示。無人機攝影測量技術(shù)的引入，在提供形變監(jiān)測信息的同時，還提供了地形地貌、建筑分布、交通情況等諸多信息要素，直觀地提供滑坡監(jiān)測點的位置信息，有利于災(zāi)害防治工作。

圖7 數(shù)字高程地圖

圖8 數(shù)字水平地圖

以大疆精靈Phantom無人機攝影測量為例，像元大小為2.41μm，相機焦距為8.8mm，飛行高度約為36.5m時，可以滿足客戶1cm的監(jiān)測精度要求，如圖9所示。考慮到在平穩(wěn)形變階段TUAV≈30d≈15T，所以在進行數(shù)據(jù)融合計算時,由無人機攝影測量得到的形變信息dLUAV的精度在1mm內(nèi)，在尺度與分辨率方面與土工帶傳感器形變量dLtu和RTK測量形變量dLRTKmm級精度保持一致。

圖9 無人機攝影測量精度分析

5 結(jié)果分析

本文按照已知的數(shù)據(jù)內(nèi)容和格式，綜合考慮典型案例中山體滑坡發(fā)生前后的形變規(guī)律，對數(shù)據(jù)進行模擬設(shè)置。基于模擬設(shè)置的警報級別滑坡監(jiān)測數(shù)據(jù)對本文的數(shù)據(jù)處理算法進行仿真分析。最后通過對比模擬設(shè)置數(shù)據(jù)與算法處理結(jié)果，檢驗山體滑坡形變監(jiān)測數(shù)據(jù)處理算法的功能。詳細(xì)分析內(nèi)容如下，主要包括警報級別數(shù)據(jù)模擬設(shè)置、抗差自適應(yīng)Kalman濾波融合結(jié)果分析、基于灰色預(yù)測理論和蠕變切線角判據(jù)的預(yù)測結(jié)果分析、抗差自適應(yīng)Kalman濾波算法對預(yù)測結(jié)果的影響分析4個方面的內(nèi)容。

5.1 實驗數(shù)據(jù)模擬設(shè)置

5.1.1 Kalman濾波效果展示數(shù)據(jù)模擬設(shè)置

為直觀清晰地展示通過抗差自適應(yīng)Kalman濾波算法進行數(shù)據(jù)融合的效果，本文模擬設(shè)置某期采集到的500條土工帶傳感器形變數(shù)據(jù)和500條RTK接收機形變數(shù)據(jù)，模擬設(shè)置無人機攝影測量到的滑坡監(jiān)測點的加速度信息，將模擬數(shù)據(jù)代入抗差自適應(yīng)Kalman濾波算法中進行數(shù)據(jù)融合。數(shù)據(jù)模擬設(shè)置如表4所示。

表4 抗差自適應(yīng)Kalman濾波數(shù)據(jù)模擬設(shè)置說明

5.1.2 警報級別數(shù)據(jù)模擬設(shè)置

為了檢驗滑坡監(jiān)測算法在預(yù)警等級為警報等級時的監(jiān)測效果，按照處理算法的結(jié)構(gòu)并綜合考慮典型案例中山體滑坡發(fā)生前后的形變規(guī)律模擬出12期的數(shù)據(jù)進行仿真分析，前8期的數(shù)據(jù)為測量模擬數(shù)據(jù)，實現(xiàn)對第9、10、11、12期數(shù)據(jù)的預(yù)測和危險等級預(yù)警功能。警報數(shù)據(jù)模擬設(shè)置說明如表5所示，數(shù)據(jù)模擬結(jié)果如表6所示。

表6 警報級別數(shù)據(jù)模擬設(shè)置結(jié)果

5.2 抗差自適應(yīng)Kalman濾波融合結(jié)果

上述內(nèi)容已經(jīng)說明每隔2d采集500條土工帶傳感器和RTK接收機的監(jiān)測數(shù)據(jù)。針對每一期的大量監(jiān)測數(shù)據(jù)，本文結(jié)合無人機近景攝影測量計算得到的監(jiān)測點加速度信息，使用抗差自適應(yīng)Kalman濾波進行數(shù)據(jù)融合，最終得到監(jiān)測點處精確的形變量，為實現(xiàn)后續(xù)的山體滑坡預(yù)測預(yù)警功能做數(shù)據(jù)準(zhǔn)備。抗差自適應(yīng)Kalman濾波算法的數(shù)據(jù)融合效果主要通過某一期的數(shù)據(jù)融合結(jié)果進行展示。

抗差自適應(yīng) Kalman濾波前后的數(shù)據(jù)對比如圖10所示。圖10中,藍(lán)色實線為經(jīng)抗差自適應(yīng)Kalman濾波之后的滑坡監(jiān)測點的形變量，紅色點代表RTK定位的形變監(jiān)測數(shù)據(jù)，綠色點代表土工帶傳感器的形變監(jiān)測數(shù)據(jù)，黑色實線表示模擬設(shè)置的形變真值。調(diào)整無人機攝影測量的形變加速度信息，使之與土工帶數(shù)據(jù)和RTK數(shù)據(jù)相匹配，可以看出濾波融合效果良好，在真值附近收斂，本次數(shù)據(jù)融合效果展示中無人機近景攝影測量的形變加速度a=0.0015(mm·d-2)。

圖10 抗差自適應(yīng)Kalman濾波前后數(shù)據(jù)對比

從圖10可以發(fā)現(xiàn)，濾波收斂并穩(wěn)定之后數(shù)據(jù)的最大值為5.984mm，最小值為5.570mm，與最初設(shè)置的真值5.8mm之間的誤差不超過0.23mm，滿足mm級監(jiān)測精度要求。需要說明的是，由于所選取的用于展示抗差自適應(yīng)Kalman濾波效果的這一期數(shù)據(jù)具有普適性，其監(jiān)測精度即可反映普遍數(shù)據(jù)經(jīng)抗差自適應(yīng)Kalman濾波之后的監(jiān)測精度，滿足本文數(shù)據(jù)處理精度要求，為之后的山體滑坡預(yù)測預(yù)警算法提供了精確的數(shù)據(jù)支撐。

5.3 基于灰色預(yù)測理論和蠕變切線角判據(jù)的預(yù)測結(jié)果

灰色預(yù)測理論形變預(yù)測結(jié)果與原模擬數(shù)據(jù)對比如圖11所示。對比結(jié)果顯示警戒級別數(shù)據(jù)經(jīng)過Kalman濾波融合并經(jīng)過灰色預(yù)測理論得到的警戒數(shù)據(jù)預(yù)測形變量曲線與警戒級別數(shù)據(jù)的形變量曲線基本重合，說明灰色預(yù)測理論預(yù)測的形變量與模擬設(shè)置數(shù)據(jù)的形變量接近，預(yù)測效果良好。

圖11 灰色預(yù)測理論的形變預(yù)測效果圖

灰色預(yù)測理論預(yù)測等級與原數(shù)據(jù)等級對比如表7所示。對比結(jié)果顯示由灰色預(yù)測理論形變預(yù)測結(jié)果得到的預(yù)警等級與模擬設(shè)置數(shù)據(jù)的預(yù)警等級基本符合，僅在第4期數(shù)據(jù)上存在些許偏差，滑坡監(jiān)測算法具備有效的滑坡預(yù)測預(yù)警功能。

表7 灰色預(yù)測理論預(yù)測等級與原數(shù)據(jù)等級對比

警報級別數(shù)據(jù)實驗結(jié)果如表8所示，均值方差比C1為0.0296，明顯小于0.35，即均值方差比精度等級為優(yōu)；小概率誤差P1為1，明顯大于0.95，即小概率誤差精度等級為優(yōu)，所以該灰色預(yù)測模型的精度等級為優(yōu)，形變預(yù)測結(jié)果精確。經(jīng)過抗差自適應(yīng)Kalman濾波的監(jiān)測點的精確形變數(shù)據(jù)通過灰色預(yù)測理論預(yù)測之后得到的多期數(shù)據(jù)dL_Kalman_Grey_predict呈現(xiàn)上升趨勢，即每隔2d滑坡監(jiān)測點的形變量越來越大，預(yù)計第12期(8d后)采集到的形變量能達(dá)到0.43714m；形變速度呈現(xiàn)明顯上升趨勢，即滑坡監(jiān)測點做加速形變，預(yù)計第12期(8天后)的形變速率能達(dá)到0.22388m/d；蠕變切線角angle呈現(xiàn)上升趨勢，預(yù)計第12期(8天后)能夠達(dá)到88.53°；預(yù)警等級warning_results也呈現(xiàn)上升趨勢，第8期(最近數(shù)據(jù)采集時刻)的預(yù)警等級為3級，預(yù)計第9期的預(yù)警等級為3級，第10、11、12期(未來8天)的預(yù)警等級為4級，即為警報等級，未來一周內(nèi)有非常大的概率出現(xiàn)山體滑坡災(zāi)害。

表8 警報級別數(shù)據(jù)實驗結(jié)果

5.4 抗差自適應(yīng)Kalman濾波算法對預(yù)測預(yù)警的影響

為進一步說明抗差自適應(yīng)Kalman濾波算法對預(yù)測結(jié)果的影響，本文對未經(jīng)過抗差自適應(yīng)Kalman濾波融合的RTK形變監(jiān)測數(shù)據(jù)單獨通過灰色預(yù)測理論模型進行形變預(yù)測和預(yù)警等級分析，將預(yù)測預(yù)警結(jié)果與模擬設(shè)置數(shù)據(jù)的形變量和預(yù)警結(jié)果進行對比，同時與前文所述的經(jīng)過抗差自適應(yīng)Kalman濾波融合之后的滑坡預(yù)測預(yù)警效果進行對比分析。如圖12所示，未經(jīng)過抗差自適應(yīng)Kalman濾波融合的RTK形變監(jiān)測數(shù)據(jù)單獨通過灰色預(yù)測理論進行形變預(yù)測得到的結(jié)果與模擬設(shè)置數(shù)據(jù)的偏差較大，形變預(yù)測效果變差。如表9所示，未經(jīng)過抗差自適應(yīng)Kalman濾波融合的RTK形變監(jiān)測數(shù)據(jù)單獨通過灰色預(yù)測理論和蠕變切線角判據(jù)得到的預(yù)警等級與模擬設(shè)置數(shù)據(jù)的預(yù)警等級存在的偏差明顯增大，預(yù)警錯誤率增加，監(jiān)測算法的預(yù)警效果變差。由表7和表9可知，使用抗差自適應(yīng)Kalman濾波融合之后的數(shù)據(jù),進行形變預(yù)測和預(yù)警分析的效果明顯優(yōu)于使用單一形變監(jiān)測數(shù)據(jù)進行預(yù)測預(yù)警分析。可見抗差自適應(yīng)Kalman濾波融合得到精確的形變監(jiān)測數(shù)據(jù),對于形變預(yù)測和預(yù)警分析是非常重要的。

圖12 灰色預(yù)測理論的形變預(yù)測效果圖(未經(jīng)Kalman濾波融合)

表9 灰色預(yù)測理論預(yù)測等級與原數(shù)據(jù)等級對比(未經(jīng)Kalman濾波融合)

6 總結(jié)

由上述分析可知，本文設(shè)計了一種基于抗差自適應(yīng)Kalman濾波和灰色預(yù)測理論的山體滑坡監(jiān)測技術(shù)。該技術(shù)通過多種監(jiān)測方式對易發(fā)生滑坡的監(jiān)測區(qū)域進行形變監(jiān)測，通過抗差自適應(yīng)Kalman濾波技術(shù)對包括RTK定位數(shù)據(jù)、無人機攝影測量數(shù)據(jù)、土工帶傳感器數(shù)據(jù)在內(nèi)的多源數(shù)據(jù)進行數(shù)據(jù)融合，最后通過灰色預(yù)測理論實現(xiàn)對監(jiān)測區(qū)域滑坡形變程度和危險等級的預(yù)測。

與以往研究工作相比，本文所設(shè)計的滑坡監(jiān)測技術(shù)的價值主要體現(xiàn)在以下幾個方面：

1)本文所設(shè)計的滑坡監(jiān)測技術(shù)具有很強的實時性。RTK定位技術(shù)、土工帶傳感器的使用，使得山體滑坡監(jiān)測功能受天氣狀況和植被覆蓋等因素的影響大幅減弱，能夠及時發(fā)現(xiàn)和準(zhǔn)確預(yù)測到可能發(fā)生的滑坡災(zāi)害，有效預(yù)防財產(chǎn)損失和人員傷亡。

2)本文所設(shè)計的滑坡監(jiān)測技術(shù)自動化程度高。場外的設(shè)備一經(jīng)布設(shè)，不需人工定期到現(xiàn)場進行勘探。場外的設(shè)備故障能夠及時反饋到云端監(jiān)測平臺，只需有針對性地進行修復(fù)即可。與此同時，少量故障設(shè)備回傳的錯誤信息也能通過算法濾波自行消化，對預(yù)測結(jié)果造成的影響很小。

3)本文所設(shè)計的滑坡監(jiān)測技術(shù)的監(jiān)測精度高。由上述仿真分析可知，由于土工帶傳感器形變監(jiān)測、RTK靜態(tài)形變監(jiān)測和無人機攝影測量形變監(jiān)測的精度均能達(dá)到mm級，外加抗差自適應(yīng)Kalman濾波融合處理，該技術(shù)對于監(jiān)測點處的形變監(jiān)測精度可達(dá)mm級。

4)本文所設(shè)計的滑坡監(jiān)測技術(shù)的設(shè)備成本低。土工帶傳感器輔助RTK形變監(jiān)測和無人機攝影測量形變監(jiān)測，大大降低監(jiān)測成本，為大范圍的山體滑坡監(jiān)測提供了可行性。

本文目前尚存在一些不足之處，主要體現(xiàn)在受疫情影響和時間限制等多方面因素制約，缺乏實測數(shù)據(jù)進行實驗驗證，算法的分析仍停留在仿真層面，該問題待采得實測數(shù)據(jù)后即可解決。此外，本文所設(shè)計的滑坡監(jiān)測技術(shù)尚未考慮地理水文、人類活動等相關(guān)因素，后續(xù)工作計劃引入隨機森林模型對影響滑坡的諸多因素進行分析。