罗小燕,胡 振,汤文聪,刘 占
(1.江西理工大学机电工程学院,江西 赣州 341000;2.江西省矿冶机电工程技术研究中心,江西 赣州 341000)
利用矿石图像进行矿石颗粒识别是机器视觉、图像处理和选矿等领域一个重要的研究方向[1~3]。矿石的粒度分布是评价破碎效果、实现选矿自动化、对选矿设备最优控制的一项重要参数,而矿石粒度检测的关键则是对采集的矿石图像准确识别分割与粒径测量。近年来,基于矿石图像的种类识别方法和技术得到极大的发展,出现了许多关于矿石图像种类识别的方法[4,5]。但由于矿石颗粒的颜色、纹理和形状各有不同,使得同种类矿石的纹理和形状存在着明显的差异性,而不同种类矿石存在一定的相似性,使得现有的识别方法的实际识别效果不明显,这些并不能满足现代矿石粒度检测系统对于矿石检测的高准确性和快速性的要求[6,7]。同时,由于矿石形状是矿石图像的边界信息,而这种边界信息是矿石物种识别的重要特征。
原有方法并不能满足矿石形状检测的需要,因此,需要对原有的方法进行改进。由于矿石边缘的傅里叶(Fourier)描述子具有旋转、缩放和平移等不变性的良好优点,所以被广泛应用于矿石颗粒的种类识别中。然而,这些纹理特征过于单一,识别率并不十分理想[8~10]。局部二值模式(local binary pattern,LBP)具有良好的灰度不变特性,能够克服位移、光照不均等问题[11~13]。近年来,研究者为了提高矿石形状的识别率,通过将矿石颗粒的多种特征相结合进行矿石识别的研究方法,提供了研究矿石形状相应的思路[14~16]。于虎[11]提出了扩展的LBP 算法,提出使用图像梯度信息来代替原始像素,将现有的像素间相关性计算推广到梯度区域间相关性计算,有效地分析出图像区域间关联的更多信息,避免了局部像素轻微变化所引起的扰动。刘念等人[12]将LBP、灰度共生矩阵、傅里叶算子等轮廓特征相结合,提高了矿石识别率,但该方法特征计算时间和训练时间过长,影响了算法识别效率。由于提取矿石特征数据不同,是2组异构数据,为了研究2 组异构数据之间相关关系,本文采用典型相关分析(canonical correlation analysis,CCA)这种多元统计分析算法[17]。CCA 能够揭示2 组特征之间的相关性,提取出的特征在模式分类中与传统的特征提取相比更具有鉴别力[18]。但由于CCA是一种线性相关性,其局限性较大,故本文提出一种基于局部判别CCA 的矿石颗粒分类方法,并在现场采集的矿石图像数据库上进行验证。
矿石颗粒采集系统主要由CCD 摄像机、图像采集系统、挡板、光源、传送带以及计算机等组成。CCD 摄像机设在传送带上方,并位于两挡板中间,与计算机直接相连[19]。可以调整CCD相机高度和摄像头焦距,来改善拍摄范围,再通过MATLAB图像采集工具,完成矿石图像数据实时采集,如图1所示。
通过图像采集系统获得的钨矿石颗粒图像如图2(a)所示,可以看出:矿石颗粒大小不一,形状、颜色各不相同,存在明显的堆积粘连、边缘模糊等问题,因此,在进行颗粒粒径测量前,需要提前解决该影响,减小测量误差。
图2 矿石轮廓边缘线描绘过程
为了方便寻找轮廓特征,提取最终效果图中的边缘轮廓线,得到二值化边缘轮廓图像,如图2(b)。这时可以发现,部分连通域又存在一些极小的子连通域,经过多次试验测量可知,大部分子连通域面积不会超过50 像素,因此,设定连通域面积范围,将其去除,最后标记处理的连通域面积与周长。
通过对破碎实验室的大量矿石进行随机抽样,统计各形态参数,并根据其水平投影的形状特征将其分为两大类,种类1:圆矿石,表示投影形状近似正圆;种类2:椭圆矿石,表示投影形状近似椭圆,如图3所示。
然而,矿石颗粒图像处理系统实现智能分类,则需要根据形状特征引入相关评价参数--形状复杂度因子。形状复杂度因子是描述目标区域的边界复杂度,根据各个矿石形状轮廓、纹理特征得到各个矿石参数值,再以该值进行矿石类别分离,实现矿堆中单一矿石粒径测量。常用的特征描述较多,这里仅介绍傅里叶形状特征描述子、LBP纹理特征描述子两种。
傅里叶描述子是一种图像特征,是一个用来描述目标形状轮廓的特征参数。通过傅里叶变换将物体边界信息转化为形状特征,转换轮廓特征的空间域,提取频域信息特征向量。即用1个向量代表1个轮廓,将轮廓数字化,只需少量的描述子,即可大致代表整个轮廓,从而能更好地区分不同的轮廓,进而达到识别物体的目的。
由于傅里叶描述子具有旋转平移、尺度不变性等特点,为了方便计算,进行归一化处理,其图像轮廓的角度、位置、及轮廓缩放等参数几乎不变,是一个稳定性高的图像特征。因此,在矿石图像中运用傅里叶描述子时,进行如下操作:
1)求取目标矿石图像的梯度幅值大小和方向。选用Sobel算子分别沿x,y方向计算平滑图像的偏导数(Gx,Gy),表达式如下
幅值G(x,y)与梯度方向θ(x,y)表达式如下
2)确定边缘点。遍历整幅图像,从4 邻域梯度方向寻找像素梯度方向的相邻像素,得到边缘点集。计算边缘中心距离。边缘的中心坐标为
以(xc,yc)为极点(xi,yi)(i=0,1,2,…,n-1)将角点直角坐标系转换为极坐标
式中ri(i=0,1,2,…,n-1)为中心角点距离。
3)按照θi的升序对ri(i=0,1,2,…,n-1)排序,得到一组中心角点距离数组,记为V=[r0,r1,…,rn-1],对V进行K点傅里叶变换
由上式构造傅里叶描述子
式中 |g|为傅里叶频谱。
设图像中任意像素点为f(xc,yc),以uc为中心点的3 ×3窗口内的8 个点为u0,u1,…,u7,定义H=h(u0,u1,…,u7)为该局部区域的纹理,对窗口内的非中心像素点做二值化处理,其中,阈值根据中心像素点的灰度值确定。即将图像中3 ×3 窗口内的像素点的灰度值与中间值作比较,其中,小于中心点灰度值的像素用0 来表示,反之则用1 表示,如下所示
按照顺时针方向可以获取一个8 位的二进制数,则该二进制数对应一个二进制模式。中心像素点的特征值则可以由此二进制模式表征。由此每个符号函数可以通过以下公式将其转换为十进制数
用LBP算子[20,21]扫描要进行识别的矿石目标区域,并以直方图的方式来描述目标区域的纹理特征。图像的LBP纹理向量图定义如下所示
式中n为LBP码的数量,3 ×3 模式的LBP码是8 位的二进制数,级数为256。由此获得的LBP 编码包含矿石图像局部模数的边缘和斑点等纹理信息,并将每个LBP码称之为一个微模式。这样可以将一幅矿石图像中的256 个LBP码构成的向量可看成LBP纹理特征。
分别提取矿石图像的傅里叶特征和LBP 特征,训练出2组异构数据。并归一化为X=[x1,x2,…,xn]∈Rp×n,Y=[y1,y2,…,yn]∈Rq×n。利用CCA寻找2 组投影矩阵ωx∈Rp×n和ωy∈Rq×n,使得投影后的矩阵X和Y之间的相关性最大。其目标函数可以表示为
式中Sxx=XXT,Syy=YYT分别为X、Y的协方差;Sxy=XYT为X、Y的互协方差;ωx,ωy称为映射矩阵。
利用CCA求映射矩阵ωx,ωy再用CCA随机测试一幅归一化后的矿石特征图,得到融合后的特征向量。
设计K-mean均值的代价函数为
式中uc(i)为与x(i)最近的聚类中心点。
K-means的另一关键在于K 类别的选择,已经确定矿石颗粒样本类别,则利用K-means 分类器根据其边缘纹理特征将其分为2类(K=2),以便于后续粒径测量。
通过采集平台,本文在16~20,20~30,30~50mm 矿堆中采集矿石样本,为方便计算,3 个粒级各随机抽取8 个矿石,如图4所示。分别计算傅里叶形状特征描述子和LBP纹理特征描述子,并进行特征融合,然后再归一化,将其输入到K-means聚类器,最后可得2组矿石种类。如图5所示。
图5 样本矿石聚类
将特征向量分类结果反馈到矿石分类模型,使矿石图像检测系统扫描样本矿堆,自动识别分类,获取圆矿石、椭圆矿石两大类。最后逐一标记矿石颗粒,并自动描绘边缘轮廓线,其处理结果如图6。
图6 矿石分类
由图6 可知,图中18 个矿石颗粒,大部分矿石颗粒边缘轮廓被描绘为椭圆矿石,识别率为88.88%,仅有2 个被描绘为圆矿石,识别率为12.12%。结合实际生产,破碎后的矿石颗粒形态各异,无明显形状特征,因此,破碎产生的椭圆矿石相对较少,往往是以不规则、类似椭圆的形状存在,导致椭圆矿石个数相对更多。说明分类结果鲁棒性高,具有实际工程意义。
根据不同类别的矿石,可计算各组矿石粒径并统计数据,比较单独使用两种测量方法和混合叠加的方法,从而筛选出最优粒径测量方法,为检测粒度分布做准备。
本文根据矿石颗粒形态参数,包括矿石颗粒大小和形状特征,对采集的矿石图像进行二值化、去分割残留等处理,并描绘出边缘轮廓线,标记矿石个数,为矿石分类奠定基础。建立了颗粒分类模型,引入傅里叶形状特征描述子和LBP纹理特征描述子转化为特征向量并标准化,输入到K-means分类器,达到矿石种类识别的目的。由实验结果表明:将20 个样本矿石颗粒分为2 个圆矿石、18 个椭圆矿石两大类,从而使检测系统根据类别分别测量粒度信息,并自动输出综合粒径参数,实现了自动化测量颗粒参数的目的。