張 思 煒
(南京信息工程大學江蘇省氣象探測與信息處理重點實驗室 江蘇 南京 210044) (南京信息工程大學濱江學院 江蘇 無錫 214105) (南京信息工程大學江蘇省大氣環境與設備技術協同創新中心 江蘇 南京 210044)
冰雹是一種不規則的固態降水物,它是一種嚴重的自然災害。冰雹災害具有局地性顯著、突發性強和持續時間短等特征。冰雹災害的發生不僅給人民生命財產造成了巨大的損失,也給該地區的農業、建筑財產、通信設備、電力設施、交通設施等帶來極大的危害[1]。目前冰雹災害等級可以通過人工測量冰雹顆粒的特征參數來評估,但是在實際測量當中,過多的冰雹勢必會相互接觸粘連,為了完成測量,手動分散費時費力而且冰雹多為類圓形或不規則圖形,在測量時勢必會造成精度誤差。為了提高測量的效率,實現冰雹災害等級評估的自動化,通過計算機智能來獲取自然環境中冰雹顆粒的特征參數成為一個值得研究的課題。
現有冰雹災害的研究大多通過對冰雹云的多普勒雷達回波特征分析,來實現冰雹災害的預警與識別[2-4],但是在自然環境下對冰雹特征參數測量的研究尚屬空白。文獻[5-9]為測量物體特征參數的一些方法。董堯培等[5]結合HSV顏色空間和Hough變換來測量鋁錠的厚度。李秋潔等[6]通過分析MLS測量系統分辨率建立變尺度格網來測量靶標葉面積,效果好、精準度高。馮青春等[7]采用光度體視覺技術對蔬菜秧苗葉片曲面形態進行測量,測量結果與人工測量結果相比,差異較小。孫玉婷等[8]將葉子的幾何特征參數的線性回歸值作為輸入,構建了SVM回歸模型,對葉面幾何特征參數進行了估計,結果均方根誤差以及平均相對誤差都比較低。李文勇等[9]運用凸包理論結合分水嶺算法對自然場景下未成熟蘋果進行直徑的測量,測量誤差較小。
為了能夠從自然環境下獲得冰雹顆粒的特征參數,分析了冰雹和背景的顏色特征,本文研究并且提出HSI模型分割方法將冰雹與背景完全分離開來,針對粘連情況運用基于分水嶺算法將其分割,之后提取冰雹特征參數。該方法可以客觀地反映災害的發生程度并且能有效地提高冰雹災害等級劃分的自動化水平。
HSI顏色空間來源于人類的視覺系統,分別通過色調(H)、飽和度(S)和亮度(I)對色彩進行描述。相比于RGB顏色空間,HSI各個分量相互獨立能將色彩變化描述更為清楚,因此HSI顏色空間是處理彩色圖像理想的顏色空間。RGB顏色空間轉換到HSI顏色空間可由如下轉換公式實現[10]:
(1)
(2)

(3)
(4)
式中:R、G、B、H、S、I表示冰雹圖片中的像素值。通過以上的轉換公式,可以將RGB顏色空間轉換成HSI顏色空間。
為了精準地獲得冰雹顆粒的顏色范圍,這里準備了50幅冰雹圖片,這50幅冰雹的圖片僅僅包含了冰雹和單一顏色的背景。本實驗考慮為測量結果提供較為準確的測量標準,通過圓形的模具來模擬冰雹,通過圓形模具可獲得測量參數的真實值。冰雹的顏色通過控制白色顏料濃度的大小來實現,統計50幅二值化冰雹圖片的像素值,計算得到單位像素面積比例尺0.006 59 cm2、單位像素周長比例尺為0.007 73 cm,單位像素長度比例尺為0.007 64 cm。通過分析冰雹圖片中的像素分布并且統計顏色直方圖,可以得出冰雹的顏色范圍。圖1為待提取的冰雹圖片,圖2為冰雹圖片S分量的直方圖。

圖1 待提取冰雹圖片

圖2 統計S分量的直方圖
之后再分別統計H、S、I三個分量的范圍。經過實驗發現冰雹圖片中的H分量無法對冰雹的分割做出貢獻,所以這里分別統計S分量和I分量的顏色范圍從而將冰雹圖片二值化。
分水嶺算法[11]是一種經典的區域分割方法,在分離粘連物體上有著廣泛的應用。基本思想是將整幅圖像看成一片崎嶇不平的山地模型。其中各個像素點的灰度值代表著該點的海拔高度,圖像中局部極小值區域對應地貌中的匯水盆地。在模擬注水的過程中,水面慢慢地從山底向山頂擴展,隨著水位的上升,匯水盆地將會被水淹沒合并成同一個區域。為了阻止區域的合并,在各個匯水盆地之間修建“堤壩”形成分水嶺,從而將整個拓撲地貌分割成若干區域。分水嶺算法分割的具體步驟:令M1,M2,…,Mi為梯度圖像h(x,y)的局部極小值,C1,C2,…,Cj為相對應匯水盆地中所有坐標集合。Cj[N]為n階段Cj的集合,C[N]表示所有Cj[N]的集合。hmin和hmax代表著圖像h(x,y)中灰度值的最小值與最大值。則有:
(5)
令Cj[N]為所有匯水盆地的集合,則有:
(6)
隨著水位以整數從n=hmin+1到n=hmax+1遞增,若設g[n]為位于平面h(x,y)=n以下點的集合,根據算法g[n]中的坐標點將會被標記成黑色,其他的坐標點將標記成白色形成一幅二值圖片。算法設定選取C[hmax+1]=g[hmin+1],然后進入遞歸調用的階段,根據n階段的C[n-1]來求解C[n]。設M為g[n]連通分量,則對于連通分量M,有三種情況:
(1)M∩C[n-1]=?。
(2)M∩C[n-1]包含C[n-1]中的一個連通分量。
(3)M∩C[n-1]包含C[n-1]中多個連通分量。
當新的極小值產生并滿足(1)時將C[n]與C[n-1]合并;若滿足(2)則將C[n]與C[n-1]相應的部分合并;若滿足(3),則在C[n]中構建分水嶺,阻止匯水盆地的結合。
傳統的分水嶺算法對噪聲相當敏感,若直接將其應用到冰雹圖像上面則會產生嚴重的過分割現象,分割結果如圖3(a)所示。為了能夠分離粘連的冰雹,本文提出一種改進的分水嶺算法。傳統的算法是對形態學重建后的灰度圖直接進行二值化,分割結果往往會保留大量不規則區域同時形成虛假邊緣,為了保證分割效果,本文直接對冰雹原圖運用HSI模型分割來濾除復雜的背景并運用數字形態學除去分割后多余的噪聲,保證了冰雹圖像邊緣的完整性,分割結果如圖3(b)所示。

圖3 改進分水嶺分割過程
接著再對二值化圖像進行歐氏距離變換[12],距離變換是對圖像中每個像素點到最近非零像素點之間距離的運算。歐氏距離計算公式為:

(7)
式中:p(x1,y1)和q(x1,y2)分別為前后景像素點的坐標。將圖像中每個連通域的內部像素點到非內部像素點最短距離構成的集合定義為S(x,y),Min、Max為集合S(x,y)中的最大值與最小值,G(x,y)為連通域中每個內部像素點經過距離變換所對應的灰度值,則歐氏距離變換的公式為:
(8)
將二值圖像中每個相應的像素點按照式(7)式(8)計算,便可以把二值圖轉換成距離變換圖,變換結果如圖3(c)所示。在距離圖中連通域內的像素點距離連通域邊界越遠,像素點所在位置的灰度值就越高,亮度也就更強。
對于距離變換圖而言依然可能會存在由暗噪聲導致的偽局部極小值點。若此時進行分水嶺變換仍會導致過分割現象。為了得到滿意的分割效果,本文采用形態學擴展極小值變換來抑制多余的偽極小值點。擴展極小值變換的技術原理是通過人為設定一個閾值將匯水盆地深度低于該閾值的局部極小值除去。匯水盆地深度h定義為圖像中某一連通域邊界處的梯度值與連通域內部最小梯度值之差,抑制極小值深度h小于所選閾值的過程如下:
(9)

(10)

經過上述算法處理,可以發現存在粘連的冰雹被分成一個個連通的區域。通過對這些連通區域的測量即可得到冰雹的特征參數。
(1) 冰雹的面積測量。采用統計像素的方法來計算冰雹的面積。通過計算每個連通域冰雹的像素個數就可以較為精準地計算每個連通域的面積。
(11)
式中:S0為實際冰雹的面積;S1為參照物的實際面積;P0為冰雹圖像所包含的像素點;P1為參照物圖像所包含的像素點。

(3) 冰雹的直徑測量。最小外接矩形(MER)[14]可以用來對一個不規則物體進行矩形計算。一個冰雹顆粒可以有許多個不同的外接矩形,其中外接矩形面積最小的被稱作是最小外接矩形。考慮到冰雹顆粒的形狀多為類圓形或者不規則多邊形,不能直接將最小外接矩形的長和寬作為冰雹的直徑。在實際測量中需要對冰雹從多個方向上進行直徑的測量從而獲得較為精準的等效直徑。為了符合現實冰雹直徑的測量要求,這里提出一種改進的最小外接矩形測量方法。測量原理如圖4所示。

圖4 冰雹直徑測量原理
(12)
(13)
(14)
式中:Nn、Em、Sg、Wh分別為最小外接矩形與冰雹上下左右四個邊界交點的集合;(xupi,yupi)和(xdownj,ydownj)為上下邊界交點的像素坐標;(xlefti,ylefti)和(xrightj,yrightj)為左右邊界交點的像素坐標;n、m分別代表上下邊界頂點的個數;g、h分別代表左右邊界頂點的個數。圖4中每個類型的線型表示邊界上不同的點,由點引出的線段代表著兩點之間的距離。按照式(12)從N1到Nn依次計算到S1,S2,…,Sm的距離取平均可得上下邊界冰雹等效直徑d1,左右邊界按照式(13)從W1到Wn依次計算到E1,E2,…,Em的距離取平均即可得左右邊界冰雹等效直徑d2,最后再按照式(14)計算即為最終的冰雹直徑d。冰雹特征參數提取整體算法的具體流程如圖5所示。
實驗運行環境:Windows 10操作系統,處理器為Intel?CoreTMi5 CPU 8400@2.8 GHz,內存為8 GB。圖像處理軟件MATLAB 2018b。
本實驗一共準備了100組冰雹圖片,因為篇幅有限選擇三幅不同背景的冰雹圖片作為實驗對象展現實驗過程。冰雹圖片如圖6(a)所示。從左到右的背景分別為沙地、草地和操場。首先將三幅冰雹圖片從RGB顏色空間轉換到HSI顏色空間,然后通過先驗的顏色范圍將冰雹顆粒從背景中分割出來,其效果如圖6(b)所示。從三幅圖片的分割結果可以看出大部分冰雹顆粒已經被分割出來,但還是存在少量噪聲,需要進行進一步的處理。第一幅圖片整體效果較好,僅存在少數多余的連通域以及噪聲;第二幅圖雖然冰雹顆粒整體已經被分割出來但是一些青草也被分割出來且存在青草的遮擋問題,使右下方的冰雹無法完整地分割出來;第三幅圖的背景存在一些細小的噪聲并且有一部分的冰雹顆粒因為過于透明而被濾除,針對以上出現的問題首先濾除一定閾值的連通域,之后構建圓形結構元素運用形態學開閉的操作將這些噪聲濾除。去除噪聲后的效果如圖6(c)所示。通過濾除噪聲,雖然有極個別的冰雹顆粒的形狀受到影響,但是大部分冰雹的輪廓清晰,外形結構無缺失。從濾除噪音的冰雹二值圖像中可以看出部分從背景中分割出來的冰雹之間相互接觸,相互接觸的冰雹會形成一大片連通區域,若不進行處理測量時冰雹特征參數的測量值會偏大,從而影響測量的精準性。因此在進行冰雹特征參數測量時需要將相接觸顆粒分離。下面使用1.4節提到的分水嶺算法對濾波后的冰雹圖片進行分割。圖6(d)為分割結果,相互接觸的冰雹顆粒已經被分水嶺線分隔開,每一個冰雹對應一個連通區域。

圖6 冰雹分離過程
選取其中一幅冰雹二值圖像,首先對每個連通域進行標記并且做出最小外接矩形,結果如圖7(a)和圖7(b)所示。可以看出共有26個冰雹而且每一個分離的冰雹都具有一個最小外接矩形。通過統計像素法計算冰雹顆粒的面積,為了實現計算單個冰雹的面積首先將每個標記的冰雹分離,部分冰雹如圖7(c)和圖7(d)所示。可以看出分離的冰雹形狀完整、邊緣清晰。對分離后的冰雹顆粒使用統計像素法來計算其面積,并將其測量結果與真實值相比較,結果如圖8(a)所示。圖8(c)中的橫坐標為每個冰雹的編號,縱坐標為本文方法的計算值,中間的直線為人工測量值,當測量值的散點離直線越遠,偏差越大;當測量值的散點離直線越近,偏差越小。可以看出一部分測量值位于直線下方,說明測量結果略微偏小。


圖7 實驗結果

圖8 冰雹特征參數的測量值與計算值

對分離的冰雹運用最小外接矩形法同時提取邊緣像素坐標,結果如圖7(e)和圖7(f)所示,可以看出單個冰雹的邊界被完整地提取出來,找出冰雹在上、下、左、右四個方向上與最小外接矩形交點的像素坐標來計算冰雹的直徑,但在實際測量中最小外接矩形與冰雹圖像邊緣的交點可能存在多個,針對這種情況運用式(12)-式(14)來計算冰雹的等效直徑。計算結束后,將計算值與真實值進行比較,其結果如圖8(c)示。表示測量值的散點圍繞直線均勻上下分布,說明測量值與真實值接近。
本實驗共準備100組冰雹圖片進行測量,分別計算出冰雹樣本的平均測量值、平均誤差、均方根誤差RMSE、決定系數R2[15]來衡量測量值的精準度。同時將本文算法與相關的算法進行精準度的對比。
(1) 均方根誤差RMSE(Root Mean Square Error):測量值與真實值誤差的平方與觀測次數n比值的平方根,計算公式為:
(15)
式中:xi為本實驗算法取得的樣本i的測量值;yi為樣本i的測量值;n為樣本數。
(2) 決定系數R2:相關系數的平方,能夠直觀地反映測量值與真實值的擬合程度。計算公式為:
(16)

冰雹特征參數的測量結果如表1所示。測量結果的RMSE越小,表明測量值與真實值越接近,擬合得越好。R2的值越接近1表明數據相關性越好、誤差越小。在冰雹面積的測量上,分別采用坐標紙法和長寬系數法與本實驗所用的統計像素法進行精準度的比較。從表1可以看出三種方法的平均誤差分別為0.184 9、0.677 3和0.615 7,且從RMSE和R2的值上來看統計像素法的測量結果更加精準。坐標紙法需要人為地將冰雹的輪廓投映到坐標紙上,其過程難免會產生誤差;長寬系數法對形狀規則的物體測量效果較好但對冰雹這類形狀多為類圓形或不規則物體測量精度較差。

表1 不同冰雹特征參數測量方法的比較分析
在冰雹周長的測量過程中分別應用Freeman鏈碼和Canny算子追蹤冰雹邊緣,通過統計像素之間距離來計算周長。兩種方法的RMSE為0.208 1和0.301 6,R2為0.881 4和0.754 5。Canny算子運算時間較長同時提取圖像的邊緣為多層像素,冗余的像素會導致測量結果存在誤差;Freeman鏈碼檢測到的圖像邊緣均為單一像素同時考慮到了像素之間傾斜方向的距離問題,測量效果更好。
最后在冰雹直徑的測量上分別采用改進的最小外接矩形法與游標卡尺法來進行實驗對比,直徑平均測量值分別為2.895 5和2.951 7,平均誤差分別為0.084 3和0.140 5。游標卡尺法測量較為精準但是逐個測量費時費力而且在讀數時容易導致測量誤差,相比之下改進的最小外接矩形法測量效率更高,測量結果更為精準。本文方法的測量值誤差均在2%到6%之間,可精準地提取冰雹的特征參數。
本文以自然環境下冰雹的分割和特征參數的測量作為研究對象,選用基于HSI顏色空間的閾值法將冰雹從背景中分割,運用改進的分水嶺方法分割粘連冰雹邊界,提取的冰雹輪廓清晰、完整性好,避免了過分割的情況。冰雹面積的計算選用統計像素法,該方法計算過程簡單、可靠性好;周長算法使用Freeman鏈碼,根據鏈碼值累加得到冰雹的周長;冰雹的直徑通過改進的最小外接矩形來測量,改進后的算法適用性廣,可用于類圓形或不規則圖形的測量。從冰雹特征參數的測量結果來看本文算法精準度較高,RMSE和R2的值分別趨近于0與1,平均誤差在2%到6%之間。下一步將研究在冰雹大量堆疊的情況下如何測量冰雹的特征參數,進一步提高算法的適用范圍。