黃河內(nèi)蒙古河段封凍期冰晶結(jié)構(gòu)特征分析

王 娟，黃 樾，李志軍，鄧 宇，張邀丹

（1.鄭州大學(xué) 水利科學(xué)與工程學(xué)院，河南 鄭州450001；2.大連理工大學(xué) 海岸和近海工程國(guó)家重點(diǎn)實(shí)驗(yàn)室，遼寧 大連116042；3.黃河水利科學(xué)研究院，河南 鄭州450003）

每年冬季，我國(guó)北方的河流、湖泊、水庫(kù)都會(huì)發(fā)生冰凌現(xiàn)象。冰凌會(huì)阻礙正常的航運(yùn)、發(fā)電、供水，破壞大壩、護(hù)坡、橋梁等結(jié)構(gòu)，還可能產(chǎn)生冰塞、冰壩等凌汛災(zāi)害，對(duì)人民的生命財(cái)產(chǎn)造成威脅。其中，針對(duì)黃河冰凌理論的研究，以河冰的生消演變和冰塞冰壩的形成機(jī)理為主［1-3］。河冰作為一種天然復(fù)合材料，不同地區(qū)、不同時(shí)間形成的河冰性質(zhì)存在差異，這主要是由時(shí)空的差異造成了河冰細(xì)觀結(jié)構(gòu)的差別，而冰晶體結(jié)構(gòu)很大程度上決定了河冰的物理力學(xué)性質(zhì)［4-5］。研究表明，粒狀冰的脆性強(qiáng)度與其晶體直徑的平方根成反比［6］，多晶冰的極限抗壓強(qiáng)度與粒徑平方根的倒數(shù)成正比例關(guān)系［7］。因此，對(duì)晶體結(jié)構(gòu)進(jìn)行研究，有利于了解冰的物理性質(zhì)，揭示冰的內(nèi)在破壞機(jī)理［8］。

相關(guān)學(xué)者嘗試過(guò)許多方法對(duì)冰晶體結(jié)構(gòu)進(jìn)行觀測(cè)。19世紀(jì)中期，Seligman使用了在紙上制作鉛筆拓片的方法［9］，使風(fēng)化的晶界清晰地顯示出來(lái)。近年來(lái)，基于偏振光的原理觀測(cè)晶體結(jié)構(gòu)的相關(guān)設(shè)備不斷發(fā)展，如自動(dòng)冰組構(gòu)分析儀［10］、費(fèi)氏臺(tái)［11］等。2006年，Kipfstuhl等［12］使用電子顯微鏡和CCD攝像機(jī)對(duì)晶體結(jié)構(gòu)進(jìn)行了觀察，通過(guò)顯微結(jié)構(gòu)映射法得到了晶體的高分辨率圖片。其中，費(fèi)氏臺(tái)因操作方便、便于攜帶、試驗(yàn)條件要求低等優(yōu)點(diǎn)，成為觀測(cè)冰晶結(jié)構(gòu)常采用的儀器。

冰晶觀測(cè)技術(shù)一直在提高，但對(duì)于晶體圖片的后續(xù)處理技術(shù)仍停留在人工操作的階段。目前常用的方法是用Adobe Photo Shop CS軟件配合手工勾選出晶體圖片中晶體的邊界，然后使用圖像處理軟件計(jì)算每個(gè)冰晶粒的尺寸［13］。該方法需要花費(fèi)大量的時(shí)間，同時(shí)手工勾選邊界存在誤差?；谝陨媳尘?，本文基于MATLAB的數(shù)字圖像分析能力，提出使用邊界提取功能分割冰晶圖像并計(jì)算晶粒尺寸的方法，分析黃河冰晶的結(jié)構(gòu)分布特征。

1 冰晶圖像的獲取

2019年2月，在黃河內(nèi)蒙古典型河段什四份彎道處（見圖1（a））采集了4組冰胚。什四份彎道位于黃河內(nèi)蒙古段緯度較高的位置，海拔990 m，每年冰期持續(xù)100 d左右。冰胚的取樣點(diǎn)選取在黃河北岸距岸邊60～150 m處。該處為黃河冬季主河道，冰層表面平整無(wú)積雪，冰下水深6～7 m，冰蓋泥沙含量較低，排除冰凌堆積成冰。采集冰坯的具體位置見表1。采冰時(shí)首先利用GPS確定正北方向，按照約30 cm寬、30 cm長(zhǎng)在冰面上切割冰塊，使用電鋸和人工板鋸將冰胚切斷取出。

圖1 黃河冰樣提取及冰晶圖像獲取

取樣完成后，將冰樣沿垂直方向分層，用鋸骨機(jī)切成厚度相同的小塊，并使用刨刀將冰塊的一面修平，貼在溫?zé)岬牟Ａ?，使冰塊與玻璃片充分黏結(jié)。隨后，再使用刨刀將冰塊削成厚度約1 mm的薄冰片，在費(fèi)氏臺(tái)偏光鏡（見圖1（b））下觀測(cè)冰薄片，拍攝冰晶圖片，并記錄每張圖片對(duì)應(yīng)的實(shí)際尺寸。

表1 取樣點(diǎn)坐標(biāo)

2 冰晶圖像處理

根據(jù)黃河冰晶結(jié)構(gòu)觀測(cè)結(jié)果，截取黃河冰晶圖片的有效部分并轉(zhuǎn)化為灰度圖像，對(duì)圖像進(jìn)行濾波除噪和圖像分割?；谶吘墮z測(cè)的方法提取冰晶體的邊界，對(duì)部分邊界特征不明顯的區(qū)域進(jìn)行校正補(bǔ)充。提取出完整邊界之后，通過(guò)圖像分析獲取黃河冰晶的結(jié)構(gòu)特征。

2.1 圖像前處理

黃河冰的晶體結(jié)構(gòu)原始圖見圖2（a），在原圖中截取出最大內(nèi)接矩形（見圖2（b）），同時(shí)對(duì)照刻度尺記錄矩形圖片的實(shí)際尺寸，用于后期的單位換算。將圖片轉(zhuǎn)化為灰度圖，提高圖片處理的效率。使用偏光鏡拍攝的冰晶圖片具有較大的噪聲，常見的空間域去噪方法有鄰域平均法、選擇平均法、中值濾波法、空間低通濾波法［14］。通過(guò)效果對(duì)比，選用中值濾波法對(duì)圖像進(jìn)行降噪處理。中值濾波的原理是在每一個(gè)像素點(diǎn)的周圍取一個(gè)鄰域，將該點(diǎn)的像素值取為鄰域的中值，從而消除與周圍像素相差較大的噪聲?；谥兄禐V波的原理，分別使用邊長(zhǎng)為4像素、8像素、16像素、24像素、30像素的鄰域?qū)ΡЩ叶葓D進(jìn)行降噪處理。通過(guò)比較，最終選用邊長(zhǎng)為24像素的鄰域進(jìn)行濾波，見圖2（c）。

圖2 冰晶圖像預(yù)處理

傳統(tǒng)的圖像分割方法有邊緣檢測(cè)法［15］、閾值分割法［16］、區(qū)域法［17］等。傳統(tǒng)的閾值法適用于將圖像分割為有限的兩個(gè)或多個(gè)區(qū)域，不適用于晶粒較多的冰晶圖像的分割。同時(shí)，冰晶圖像噪聲大，部分冰晶體的特征不明顯，存在像素值接近、邊界不明顯的相鄰晶粒。這些現(xiàn)象給使用區(qū)域法和聚類算法帶來(lái)較大難度?；谏鲜銮闆r，最終選擇邊緣檢測(cè)方法對(duì)冰晶圖像進(jìn)行分割。常見的邊緣檢測(cè)算子有Robert算子、Sobel算子、Prewitt算子、Laplacian算子、Canny算子［18］等。前4種算子都屬于局域窗口梯度算子，對(duì)噪聲的抵抗性差，而Canny算子對(duì)噪聲的抵抗性較強(qiáng)，故采用Canny算子對(duì)冰晶圖像進(jìn)行邊界提取。從圖3（a）可看出，大部分冰晶邊界已經(jīng)被提取出來(lái)，但仍存在少量邊界的缺失和少量噪音。進(jìn)一步處理后，邊界清晰連續(xù)而且沒(méi)有明顯錯(cuò)位的現(xiàn)象（見圖3（b））。

圖3 冰晶圖像分割

2.2 冰晶圖像分析

采用MATLAB，通過(guò)編程將圖3（b）中的連通區(qū)域進(jìn)行編號(hào)并賦予不同的灰度值，每一個(gè)連通區(qū)域代表一個(gè)冰晶粒。通過(guò)灰度值的不同，計(jì)算每個(gè)連通區(qū)域像素點(diǎn)的數(shù)量以及其邊界的像素點(diǎn)數(shù)量，根據(jù)記錄的圖像寬度，將像素面積和像素周長(zhǎng)換算為實(shí)際面積和周長(zhǎng)。由于冰晶粒形狀不規(guī)則，因此用與晶粒面積相等的等效圓直徑來(lái)描述冰晶粒的直徑，面積和直徑的計(jì)算公式如下:

式中:S為實(shí)際面積；D為等效直徑；S n為區(qū)域像素點(diǎn)數(shù)量；k為每個(gè)像素點(diǎn)對(duì)應(yīng)的實(shí)際長(zhǎng)度。

由于直線邊界和斜線邊界每個(gè)像素點(diǎn)對(duì)應(yīng)的實(shí)際長(zhǎng)度不同，因此在計(jì)算實(shí)際周長(zhǎng)L時(shí)需要將斜邊界上的像素點(diǎn)額外乘以系數(shù)見式（3），其中L直和L斜分別為直線邊界和斜線邊界的像素點(diǎn)數(shù)量。

3 結(jié)果分析

3.1 冰晶?？傮w分布規(guī)律

對(duì)共計(jì)50張黃河冰晶圖片進(jìn)行處理，獲取7 000多個(gè)冰晶粒尺寸數(shù)據(jù)。等效直徑小于3 mm的冰晶粒占12.99%，等效直徑為3～6 mm的冰晶粒最多（見圖4），占計(jì)算總量的31.94%。隨著等效直徑的增大，冰晶粒數(shù)量逐漸減少，等效直徑大于21 mm的冰晶粒僅占計(jì)算總量的5.67%。

圖4 冰晶粒尺寸總體分布

冰晶體類型可分為粒狀冰和柱狀冰兩種，對(duì)兩種不同的冰晶結(jié)構(gòu)進(jìn)行分析，見圖5。粒狀冰的晶粒分布與圖4中總體分布的規(guī)律相似，其中:等效直徑0～3 mm的晶粒數(shù)量小于3～6 mm的晶粒數(shù)量，等效直徑大于3 mm的晶粒占比隨著晶粒尺寸的增大逐漸減小。粒狀冰的大尺寸晶粒較少，等效直徑大于12 mm的晶粒僅占計(jì)算總量的6.57%。柱狀冰尺寸偏大，小于3 mm的冰晶數(shù)量較少（僅占計(jì)算總量的6.57%），等效直徑3～12 mm晶粒數(shù)量隨晶粒等效直徑的增大而減少，但大于12 mm的晶粒占比比粒狀冰的多，晶粒數(shù)量下降的幅度相對(duì)平緩。這說(shuō)明柱狀冰的大尺寸冰晶粒數(shù)量更多，且尺寸分布較為均勻。

圖5 粒狀冰與柱狀冰晶粒尺寸分布

3.2 冰晶粒水平分布規(guī)律

圖6 展示了兩組代表性冰胚的冰晶水平圖像處理結(jié)果，左邊冰樣從上到下分別是1號(hào)冰胚8、16、24、32、40、48、56、63 cm深度的水平切片，右邊冰樣從上到下分別是2號(hào)冰胚5、15、25、35、45、55、68 cm深度的水平切片。從圖6可看出，晶粒分割圖與原圖基本一致，邊界連續(xù)完整，誤差很小。圖7（a）為計(jì)算得到的1號(hào)冰胚每層水平切片的等效直徑分布，可以看出，8 cm和24 cm深度處的晶粒等效直徑分布集中在0～3 mm和3～6 mm區(qū)間，大于6 mm的晶粒占比不到10%，原因是1號(hào)冰胚0～24 cm深度是粒狀冰，晶粒較小。40～63 cm之間的冰晶體是柱狀冰，隨著深度的增加，大尺寸冰晶粒呈現(xiàn)逐漸增加的趨勢(shì)。63 cm深度處的晶粒異常的原因是晶粒等效直徑大，圖片拍攝的晶粒數(shù)量少，統(tǒng)計(jì)規(guī)律不明顯。圖7（b）為計(jì)算得到的2號(hào)冰胚每層水平切片的等效直徑分布，2號(hào)冰胚0～3 cm和20～35 cm深度是粒狀冰，粒狀冰層等效直徑小于6 mm的晶粒占比達(dá)到了70%。隨著深度的增加，柱狀冰層晶粒占比增大。

圖6 水平方向冰晶體圖像識(shí)別

通過(guò)對(duì)兩組冰晶體水平方向晶粒等效直徑分布的對(duì)比可知，所有粒狀冰層的晶粒占比均在3～6 mm區(qū)間達(dá)到最大，隨后大幅度減??；所有柱狀冰層的晶粒占比分布較為均勻，差異并不顯著。

3.3 冰晶粒豎直分布規(guī)律

圖7 不同深度水平方向冰晶粒尺寸分布

兩組代表性的冰晶體垂直分布結(jié)果見圖8。圖8（a）為1號(hào)冰胚垂直切片的等效直徑分布。深度0～24 cm是粒狀冰層，晶粒等效直徑維持在3～4 mm之間；深度24～63 cm是柱狀冰層，晶粒等效直徑隨著深度的增加而增大，并在56 cm深度達(dá)到最大值11.5 mm；在深度63 cm處晶粒等效直徑出現(xiàn)了小幅的減小，推測(cè)是冰花凝結(jié)的作用導(dǎo)致。圖8（b）為2號(hào)冰胚垂直切片的等效直徑分布，晶粒等效直徑在5 cm深度達(dá)到極大值6.2 mm，隨后又在25 cm深度減小至3.38 mm，再往下晶粒等效直徑逐漸增大，并在68 cm深度達(dá)到最大值12.2 mm。冰晶粒等效直徑的變化規(guī)律整體反映了冰晶結(jié)構(gòu)粒狀—柱狀—粒狀—柱狀的交替變化趨勢(shì)。

圖8 垂直方向不同位置冰晶粒尺寸分布

4 結(jié) 論

通過(guò)濾波除噪、邊界提取等方法對(duì)黃河內(nèi)蒙古典型河段什四份彎道處所取的冰胚進(jìn)行了處理，得到了黃河冰晶粒的完整邊界，并通過(guò)統(tǒng)計(jì)像素點(diǎn)的方法計(jì)算了黃河冰晶粒的尺寸、面積和周長(zhǎng)，分析了黃河內(nèi)蒙古河段凌汛期冰晶粒的分布特征。

（1）等效直徑3～6 mm的冰晶粒占比最高，對(duì)于等效直徑大于3 mm的冰晶粒，其占比與等效直徑成反比。粒狀冰的平均粒徑小于柱狀冰，大尺寸冰晶粒相對(duì)較少，而柱狀冰的大尺寸晶粒相對(duì)較多，且分布較為均勻。

（2）隨著深度的增加，冰層大尺寸冰晶粒占比逐漸增大。隨著冰晶粒等效直徑的增大，粒狀冰占比呈現(xiàn)先增大后減小的趨勢(shì)，其中冰晶粒等效直徑小于6 mm的粒狀冰占比在70%以上，柱狀冰的冰晶粒分布則較為均勻。

（3）冰晶粒的垂直分布與冰胚的晶體結(jié)構(gòu)相關(guān)。隨著深度的增加，粒狀冰層和柱狀冰層交替出現(xiàn)，晶粒尺寸也隨之出現(xiàn)不同的極值點(diǎn)。