張慶斌,趙之星,劉志奇
(中晉環(huán)境科技有限公司,山西 太原 030002)
滑坡監(jiān)測歷來是地質(zhì)災害研究領域中的重點,大量的開采、施工或長期惡劣氣候的影響,使得原始地貌形成了多種不穩(wěn)定邊坡,這些邊坡由于受力失衡會產(chǎn)生定向移動[1]。大面積滑坡會對礦區(qū)工作人員、附近居民以及其他相關作業(yè)人員的生命和財產(chǎn)造成嚴重危害。其中突發(fā)性的滑坡現(xiàn)象危險性最大,其隱蔽不易辨別的征兆和災害發(fā)生時無法應對的緊迫性,是迄今為止災害監(jiān)測預警的難題[2]。為此大量學者和機構開發(fā)了各種監(jiān)測手段和預警系統(tǒng)。上海華測導航科技有限公司推出了基于綜合技術的自動實時監(jiān)測系統(tǒng),該系統(tǒng)應用了GNSS接收機、固定式測斜儀、裂縫計、雨量計等傳感器,以及無線傳輸裝置和太陽能電池等輔助設備[3-6]。該系統(tǒng)利用多種手段獲取信息,綜合評估目標變化,得出結果,但是其設備價格昂貴,一旦發(fā)生滑坡,系統(tǒng)設備全部損壞,監(jiān)測成本較高。成都理工大學團隊[7-9]利用機載雷達和傾斜相機對災害體進行采集,通過內(nèi)業(yè)處理生成三維激光點云實景三維模型,并對多期數(shù)據(jù)進行比較,該方法能直觀顯示出測區(qū)隱患點,為決策提供準確的依據(jù),且對植被茂密災害點隱蔽的測區(qū)尤為勝任。但該方法處理過程較長,特別是三維模型的生產(chǎn)占用大量的時間,因此監(jiān)測實時性略顯滯后。文獻[10-12]利用近景攝影測量、計算機視覺等手段,通過固定雙相機,利用影像匹配技術,結合內(nèi)外方位元素計算目標的三維坐標,并分析多期的坐標變化以達到目標位移的計算。該方法不需要在災害體上安裝設備,作業(yè)成本低[13-15]。但為了得出空間坐標,需要進行較為復雜的計算,且受制于攝影測量本身精度的限制,對突發(fā)性應急存在一定的不適應[16-18]。針對滑坡的突發(fā)性和應急性,通過固定相機對滑坡體拍照,對相片的像素進行計算,設定多層篩選條件和閾值,令計算機準確分析出實際滑坡和誤差影響或分散式局部移動的情況,做出相應決策。該方法延續(xù)了近景攝影非接觸式的優(yōu)勢,同時舍棄求取三維坐標的繁雜計算過程,僅通過影像進行判別,響應速度快,自動化程度高,對滑坡災害的應急較為適用[19-20]。
固定相機按照一定時間間隔對滑坡體進行連續(xù)拍照。相機固定即像片不動,像片不動即像平面坐標系不變,當滑坡目標不發(fā)生移動時,理論上兩幀影像提取的同名點像素坐標應相等,即
受到噪聲影響,兩幀影像同名點坐標會產(chǎn)生一定的偏差,由于噪聲服從各向同性高斯分布,使得總體上點坐標偏差方向各異,偏差均值為0,即
(1)
事實上,除了隨機噪聲影響外,特征提取中存在粗差以及樣本數(shù)量有限等影響因素,導致上述坐標偏差的均值與理想值0存在偏離,因此需要定義置信閾值,當該偏離在閾值范圍內(nèi),則同名點坐標的變化僅受到噪聲的干擾,沒有發(fā)生實際的移動。結合坐標誤差的理論分布,令
式中:Δui、Δvi分別為第k張影像與第1張影像的同名點橫、縱坐標差。
則樣本方差的估值為
(2)
式中:S為標準差。
構建樣本樞軸量
(3)
圖1 誤差服從t分布的α分位Fig.1 Alpha quantiles with errors following t distribution
反之,該情況下特征點坐標偏差的均值和理論期望之間的差異是顯著的,因此認為像點受到了噪聲以外的干擾,需進行進一步分析。
若像點的偏差均值過大,則需深入分析,以判斷變化是由粗差影響還是本身發(fā)生移動造成的,具體方法如下。
計算特征點的梯度方向:
(4)
式中:θi為第1張影像相較于第k張影像同名點發(fā)生位移的梯度方向。
繪制特征點梯度方向的直方圖,橫軸為梯度方向的角度區(qū)間,在45°~135°內(nèi),每間隔5°為1個區(qū)間,共18個,縱軸為滿足每個區(qū)間點的數(shù)量,即N1、N2、N3,…,N18,將點數(shù)量按大小排列,并提取出最多的2組進行分析。假設點數(shù)量最大的一組為Ni,將Ni個像點按照橫坐標從小到大依次計算相鄰點的距離
式中:di為同一影像上相鄰像點間的距離。
以第1個點與第2個點的相鄰距離d1為基準,定義相鄰間距的權重pi為
pi=d1/di
(5)
式中:mi為第1張影像相較于第k張影像第i個同名點發(fā)生的位移。
(6)
同時根據(jù)相似關系得出一個像素對應物方的實際距離為
GSD=f/(Lδ)
(7)
式中:GSD為1個像素對應的實際距離;f為焦距;L為拍攝距離;δ為單個像素大小。
因此發(fā)生移動的距離為D=GSD(u2+v2)1/2,其中,u、v為目標點像素坐標。
根據(jù)以上分析,總結本次影像應急監(jiān)測預警的條件,如圖2所示。
圖2 滑坡應急預警決策流程
由于數(shù)量稀少,危險性等原因,對理想突發(fā)性滑坡測區(qū)的選定較為困難,考慮時間成本等因素,選擇進行模擬試驗。試驗區(qū)域選擇太原市中國冶金地質(zhì)總局新辦公樓施工區(qū)西北角的人工邊坡,邊坡坡度正切值約0.7(坡面角35°),土層15 cm以下布置1.2 cm×1.2 cm的襯墊,以便模擬滑坡移動過程。為盡可能接近滑坡體結構,在邊坡上摻雜了碎石組成的松散堆積體。采用SONY_A7Ⅱ相機,將其固定在改裝后的三腳架上,在距目標3 m處粗略整平,并保持不動。首先在坡面未動的情況下拍攝首張照片(圖3a),然后分3種情況拍攝后續(xù)照片:保持坡面不動(圖3b),人為使坡面上某點或較為分散的零星點位發(fā)生移動(圖3c),人工拉拽襯墊,坡面滑移(圖3d)。分別將圖3b、圖3c、圖3d與圖3a進行對比分析,驗證方法的正確性。
將圖3b中影像首先通過優(yōu)化后的SIFT算法進行同名點提取與匹配,利用基于本方法開發(fā)的程序進行處理分析,計算每組像片間同名點的像素坐標變化情況,驗證本方法對3種位移情況的判斷及響應效果。
2.2.1 未移動影像組處理分析
第1對影像選擇初始影像和滑坡未動影像進行軟件分析,通過同名點匹配、統(tǒng)計匹配點樣本數(shù)量,計算同名點像素坐標差值、平均值及標準差估計值,依據(jù)假設檢驗求取檢驗統(tǒng)計量的絕對值,并與其服從的分布拒絕閾的臨界點進行比較。處理結果如圖4所示。共提取特征點1 858個,匹配成功的同名點1 553個,紅色箭頭為像點移動的方向,經(jīng)比較計算得出檢驗統(tǒng)計量為1.78,α取0.05時,t分布的拒絕域臨界點為1.96,顯然統(tǒng)計量小于臨界點。
2.2.2 局部分散式移動影像組處理分析
第2對影像為初始影像與分散局部移動個別點位的影像之間的對比,如圖5所示。經(jīng)計算得出像點坐標差值組成的檢驗統(tǒng)計量已超過臨界值位于拒絕域內(nèi),需要進一步分析,在分析階段,統(tǒng)計梯度方向角度45°~135°內(nèi)每5°一個區(qū)間的方向數(shù)量,如圖6所示。
圖3 試驗序列影像Fig.3 Experimental sequence images
圖4 未移動試驗序列影像Fig.4 Unmoved experimental sequence images
圖5 分散局部移動試驗序列影像Fig.5 Scattered local moving experimental sequence images
圖6 同名點梯度方向數(shù)目Fig.6 Number of gradient directions with same name
其中選定數(shù)量最多的前2組,進行第2步分析,依據(jù)程序計算2組中特征點的平均相鄰距離和影像間同名點的移動量,計算得出各自數(shù)值見表1。
表1 梯度方向數(shù)量前2組的距離Table 1 Top two distance of gradient directions
2組方向中發(fā)生滑移的特征點相鄰距離加權平均值遠大于所設定的3 pix,而影像間發(fā)生的位移也都遠大于5 pix。
2.2.3 移動影像組處理分析
第3組影像為初始影像和人工模擬移動影像對比,匹配點如圖7所示。
經(jīng)程序計算,像點坐標差值組成的檢驗統(tǒng)計量已遠大于t分布的拒絕閾臨界點,需要進一步分析,在分析階段,統(tǒng)計得出45°~135°內(nèi)按照每5°一個區(qū)間的方向數(shù)量(圖8)。
圖8 同名點方向梯度方向數(shù)目Fig.8 Number of gradient directions with same name
選定數(shù)量最多的前2組,進行第2步分析,依據(jù)程序計算2組中特征點的平均相鄰距離和影像間同名點的移動量,計算得出各自數(shù)值見表2。
表2 梯度方向數(shù)量前2組的距離Table 2 Top two distance of gradient directions
2組方向中發(fā)生滑移的特征點相鄰距離加權平均值小于方法中規(guī)定的3 pix,而影像間發(fā)生的滑移距離也都遠大于5 pix,同時根據(jù)相似幾何關系,按照拍攝距離3 m計算得出該滑移量大約為5.1 m。
1)第1組試驗計算結果表明:影像間同名點像素坐標差平均值,其絕對值小于臨界值,說明不存在顯著性差異,造成的位移很大概率上是由于受到高斯噪聲的影響,因此可以判定未發(fā)生移動,并終止分析。
2)第2組試驗中,僅在較為分散的局部區(qū)域內(nèi)發(fā)生變動,其像素坐標差絕對值超出臨界值,經(jīng)計算得出兩影像間像點位移加權平均值超過設定閾值5 pix,而發(fā)生變化的點位相鄰間距加權平均值卻大于設定的3 pix,表明盡管點位發(fā)生了較為明顯變化,但點位之間的距離較遠,位置分散,區(qū)別于真實滑坡情況下的連片移動,故判定疑似人為所致,或滑坡前兆,發(fā)送提醒信號、不報警。
3)第3組試驗為已發(fā)生移動情況下的分析,該像素坐標差平均值遠大于臨界值,在分析階段,篩選出梯度方向數(shù)量最多的2組,計算得出變化點相鄰間距加權平均值小于設定值3 pix,同時影像間像點滑移量的加權平均值遠大于設定閾值5 pix,證明發(fā)生變化的點分布緊密且連接成片,目標位移也遠超過高斯噪聲和粗差的影響,故可判定目標發(fā)生滑坡移動,并報警。
1)提出了依靠單目序列影像進行滑坡位移的判斷算法,利用像素變化的統(tǒng)計規(guī)律,加入梯度方向、目標點相鄰距離和影像間目標點位移等參數(shù),并設定初始篩選、數(shù)據(jù)分析、相應決策3個環(huán)節(jié),提高了該算法的魯棒性。
2)該算法可對高斯噪聲影響、分散式局部移動、和大面積滑坡3種現(xiàn)象進行有效甄別,省略了計算三維空間坐標的繁瑣過程,提高了運行速度,節(jié)省生產(chǎn)成本和保證人身安全的前提下,對突發(fā)性滑坡事故應急工作能起到一定的幫助。
3)由于受試驗條件及時間限制,目前僅利用模擬試驗進行測試,且對于閾值設定的合理性和普適性還需要更多的實地場區(qū)進行驗證,這也是后期的優(yōu)化工作。