汪恩良，胡勝博，韓紅衛(wèi)，劉承前

（1.東北農(nóng)業(yè)大學 水利與土木工程學院，黑龍江 哈爾濱 150030；2.黑龍江省寒區(qū)水資源與水利工程重點實驗室，黑龍江 哈爾濱 150030）

1 研究背景

冬季封河和春季開河階段，我國北緯30°以北的河流普遍存在冰情危害［1］，隨著河流斷面冰凌分布密度的增大，冰凌撞擊水工建筑物的概率增大，且大量流凌在特定河段，如河道束窄處、淺灘、急彎、連續(xù)彎道、尚未解凍冰蓋的前緣等處受阻，形成冰凌阻水體，往往嚴重阻塞過水斷面，使上游水位顯著壅高，造成凌汛災(zāi)害［2-3］。在冰情嚴重的年份，春季流凌的撞擊會對橋墩等水工建筑物造成不同程度的破壞，每年國家和沿海各省區(qū)都要投入巨大的人力物力和財力專門用于凌汛災(zāi)害的預(yù)防和處理，故監(jiān)測河流斷面的冰凌變化信息，掌握冰凌分布密度數(shù)據(jù)是冰凌災(zāi)害預(yù)防的關(guān)鍵［4-5］。

河流斷面冰凌分布密度是在指定斷面內(nèi)，冰凌分布面積占河流單位面積的百分比，可以反映出河流斷面上不同時間冰凌流量的大小，是預(yù)測冰凌災(zāi)害的重要指標［6］。目前獲得冰凌分布密度常用的方法有人工望遠鏡目估法和衛(wèi)星遙感監(jiān)測法：（1）人工望遠鏡目估法［7］，通過人為觀測，估算河流斷面冰凌分布密度。這種方式操作簡單，但對觀測人員要求較高，需要有比較豐富的經(jīng)驗，既要保持眼睛的健康，也要保持清醒的頭腦，但無法對冰凌分布密度進行準確地量化，缺乏客觀性；（2）衛(wèi)星遙感監(jiān)測技術(shù)［8］，使用衛(wèi)星拍攝的冰凌圖像提取冰凌分布密度，由于衛(wèi)星遙感方式的時間間隔大，空間分辨率較低等因素的影響，導致冰凌分布密度監(jiān)測存在較大誤差。國內(nèi)外學者對海冰分布密度展開了大量的研究，提出了使用船載攝像機作為觀測工具，對航行途徑海域的海冰狀況進行觀測，Mu?ramoto 等［9］在1993年將攝像技術(shù)首次應(yīng)用到海冰觀測的實踐中；Hall 等［10］使用船載攝像設(shè)備對海冰進行觀測，基于閾值技術(shù)提取海冰分布密度，并使用SSM/I 的數(shù)據(jù)進行驗證；Worby 等［11］使用多種衛(wèi)星數(shù)據(jù)結(jié)合船上觀測研究南極海冰變化等；盧鵬［12］利用傳統(tǒng)圖像形態(tài)學算法，構(gòu)造出一種新的海冰邊緣檢測算法，提高拍攝圖像中單個海冰識別效率；鄧霄等［13］基于攝像機和透視變換圖像校正算法設(shè)計了冰凌密度監(jiān)測系統(tǒng)，改善了傳回圖像存在“近大遠小”的問題；攝像機冰凌監(jiān)測技術(shù)能夠在一定程度上滿足冰凌分布密度監(jiān)測的要求，且圖像校正算法能夠改善圖像存在的畸變，但由于圖像校正算法自身的限制，對攝像機布置有一定的要求，且無法徹底消除圖像畸變［14］。

隨著無人機（unmanned aerial vehicles，UAV）和傳感器技術(shù)的不斷提高，小型無人機在低空遙感測量領(lǐng)域得到了廣泛的應(yīng)用［15］。對比人工望遠鏡目估法和衛(wèi)星遙感監(jiān)測法，可搭載多傳感器的無人機平臺具有效率高、成本低、操作靈活、更適合復雜河道環(huán)境等特點，不僅克服了人工獲取影像數(shù)據(jù)的困難，還可以從根源上避免了岸邊側(cè)視拍攝導致的圖像畸變，提高原始數(shù)據(jù)精度［16-18］，為監(jiān)測大范圍冰凌變化信息提供了新的技術(shù)手段。

本文通過無人機低空遙感技術(shù)對目標河段進行俯視拍攝，將視頻逐幀提取圖片，利用matlab 仿真軟件編程，通過頂帽變換算法對圖片進行亮度均衡化處理，基于OTSU算法（最大類間方差法）提取最佳閾值，對圖像進行閾值分割，最終通過面積濾波降噪，計算河冰分布密度，繪制河冰分布密度圖，實現(xiàn)對河冰分布密度連續(xù)監(jiān)測。

2 研究區(qū)及數(shù)據(jù)預(yù)處理

2.1 研究區(qū)域概況黑龍江位于我國最北端，是中俄界河，流域處于寒溫帶地區(qū)，東亞大陸季風氣候特點突出［19］，且黑龍江上游流經(jīng)山區(qū)，河道狹窄曲折，雨量徑流比較充沛，每年流凌期，冰凌洪水以疊加組合的形式促成倒開江向下游發(fā)展，導致漠河段易發(fā)生凌汛災(zāi)害（圖1）［20］，黑龍江漠河段封凍期從11月份持續(xù)到次年4月，是典型的季節(jié)性封凍河流。

圖1 冰凌觀測地理位置示意

2.2 遙感影像獲取及預(yù)處理本文試驗數(shù)據(jù)于2021年4月25日通過無人機平臺獲取，采用大疆DJI四旋翼小型無人機御mavic 2 pro 進行定點拍攝，無人機體積為214 mm×91 mm×84 mm，可折疊機身，方便隨身攜帶，搭載哈蘇L1D-20c（CMOS），有效像素為2000 萬，可錄制4k HDR 視頻，每次航拍的天氣條件均為晴朗，風速低于5 級，飛行高度設(shè)置為400 m，起飛點與水面高程差為5 m，鏡頭垂直向下，懸停拍攝。

現(xiàn)場通過錄制冰凌運動視頻，內(nèi)業(yè)處理時每隔29.97 幀截取圖片影像，并通過圖像處理工具對影像數(shù)據(jù)進行裁剪，去除邊緣異常值。

為實現(xiàn)對河道冰凌分布圖像數(shù)據(jù)的量化分析，首先，采用無人機拍攝不同高度下A1 紙，得到無人機遙感影像拍攝高度與單位像素點所代表的實際面積的擬合公式，如圖2所示。

圖2 無人機高度與單位像素點代表實際面積擬合曲線

3 研究方法

利用冰凌錄像數(shù)據(jù)可以掌握被監(jiān)測河道的實時冰凌分布密度與流速情況。在此基礎(chǔ)上通過matlab編程對冰凌分布圖像進行量化分析，提供新的流凌分布密度監(jiān)測方式，為相關(guān)部門制定防凌方案提供理論依據(jù)和決策支持。

3.1 灰度圖像轉(zhuǎn)化通過攝像機拍攝視頻能夠?qū)崿F(xiàn)河冰分布密度數(shù)據(jù)提取，主要是因為在太陽光照射條件下，河冰與河水對可見光反射的亮度和波長不同，使得它們在攝像機拍攝出來的影像中所顯示的亮度有明顯區(qū)別，即河冰要比河水更亮［21］。若將河冰分布圖像轉(zhuǎn)化為灰度圖，這種不同也體現(xiàn)在灰度值上的差別。采用加權(quán)平均法將拍攝的彩色圖像轉(zhuǎn)化為灰度圖像，計算公式如下所示：

式中：I（x，y）為灰度圖像；R（x，y）、G（x，y）和B（x，y）分別為彩色圖像的RGB 分量值。三個加權(quán)系數(shù)是根據(jù)人的亮度感知系統(tǒng)調(diào)節(jié)得到的標準化參數(shù)。

3.2 圖像亮度均衡化處理采用灰度形態(tài)學中的組合算法，即頂帽變換算法對河冰圖像進行亮度均衡化處理［22］?；叶刃螒B(tài)學包括四種基本算法：灰度膨脹算法、灰度腐蝕算法、開運算和閉運算?；叶扰蛎涍\算［23］是對原灰度圖像進行“加長”或“變粗”；灰度腐蝕算法［24］是在灰度圖像上“收縮”或“細化”；開運算［25］將前兩種算法進行結(jié)合，對灰度圖像進行先腐蝕運算，再膨脹運算，可以除去圖像背景中比結(jié)構(gòu)元素尺寸更小的亮度明亮細節(jié)；頂帽變換［26］則是在開運算對原始圖像進行處理的基礎(chǔ)上，用原始灰度圖像再減去開運算結(jié)果的運算。這種方法能夠進一步增強圖像的對比度，適用于處理具有暗背景（河水）、亮物體（河冰）特征的圖像，頂帽變換算法推導公式為：

腐蝕運算：

膨脹運算：

開運算：

頂帽變換計算：

式中：f（x，y）為圖像M（x，y）為結(jié)構(gòu)元素，AM是M（x，y）的定義域，結(jié)構(gòu)元素M（x，y）對圖像f（x，y）的灰度膨脹記為f+M，灰度腐蝕記為f-M，開運算記為f·M，頂帽變換記為fM。

3.3 OTSU算法閾值選取是閾值分割算法的核心［27］。通過確定合理的灰度閾值，可以實現(xiàn)圖像分割，將河冰與河水區(qū)分開，計算河冰所占比例，從而得到河冰分布密度。

OTSU算法是借助最小二乘法原理在直方圖技術(shù)上推導出來，具有簡單、速度快等特點，是一種常用的閾值選取方法，適合于物體目標與背景灰度差明顯的情況［28］。該算法以灰度分布均勻性的度量單位為方差，方差值越大，說明構(gòu)成圖像的兩部分差別越大［29］。當部分目標（冰凌）錯分為背景（水）或部分背景（水）錯分為目標（冰凌）時，就會導致類間方差變小，因此使類間方差最大的分割閾值意味著計算精度變高［30］。

將原始冰凌圖像f（x，y）進行灰度化處理，得到圖像灰度范圍為{0， 1，…，l-1}級。灰度值為i的像素數(shù)設(shè)為ni，推導出灰度圖像總像素值為從而灰度值為i像素出現(xiàn)的概率為：pi=ni/N，其中選擇閾值t，將灰度圖像劃分為兩類，即D1：{0，1， 2，…，t}，D2：{t+1，t+2，t+3，…，l-1}，D1和D2出現(xiàn)的概率分別為：

D1和D2的灰度均值分別為：

圖像的總體灰度均值為：

可求出D1和D2的類間方差為：

類間方差越大，說明D1和D2之間差距越大，即冰凌區(qū)域與水分割效果越好，所以確定最佳分割閾值，就是在求兩類區(qū)域的類間方差最大值，即：

求出最佳閾值，即可實現(xiàn)圖像分割，生成二值圖像。

3.4 面積去噪在野外河道環(huán)境中，河流表面波浪反光往往是影響計算結(jié)果的重要因素。由于陽光反射產(chǎn)生的光點灰度值遠高于圖像分割閾值，導致反光點被錯分為冰凌，冰凌分布密度增大，但反光點面積小且較為密集，只存在于圖像的局部區(qū)域，不具備如高斯噪聲［31］、椒鹽噪聲［32］等常見圖像噪聲的特征，無法通過均值濾波法［33］、自適應(yīng)維納濾波法［34］、中值濾波法［35］等去除。同時，無人機能夠搭載的熱成像傳感器分辨率較低，難以滿足白天流凌監(jiān)測需求，因此通過對局部采用面積去噪法剔除反光噪聲［36］，即：目標面積S小于X個像素點時，將會被剔除，實現(xiàn)圖像分割。通過反復代入X值，比較圖像的去噪效果，確定X取值區(qū)間。

4 結(jié)果與分析

4.1 頂帽算法對比分析在野外環(huán)境中，由于存在陽光照射角度，光照不均勻等因素會導致原始冰凌圖像背景亮度不均勻，直接采用OTSU算法進行圖像分割，會導致閾值分割失敗，導致計算結(jié)果存在巨大誤差，如圖3所示。將采用頂帽變換算法前后得到的冰凌分布密度對比可以發(fā)現(xiàn)，兩者之間相差最大為19.00%，最小為0.71%。采用頂帽變換算法后，計算獲得的流凌分布密度均變大，但由于冰凌表面存在雜質(zhì)，部分碎冰過薄顏色與江面接近，導致計算結(jié)果與人工選擇冰凌得到的冰凌分布密度結(jié)果相比仍偏小。標準的光源和監(jiān)測環(huán)境在野外環(huán)境下不可能達成，難以與高效數(shù)字化監(jiān)測流程相結(jié)合，而采用頂帽變換算法可以補償不均勻的背景亮度，使圖像背景亮度均衡化，分割結(jié)果與原始冰凌圖像（圖4（a））情況一致。

圖3 三種方法冰凌分布密度對比

從原始灰度圖像三維圖（圖4（b））中可以得到圖中邊界區(qū)域參差不齊，并且從經(jīng)過開運算得到的背景灰度三維圖（圖4（c））中可以得到河道中心區(qū)域灰度值高于河道兩側(cè)，河道左上角區(qū)域相對于周邊區(qū)域灰度值略大，這是由于光照不均勻?qū)е卤尘盎叶炔痪鶆?。由于原始冰凌圖像背景中河道中間區(qū)域和河道左上角區(qū)域灰度值較高，在通過OTSU算法選定閾值后，由于閾值選取過高，導致河道兩側(cè)冰凌和灰度值較小的薄冰、碎冰被錯誤的分割為河水，即河水面積擴大，冰凌分布密度偏小，與實際情況不符。對原始冰凌圖像采用頂帽算法后，補償灰度值不均勻的背景，從經(jīng)過頂帽變換處理后的灰度三維圖（圖4（d））中可以看出，圖像亮度均衡化，冰凌與河水成功分割。通過人為選取確定冰凌分布密度為46.23%，對比頂帽變換處理前后，冰凌分布密度誤差由9.00%減小至0.97%，因此頂帽算法與OTSU算法相結(jié)合，能夠解決圖像亮度不均衡導致計算結(jié)果誤差過大的問題。

圖4 冰凌圖像處理過程

4.2 OTSU 和迭代算法對比分析閾值分割主要由四種方法：人工選擇法、直方圖技術(shù)法、OTSU算法和迭代法。人工選擇法［37］是通過肉眼識別判斷閾值，并反復代入閾值，比較圖像分割效果，縮小閾值選取區(qū)間，直至最終確定閾值，操作簡單但效率低，且缺乏客觀性，受到野外環(huán)境的影響，圖像總體亮度也會隨時間改變而改變，不適合批量處理所有冰凌圖像；直方圖技術(shù)法［38］適用于處理目標與背景灰度對比大的圖像，即具有典型雙峰直方圖特征，通過人為選擇兩處波峰之間波谷的灰度值為分割閾值，如圖5所示原始圖像灰度值不存在兩個明顯的波峰，所以本研究不采用直方圖技術(shù)法。迭代法［39］是通過選擇初始閾值，基于閾值迭代逼近思想，直至得到滿足給定的準則的最佳閾值；OTSU算法［40］按照圖像的灰度特性將圖像分為背景和目標兩部分，以目標和背景的類間方差最大為閾值選擇原則，自動確定最佳閾值。

圖5 原始冰凌灰度直方圖

本研究從2400 張原始冰凌圖像中選取流凌初期、流凌中期和流凌末期3 個階段的冰凌圖像，根據(jù)冰凌的顏色、形態(tài)，人為將冰凌從圖片中分割出來，轉(zhuǎn)化為二值圖像，將得到的冰凌分布密度作為該情況下的標準，再分別采用迭代法和OTSU算法進行閾值分割處理，冰凌分布密度處理效果如圖6所示。

圖6 冰凌分布密度計算效果對比

結(jié)果表明，對比OTSU算法和迭代算法的冰凌分布密度計算結(jié)果，流凌初期，OTSU算法和迭代算法與人工提取冰凌圖像計算結(jié)果相差分別為1.39%和0.03%，圖像分割效果較好；流凌中期，OTSU算法與人工提取冰凌圖像誤差最大為1019 m2；流凌末期，冰凌分布密度較低，由于碎冰、薄冰灰度值接近河水，對算法的計算精度要求更高，OTSU算法和迭代算法流凌分布密度計算結(jié)果相差可達3.6%，迭代算法相較于人工處理的冰凌分布密度，誤差高達4360 m2，而OTSU算法相較于迭代算法更適合處理冰凌分布密度低的圖像，OTSU算法與人工提取冰凌圖像計算結(jié)果相差僅為0.5%。且對比OTSU算法和迭代算法的運行時間，兩者相差約3.5 s，在大批量實時處理原始冰凌圖像情況下，OT?SU 算法能夠節(jié)約計算時間，且不同冰凌分布密度的圖像分割精度高。

4.3 面積去噪分析由于波浪與陽光等因素地存在，部分時間段的圖像會存在陽光反射的情況，在本研究中選取存在陽光反射的流凌圖像，如圖7（a）原始局部圖像所示，左下角存在高密集，小面積的反光點，影響冰凌分布密度計算結(jié)果。反光點的灰度值（188 ～ 255）遠大于圖像分割閾值（129），導致冰凌分布密度偏大0.15%，通過以10 個像素的速度增加X值，發(fā)現(xiàn)當X值小于40 時，白色像素點數(shù)量迅速減小，當X值大于40 時，白色像素點減小速度變緩，因此確定該高度下噪聲面積在40 個像素點以下，剔除反光噪聲效果如圖7（b）（c）所示，對比去噪前后的圖像，該方法對于陽光反射此類噪聲去除效果較好，且只針對圖像的局部區(qū)域，不存在影響冰凌圖像的形態(tài)特征。該方法實現(xiàn)簡單，去噪效率高，可根據(jù)實際情況結(jié)合OTSU算法，進一步減小誤差。

圖7 采用面積去噪處理反光效果

4.4 冰凌分布密度計算結(jié)果無人機在一定高度懸停，定點拍攝視頻，并從視頻中每隔3 s 截取圖片，同時在河岸邊設(shè)置攝像機，每隔30 s 自動拍攝。首先將河道以外部分裁剪舍去，采用matlab 仿真軟件編寫程序，對2400 張無人機拍攝圖像和2235 張攝像機拍攝圖像進行閾值分割，得到閾值分割圖像，計算結(jié)果作為每張圖像對應(yīng)的冰凌分布密度。冰凌分布密度隨時間變化情況如圖8所示，本次冰凌觀測發(fā)生在2021年4月24日至26日之間，如圖8所示，4月25日8 時前由于陽光照射角度的影響，且清晨可見光線不足，通過攝像機得到的冰凌圖像分割效果差，導致冰凌分布密度誤差大，計算結(jié)果偏小，與實際情況不符。4月25日8 時后，光線充足，均能實現(xiàn)冰凌圖像分割，但攝像機拍攝圖像存在“近大遠小”的問題，且冰凌浮出水面一定高度，遮擋了冰凌間空隙，導致計算結(jié)果相較于實際情況偏大。攝像機拍攝冰凌圖像選取的范圍尺度較小，冰凌分布密度容易出現(xiàn)急劇波動。隨著時間的推移，通過無人機采集的冰凌分布密度總體呈現(xiàn)下降趨勢，期間監(jiān)測到的冰凌分布密度最大為81.05%，且在2021年4月25日15 時45 分前，冰凌分布密度穩(wěn)定在60%之上，4月25日15 時左右，上游發(fā)生冰塞，到4月25日15 時45 分后觀測斷面冰凌分布密度迅速下降，直至4月25日16 時35 分，冰凌分布密度下降至33.2%，隨后回升至60%，符合現(xiàn)場實際情況，在4月26日8 時17 分以后，目標斷面冰凌分布密度穩(wěn)定在10%以下，無大面積冰凌存在，開江結(jié)束。

圖8 冰凌分布密度曲線

5 結(jié)論

為了實時監(jiān)測河流斷面的冰凌變化信息，預(yù)防冰凌災(zāi)害，本文探討了使用無人機低空遙感結(jié)合圖像處理技術(shù)提取冰凌分布密度，利用matlab 仿真軟件編程，最終獲得了河冰分布密度隨時間變化的動態(tài)曲線圖，得到了以下結(jié)論：

（1）閾值分割法計算量小，性能穩(wěn)定，在圖像分割算法中應(yīng)用最為廣泛。閾值選取是閾值算法的核心，通過合理確定灰度閾值，從而實現(xiàn)圖像分割，而最大類間方差法（OTSU）因其直觀性和實現(xiàn)的簡單性，類間方差越大，說明背景和目標兩部分差異越大，適合大批量提取冰凌密度圖像數(shù)據(jù)。

（2）冰凌分布密度計算結(jié)果受到光照不均勻的影響。頂帽變換算法通過對流凌觀測圖像進行亮度均衡化處理，與OTSU算法相結(jié)合，最終將與標準的誤差縮小至1.5%以內(nèi)，從而解決背景亮度不均勻條件下閾值的選擇。

（3）面積圖像分割法能夠解決由于冰凌監(jiān)測影像中波浪反光的影響，提高圖像識別的準確性，而且該方法實現(xiàn)簡單，去噪效率高。

（4）通過無人機低空遙感結(jié)合圖像處理技術(shù)，本次觀測流凌發(fā)生在2021年4月24日至26日，期間監(jiān)測到的冰凌分布密度最大為81.05%，且在4月25日15 時45 分前，冰凌分布密度穩(wěn)定在60%之上，在4月25日15 時左右，上游發(fā)生冰塞，從監(jiān)測結(jié)果可以發(fā)現(xiàn)冰凌分布密度迅速下降，符合現(xiàn)場實際情況，在4月25日16 時35 分左右，上游冰塞解除，冰凌分布密度回升至60%，在4月26日8 時17 分以后，觀測斷面流凌分布密度穩(wěn)定在10%以下，開江結(jié)束。