吳明清,弋曉康,羅華平,李傳峰,唐曉燕,陳坤杰
(1.南京農業大學工學院,南京210031;2.塔里木大學機械電氣化工程學院,阿拉爾843300;3.江蘇省農業農村廳綠色食品辦公室,南京210036)
新疆因光照豐富,降雨量少,晝夜溫差大的獨特環境而出產優質紅棗[1]。近年來,新疆紅棗產量隨種植面積逐年擴大,紅棗已經成為新疆特色農業經濟的支柱產業,紅棗產業的快速發展對紅棗采后加工的技術和設備提出了更高的要求[2]。分級是紅棗采收后貯藏、加工及銷售過程中的重要環節,同時也是提高紅棗質量和商品化價值的重要手段。目前紅棗主要的分級方法有人工、機械和機器視覺技術。人工分級的效率和精度都比較低;機械分級克服了人工分級缺點,極大的提高了分級效率和質量,但機械分級只能根據紅棗的長短徑進行分級,分級指標比較單一,不能夠真實反映紅棗的形狀[3]。形狀是影響紅棗等級和價格的重要因素。機器視覺技術利用圖像傳感器,獲取紅棗的大小、形狀、顏色、紋理和表面缺陷等信息,實現紅棗的分級[4-15]。新疆干制紅棗的品種有駿棗和灰棗,其中駿棗長徑、短徑和長度各不相同。由于二維圖像只能提供紅棗的兩維信息,而且還因視角變化,造成幾何形狀和投影面積的差異[16],會降低紅棗的分級精度。
體積和表面積都是農產品重要的物理特征,準確獲取體積和表面積對農產品的大小分級有重要的意義[17]。但傳統類球體體積的測量多采用排水法測量,表面積采用削皮或切片的方法進行測量[18]。效率低、無法實現實時測量,無法用于農產品的在線分級。近年來,已有不少文獻研究不規則農產品的體積和表面積的測量方法。李景彬等[10]將紅棗近似成類圓柱體,提取紅棗二維圖像特征參數得出類圓柱體體積,建立圖像體積與質量的線性回歸方程。Clayton等[19]分別以球體和橢球體為模型估算4種蘋果的表面積,采用削皮測量蘋果的表面積,兩種模型分別低估了實際表面積15%和18%。Uyar 等[20]通過三維激光掃描儀獲取農產品三維數字圖像的體積和表面積,分別計算雞蛋和球體的體積和表面積,其中球體體積和表面積的計算誤差小于1%,雞蛋測體積計算誤差在1%~3.15%,雞蛋測表面積計算誤差小于1%。熊妮娜等[21]應用三維激光掃描儀系統,通過掃描獲取單株樹木的三維空間的點云數據,將樹冠近似為多個圓臺體,求他們的體積之和來計算樹冠體積。龔愛平等[22]利用三維線框模型圖像處理方法測量不規則球類形狀的農產品的體積和表面積,用所提及方法、排水法和掃描圖像處理法對不規則的柑橘、蘋果和梨進行體積和表面積測量,結果顯示3 種方法的相關系數大于0.98;在測量精度方面,當圖像數量大于10幅時,測量精度大于98%,測量數量小于7幅時則精度下降明顯,測量系統能在5 s內完成被測物體的體積和表面積的測量。但是,目前國內針對紅棗體積和表面積檢測的研究幾乎處于空白,根據體積進行紅棗分級的研究還沒有相關報道。
本文搭建紅棗圖像采集的簡易裝置,利用圖像處理技術進行紅棗圖像的實時采集和特征提取,再通過特征參數構建紅棗的三維模型并計算模型的體積和表面積,為將來開發基于圖像技術和體積參數的紅棗分級設備提供理論和技術基礎。
三維輪廓模型是表示物體幾何形體的模型,是計算機中對物體的形狀可視化表現的一種方式,模型由對象輪廓的邊緣點按映射關系轉換得到,其特點是存儲量小、速度快、操作靈活和處理時間較短[22]。要建立被測物體的三維輪廓模型,需要對旋轉過程的物體規定時間內按等時間間隔測量采集N 副圖像,完成等間隔角度圖像的采集,對采集的到圖像進行處理,提取圖像邊緣輪廓點的坐標,按以下原理進行映射轉換。
在二維空間,假設向量p,q 圍繞原點旋轉,角度為θ,順時針為正,逆時針為負。圖1 展示向量p,q 繞原點旋轉,得到新向量p',q',構造矩陣如式(1)所示:

圖1二維向量原點旋轉圖Fig.1 Two-dimensional vector origin rotation diagram
圖像輪廓點的坐標點圍繞著坐標系z 軸進行旋轉,同一點旋轉前后的關系式如(2)所示:

式中(x,y)為旋轉前圖像坐標的位置;θ 為圖像的旋轉角,(°);(x',y')為旋轉后圖像坐標的位置。設坐標系(xoy)內建立空間直角坐標系xyz,空間直角坐標系xyz 轉動后形成新的直角坐標系x'y'z',如圖2 所示,轉換公式如(3)所示。

圖2 三維旋轉坐標系Fig.2 Three-dimensional rotational coordinates

式中(x,y,z)T為空間直角坐標系中任意一點;(x',y',z')T為空間直角坐標系x'y'z'中的坐標,(△x,△y,△z)T為空間直角坐標系xyz 到空間直角坐標系x'y'z'的平移參數。其中Rx(φ),Ry(α),Rz(θ)為空間直角坐標系xyz到空間直角坐標系x'y'z'的旋轉矩陣,φ,α,θ 表示紅棗輪廓將其圍繞x,y,z,軸旋轉角度(°),如(4)、(5)、(6)所示:

在公式(4)中,假設當坐標點圍繞x 軸進行旋轉,物體不發生位移,(△x,△y,△z)T=(0,0,0)T,其中(x,y,z)T表示坐標點在平面直角坐標系zoy 上的坐標,(x',y',z')T表示輪廓上的像素點在平面直角坐標系z'o'y'上的坐標。
多輪廓三維點云模型如圖3 所示,兩個相鄰輪廓的間角為△θ,輪廓個數用N,輪廓的個數N=180°/△θ,假設輪廓個數N=20,則輪廓在z軸x心上的投影如圖3a、b所示,ri和ri+1為投影面輪廓邊緣點和輪廓軸心點之間的距離,投影面輪廓邊緣相鄰兩點A,P 和軸心點Z0構成三角形A=△Z0PA,相鄰輪廓點和軸心點Z0形成的投影面積可以用式(7)表示:

式中S 為輪廓點的投影面積,cm2,Ai,Ai+1,Ai+2…,Ai+n為投影面輪廓邊緣相鄰兩點和軸心Z0點組成三角形的面積,cm2;n 為輪廓點形成三角的個數。在圖3c 中,在Z 軸方向,相鄰輪廓投影面形成一個橢圓臺,多輪廓模型的體積可以用式(8)表示。

式中V 為多輪廓模型的近似體積,cm3;ab 為多輪廓模型縱軸的長度,cm;Si,Si+1為第i,i+1 層投影面的面積,cm2;△h為相鄰投影面的高度,cm。
在圖3b中,輪廓投影面的邊緣點的周長可以用式(9)表示。

式中C為投影面積的周長;Li,Li+1,Li+2…,Li+n為輪廓邊緣相鄰兩點的距離,cm,n 為邊緣輪廓點的個數。在圖3c 中,在z 軸方向,相鄰輪廓投影面形成一個橢圓臺,多輪廓模型的近似表面積可以用式(10)表示。

式中P為多輪廓模型的近似表面積,cm2;Ci,Ci+1為第i,i+1層投影面輪廓周長,cm,△h 為相鄰投影面的高度,近似等于橢圓臺的母線的長度,cm。

圖3體積和表面積計算示意圖Fig.3 Schematic diagram of volume and surface area calculation
本文在VS2010 開發平臺下采用VC++編程語言,結合PCL(Point Cloud Library1.6 )點云庫和OpenCV 跨平臺計算機視覺函數庫(Open Source Computer Vision Library3.0)對圖像進行處理。圖像分割的目是將目標從圖像中精確的分割出來,對提取的結果進行濾波去噪,提高測量的精度[23]。原始圖像如4a 所示。相機拍攝過程受到白色背景,燈光的影響,圖像存在噪音,需要對灰度化圖像進行濾波處理,本文應用加權平均法將彩色圖像轉換為單色圖像方法進行灰度化[24],計算公式為Gray=0.299R+0.587G+0.114B,式中R、G、B 分別代表彩色圖像3 個基本單色光紅色、綠色、藍色的強度等級,其取值范圍均為0~255,灰度化效果如圖4b 所示。均值濾波法[25]是一種把相鄰像素的相應分量值的平均值作為中心點像素相應分量值的方法,本文采用均值濾波(5×5算子)的方法對紅棗圖像進行濾波操作,濾波效果如圖4c 所示。圖像的二值化就是將圖像的像素點灰度值設置為0 或255,將整幅圖呈現出明顯的黑白效果,本文采用Otsu 算法,該算法進行二值化具有簡單,時效性高,能夠兼顧邊緣細節和分割效果的優點[26],二值化效果如圖4d 所示。二值化分割后,圖像的邊緣存在不連續現象,影響邊緣輪廓的提取,應用開運算法使圖像邊緣光滑,應用腐蝕算法增加連通區域,效果如圖4e f 所示。輪廓邊緣坐標的提取方法分為兩步:1)首先應用Opencv 中的findContours()函數對圖像輪廓進行檢查,檢索模式為RETR_EXTERNAL,該模式只檢索最外層輪廓;模式標識符選CHAIN_APP_ROX_NONE,其作用獲取輪廓邊緣每一個像素的坐標位置;2)應用Opencv中的drawContours()函數繪制外輪廓,函數中參數thickness 來控制輪廓線的粗細,一般設置為1,效果如圖4g所示。輪廓邊緣像素點是二維坐標,需要將二維坐標轉換為三維坐標[26];具體操作如下,構建一個三維結構體CvPoint3D32f(float x,float y,float z),將二維坐標存儲到三維結構體中;再將三維結構體CvPoint3D32f 的數據類型轉化為PointXYZ 點云類型,存儲為.pcd 文件格式[27],效果如圖4h 所示。由于圖像采集過程,不同間隔角度的輪廓坐標都有一定的差異,點云的坐標系不同,需要對點云進行平移到統一的坐標系中,即采用getMin-Max3D 最小包圍盒函數計算輪廓坐標中心(x,y,z),輪廓點坐標中心(x,y,z)與原點坐標(0,0,0)位置相減就實現輪廓平移,效果如圖4i 所示。輪廓以Z 軸為旋轉中心按二維輪廓與三維輪廓空間的映射關系轉換為多輪廓模型[28],效果如圖4j 所示,模型中點云十分密集,存在冗余,無法有效運算,降低了輪廓模型體積和表面積的計算速度,需要對多輪廓點云模型進行下采樣操作。本文采用體素柵格[29]的點云下采樣方法,設置素網格立方體的大小為0.1 cm×0.1 cm×0.1 cm,如圖4k 所示。

圖4圖像分割和處理過程Fig.4 Image segmentation and processing
試驗樣本為阿克蘇駿棗,購買于2018年10月,采用排水法測量紅棗的體積,將紅棗按體積大小劃分成5 個等級,如表1 所示,每個等級取20 個紅棗,共計100 個紅棗樣本備用。
電子數顯游標卡尺(中國桂林量具刃具有限責任公司,精度:0.01mm),直嘴溢水燒杯,量筒(量程:100 mL,精度:1 mL),錫紙,剪刀,記號筆,密封袋,標簽紙,鑷子等。
自主開發一套機器視覺圖像采集裝置,如圖5 所示,試驗裝置由工業相機,計算機,旋轉臺,光源,暗箱等組成。相機采用維視MV-EM510C 型工業相機及4.0 mm 焦距鏡頭,像素尺寸3.45μm×2.2 μm,幀率15 幀/s,曝光時間為30~5 000 000 μs,以太網與計算機相連,相機通過云臺夾固定在暗箱邊框,水平朝向旋轉臺中心,距離被測物體約300 mm;照明光源采用兩個LED 燈;旋轉臺采用寶康隆NA250 電動旋轉臺,轉速60 s/圈,旋轉角度360°,旋轉臺半徑為125 mm,在旋轉臺的中心固定一根長度為160 mm 的鋼針,鋼針表面涂為白色。

圖5圖像采集裝置Fig.5 Image acquisition device
用戶界面用MFC(Microsoft Foundation Class)進行界面設計,該用戶操作采用軟觸發,定時器控制,在紅棗的背后放置白板,減少復雜環境的干擾。鏡頭到紅棗的距離150 mm,主光軸高度300 mm,在圖像采集時,首先將紅棗樣本果蒂部位刺插在針尖上,啟動旋轉按鈕,調試光源和焦距,當樣本清晰呈現在相機鏡頭正前方,觸發相機。30 s內分別得到樣本的圖像的幀數為12,15,20,30,45;間隔角度分別為15°,12°,9°,6°,4°;采集角度范圍為0°~180°,連續采集圖像,圖像的存儲格式為JPG。
本文采用排水法測量紅棗體積,操作過程如下:將蒸餾水倒入直嘴溢口燒杯,直至直口嘴溢口燒杯不漏水為止;用鑷子夾住紅棗浸入燒杯中使水溢出到100 mL量筒中,讀出量筒中水的體積,即紅棗的體積,每個樣本測量3 次,取平均值;再用平均值減去鑷子的體積即為紅棗的體積。本文圖像的處理的方法是用錫紙包裹紅棗表面,用剪刀剪切錫紙包裹的多余部分,將錫紙剝離壓平放置相機鏡頭的黑色A4 紙上,采集俯視圖如圖6所示,圖中白色區域為紅棗表面錫紙圖像。在圖像采集過程中,相機實際成像與理想成像之間存在非線性光學變形,為了減少系統誤差,提高系統測量的精度,需要對系統進行標定。本文參照[30-31]提出的標定方法,提取圖像的特征向量為像素單位,將其轉換為圖像的實際值。在系統硬件參數、物距確定的情況下,圖像的物理尺寸和像素尺寸比值固定不變[32],如式(11)所示。

圖6紅棗表面錫紙圖像采集Fig.6 Tin paper on red jujube surface image acquisition

式中A為比例系數,Lr為物理尺寸,Lp為像素尺寸;實際面積=0.005 6×像素尺寸,實際體積=0.000 42×像素尺寸,應用Matlab中的bwarea函數統計白色區域面積。
將紅棗按體積大小劃分為5 個等級,每個等級20 個紅棗,計算各等級的最小值、平均值、最大值、標準偏差和變異系數。不同等級紅棗的實際體積和表面積統計分析結果如表1所示,體積的標準偏差在1.09~2.05之間,變異系數在6%~12%之間;表面積的標準偏差在2.56~2.88 之間,變異系數在4%~10%之間。

表1 不同等級紅棗實際值統計分析Table 1 Different grades of red jujube actual value statistical analysis
球體的圖像如圖7a所示;提取球體的邊緣輪廓如圖7b所示;將邊緣輪廓像素點坐標位置轉換為三維點云,如圖7c所示;將球體邊緣輪廓的中心偏移到原點位置,如圖7d所示;按映射關系轉換為多輪廓模型,如圖7e所示;將多輪廓模型進行下采樣,圖7f 所示;對多輪廓模型進行垂直Z軸中心點進行投影,如圖7g所示;提取投影面的邊緣點,如圖7h所示。

圖7 球體圖像和點云處理過程Fig.7 Process processing of spherical image and point cloud
球體樣本的的直徑為4.2,3.4,3.0,2.8和2.4 cm,按等間隔角度采集圖像,間隔角度設置為15°、12°、9°、6°和4°,對應采樣圖像幀數12、15、20、30和45,采集范圍為0°~180°。
在表2 中,設球體直徑為3.4 cm 固定值,輪廓間隔角度為15°、12°、9°、6°和4°和投影高度為0.1,0.2,0.3,0.4和0.5 cm。為分析投影高度和輪廓間角度對多輪廓球體模型體積和表面積的計算結果有影響。根據表2,找出多輪廓模型實際值和測量值之間相對誤差的最大值和最小值,來設定合適的投影高度和輪廓間角度,當輪廓間角度為4°,投影高度為0.1 cm 時,體積的最小的相對誤差為6.0%;當輪廓間角為15°時,表面積最小相對誤差為1.0%。在本試驗中,多輪廓球體模型的體積隨輪廓間角和投影高度的相對誤差的增大而增大,表面積的相對誤差隨輪廓間角和投影高度的增大而減?。浑S球體模型輪廓間角的增大,耗時和內存隨之減少,間隔角度4°、6°、9°、12°和15°相對應的消耗時間為180,120,100,60和30 s。

表2 不同投影高度和間隔角度對球體模型的測量誤差Table 2 Measurement error of spherical model with different projection height and contour angle
為分析不同直徑對球體模型體積和表面積的計算結果的影響,在表3 中,輪廓間角和投影高度為固定值,分別為9°和0.1、0.4 cm,直徑為4.2,3.4,3.0,2.8和2.4 cm。觀察表3,多輪廓體積模型體積和表面積的相對誤差隨直徑增大變化不明顯,但球體直徑越小,誤差越大;計算多輪廓球體模型體積和表面積實際值和測量值之間相對誤差的均值,分別為9.1%和4.34%。
通過對多輪廓球體模型的分析,分別確定紅棗模型的投影高度分別為0.1 cm 和0.5 cm,輪廓間角為9°,每個等級紅棗數量為20 個。計算不同等級紅棗模型的體積和表面積。通過觀察表4,多輪廓紅棗模型體積和表面積的均方根誤差隨著等級變化不明顯;體積平均相對誤差隨等級增大而增大,表面積的平均相對誤差隨著等級變化不明顯;不同等級紅棗體積模型的均方根誤差和平均相對誤差的均值分別為2.45 cm3和10.02%;表面積均方根誤差和平均相對誤差的均值為3.65 cm2和7.09%。

表3 不同直徑球體的測量誤差Table 3 Measurement errors of spheres with different diameter

表4 不同等級紅棗體積和表面積的測量誤差Table 4 Measurement error of volume and surface area of different grades of red jujube
本文搭建圖像采集的簡易裝置,按不同角度采集紅棗的圖像,編寫圖像處理軟件,提取紅棗二維圖像輪廓構建多輪廓三維模型,分析在不同輪廓角度,不同投影高度和不同直徑下球體和紅棗三維模型的體積和表面積。
1)當球體直徑為3.4 cm 時,球體的輪廓間角和輪廓投影高度分別為15°,12°,9°,6°,4°和0.1,0.2,0.3,0.4,0.5 cm 時,球體體積的測量值隨輪廓間角和投影高度增大增大,而表面積的測量值輪廓間角和投影高度增大減小,體積和表面積的最大相對誤差分別為6.0%和1.0%。
2)當球體輪廓間角為9°,球體體積和表面積的投影高度為0.1 cm 和0.4 cm,直徑為4.2,3.4,3.0,2.8,24 cm時,不同球體的體積和表面積隨直徑變化不明顯,但直徑越小,體積和表面積的誤差越大,體積和表面積的相對誤差的均值為9.1%和4.34%。
3)當紅棗輪廓間角為9°,體積和表面積的投影高度為0.1 cm和0.5 cm時,紅棗的體積測量平均相對誤差隨等級的增大而增大,表面積隨等級變化不明顯,體積的均方根誤差和平均相對誤差的均值為2.45 cm3和10.2%;表面積的的均方根誤差和平均相對誤差的均值為3.65 cm2和7.09%。
本文采用多輪廓模型測量紅棗的體積和表面積,判斷紅棗過程有一定偏差,需要改進算法來提高測量的精度。上述結論對測量紅棗的體積的表面積的測量提供一種無損的測量方法,為紅棗三維特征的分級裝備開發提供技術參考。