随着数字图像处理技术和模式识别技术的发展，机器视觉已作为一种有效的缺陷检测手段应用于工业生产中。产品表面缺陷检测是机器视觉检测的一个重要部分，基于机器视觉的人造板表面缺陷检测系统能够实现对板材表面缺陷的快速、准确定位和识别。在缺陷检测过程中，由于缺陷面积微小，为了提高检测精度，系统会使用像素在2 000万以上的多个相机对板材表面进行图像采集。刨花板是一种常用的人造板材，是由木材或其他木质纤维材料的碎料胶合而成的，其表面存在细小的刨花，在缺陷检测系统采集到的高清晰度、高分辨率图像中，细小刨花会使图像呈现一定纹理，对图像分割造成干扰，同时噪声对图像分割的影响也更加显著，这些因素常导致无法将缺陷准确地从背景图像中提取出来。

针对复杂纹理背景下目标区域难以提取的问题，基于纹理特征的聚类方法是一种有效的解决手段(于海鹏等，2005；Wang et al., 2006；任宁等，2007；业宁等，2007；尹建新等，2011)，首先对图像各区域的纹理特征进行表征，然后使用聚类算法对具有不同纹理特征的区域进行分类，从而将目标与背景区分开，其关键在于纹理特征的描述和聚类方法的选择。灰度共生矩阵是一种描述纹理特征的常用方法，使用统计特征参数对图像不同区域的纹理进行表征，其抗噪能力强、稳定性好、受图像质量影响小(杨旭，2016)。吴东洋等(2010)使用灰度共生矩阵和K-means聚类方法对木材缺陷进行了识别，Boss等(2012)和白雪冰等(2008)将灰度共生矩阵与FCM模糊聚类方法相结合分别对乳房X光片中的病变区域和木材缺陷进行了识别，还有研究者将AP聚类算法用于基于灰度共生矩阵的纹理特征分析中(杨旭，2016)。

灰度共生矩阵的统计特征参数多达14个，且有些参数计算复杂，以往的通常做法是凭经验选择一部分参数进行研究(尹建新等，2011；吴东洋等，2010；朱蕾，2011；邹丽晖，2007)；但针对不同的纹理，每个参数的表征效果并不相同，因此在对刨花板表面图像进行分析时，需要通过科学的手段对各个特征参数的表征能力进行量化评估，挑选出最能反映图像纹理特点的参数。此外，在进行聚类时，常用的K-means、FCM模糊聚类等方法虽然可以得到较高的聚类精度，但均需要进行迭代计算，运算时间较长(Yedla et al., 2010；Boss et al., 2012)。

针对以上问题，本研究提出一种基于灰度共生矩阵和分层聚类的刨花板表面图像缺陷提取方法。首先将图像划分成若干个窗口，为每个窗口构建灰度共生矩阵，使用Fisher准则和线性相关性对灰度共生矩阵各统计特征参数的表征能力进行度量，选取出分类能力强且相关性弱的特征构成特征向量。然后运用BIRCH分层聚类方法对特征向量集进行聚类，借助均值和特征的统计直方图对聚类结果进行优化。最后根据聚类结果，对原图像各窗口进行标记，提取出缺陷所在区域。

1 灰度共生矩阵的构建 1.1 灰度共生矩阵及统计特征参数

灰度共生矩阵是一种通过研究灰度的空间相关特性来描述纹理的常用方法，反映的是图像灰度分布关于方向、局部邻域和变化幅度的综合信息(Haralick et al., 1973)。灰度共生矩阵共有14个统计特征参数(Haralick，1979)，包括角二阶矩(angular second moment，ASM)、对比度(contrast，CON)、相关(corrlation，COR)、熵(entropy，ENT)、方差(variance，VAR)、均值和(sum of average，SOA)、方差和(sum of variance，SOV)、逆差矩(inverse difference moment，IDM)、差方差(difference of variance, DOV)、和熵(sum of entropy，SOE)、差熵(diffenence of entropy，DOE)、最大相关系数(maximal correlation coefficient，MCC)和2个相关性度量参数(information measures of correlation，IMOC)，其中后3个参数是对前11个参数的进一步描述，并非直接对纹理进行表征，且计算过程复杂，因此后3个参数很少被采用(Ulaby et al., 1986)。对图像进行纹理分析时，首先将图像划分成若干个窗口，为每个窗口分别构建灰度共生矩阵并计算各统计特征参数值，然后通过对各窗口的特征参数值进行比较和分析，将图像分成不同的纹理区域。

1.2 构造因子的取值

在构建灰度共生矩阵时，要根据所要表征的纹理特点，合理设置各构造因子的取值，包括步长d、方向参数θ、灰度级Ng和窗口大小S。刨花板表面图像背景部分纹理明显，非常细小，对比度强，杂乱无章，没有方向性和周期性；而缺陷部分纹理不明显，灰度值比较均一，对比度差。本研究在对刨花板表面图像纹理特点进行分析的基础上，结合试验确定各构造因子的取值。

1.2.1 方向参数θ的选择

方向参数θ表示2个像素点之间的角度关系，当纹理显示出一定的方向性时，主方向上的灰度共生矩阵与其他方向具有较大差异，常用的θ取值包括0°、45°、90°和135°等。由于刨花板表面图像纹理不具有方向性，因此参数θ影响较小。为了减小计算量，缩短计算时间，本研究取θ为0°。

1.2.2 窗口大小S和灰度级Ng的选择

由于刨花板表面图像纹理细小且没有明显的周期性，因此为了能够更好地保留细节、提高精度，应尽可能减小窗口大小S，保留较高的灰度级Ng。但S减小会使图像中窗口数量增多，Ng增加会增大灰度共生矩阵维数，从而引起运算时间增加，故要根据系统的实时性要求在S、Ng取值与计算时间上进行权衡。本研究通过试验来确定S和Ng的取值，在分辨率为512像素×512像素的刨花板表面图像中，当S和Ng取不同值时计算0°灰度共生矩阵的11个特征参数，整幅图像的矩阵构造和特征参数计算时间如表 1所示。对表 1进行分析，综合窗口大小S、灰度级Ng和计算时间3方面因素，本研究选择窗口大小S为16×16，灰度级Ng为16。

1.2.3 步长d的选择

d的取值会影响特征参数的计算结果，本研究通过试验考察步长d对各特征参数的影响情况，最终选择出合适的步长。选择大小为16×16、灰度级为16的正常窗口和带有胶斑、松软、油污、杂物、大刨花5种缺陷的窗口，在步长分别为1~8的情况下构建灰度共生矩阵，计算角二阶矩、熵、对比度、逆差矩、均值和、方差和6个常用特征参数，结果如图 1所示。可以看出，当d小于4时，各特征参数随d的变化比较明显，当d大于4时，大部分特征参数变化比较平稳，且此时各特征参数的离散度较好，有助于进行区分，因此本研究取步长d为4。

2 特征参数的评价和选择

灰度共生矩阵共有14个统计特征参数，窗口的特征向量维度较高，存在信息冗杂和重复表述的问题(余丽萍等，2010)，通常凭经验选择4~6个特征参数作为研究对象；但每个参数对不同纹理的表征效果并不相同，在对刨花板表面图像进行分析时，要根据其纹理特点对各统计特征参数的表征能力进行度量。本研究提出一种基于Fisher准则和线性相关性分析的特征参数选取方法，使用Fisher准则对参数的分类能力进行量化，计算线性相关系数对参数间的相关程度进行判断，最终结合二者结果选出分类性能好、独立性强的特征参数。

2.1 特征参数分类能力评价

Fisher准则是特征选择的有效方法之一(Lu et al., 2014)，其基本思想是鉴别能力较强的特征会使分开的2类类间距离尽可能大，而每类的类内距离尽可能小，Fisher比值越大，说明该特征的分类能力越强。本研究采用单个特征参数的Fisher比值作为判别标准，对参数的分类能力进行排序，筛选出分类能力较强的特征参数。

设将n个样本分为C类w₁，w₂，w₃…w_C，其中第i类w_i包含n_i个样本，x^k、m_i^k和m^k分别表示在第k维上样本x的值、第i类样本均值和所有样本均值。第k维特征Fisher比值的计算公式为：

$ {J_{{\rm{Fisher}}}}\left(k \right) = S_{\rm{B}}^{(k)}/S_{\rm{W}}^{(k)}。$

(1)

其中，$S_{\rm{B}}^{(k)} $和$ S_{\rm{W}}^{(k)}$分别表示第k维特征在样本集上的类间方差和类内方差，计算公式为：

$ S_{\rm{B}}^{(k)} = \sum\limits_{i = 1}^C {} {n_i}n{(m_i^{(k)} - {m^{(k)}})^2}; $

(2)

$ S_{\rm{W}}^{(k)} = \frac{1}{n}\sum\limits_{i = 1}^C {} \sum\limits_{x \in {w_i}} {} {({x^{(k)}} - m_i^{(k)})^2}。$

(3)

选取正常窗口和包含缺陷的窗口共121个，分别按上述确定的构造因子构建灰度共生矩阵。在14个统计特征参数中，由于后3个参数是对前11个参数的进一步描述，因此本研究使用各窗口前11个特征参数作为样本集，计算这11个特征参数的Fisher比值，结果如表 2所示。

2.2 特征参数相关性评价

线性相关系数是一种比较常用的度量特征之间相关性的方法，其反映了特征之间关系的密切程度(谢娟英等，2014)。本研究通过计算线性相关系数来度量各个特征之间的相关程度，筛选出相关系数低的特征。

设ε和η是2个随机变量，Conv(ε, η)表示2个变量的协方差，D_ε和D_η表示2个变量的方差，若Conv(ε, η)存在且D_ε>0、D_η>0，则ε和η的线性相关系数ρ为：

$ \rho = \frac{{{\rm{Conv}}\left( {\varepsilon ,\eta } \right)}}{{\sqrt {{D_\varepsilon } \times {D_\eta }} }}。$

(4)

ρ值反映了2个变量的线性相关程度，ρ的绝对值越接近于1, 则ε和η的相关性越强，当ρ为0时，2个变量无线性关系。

采用上述相同的方法构造样本集，计算每2个特征参数之间的线性相关系数并取绝对值，得到线性相关系数矩阵如表 3所示。

本研究将线性相关系数阈值设为0.85，即2个特征之间的线性相关系数大于0.85则认为有较强的相关性，否则认为相关性较低。分析表 3数据，可将11个特征参数分为5组，每组内各参数相关性较强，而各组间相关性较弱。第1组包括角二阶矩、熵、和熵、差熵、逆差矩，第2组包括对比度、差方差，第3组包括方差、均值和，第4组为方差和，第5组为相关。结合表 2各特征参数的Fisher比值，在每组中选取Fisher比值最高的特征。由于第5组相关特征的Fisher比值过低，因此只从前4组中进行选取，最终选择结果为角二阶矩、对比度、均值和、方差和4个特征参数。

3 基于分层聚类的纹理特征分类

图像中各窗口灰度共生矩阵构建完成后，计算角二阶矩、对比度、均值和、方差和4个特征参数的值作为该窗口的特征向量，所有窗口的特征向量构成样本集，通过对样本集进行分类将正常窗口和含有缺陷的窗口区分开，达到缺陷提取的目的。

聚类分析是一种常用的无监督学习方法，将样本集按某种准则划分成若干类，使相似的样本尽可能归为一类，不相似的样本尽量划分到不同的类中。鉴于人造板表面缺陷检测系统的实时性要求，本研究使用计算复杂度较低、聚类速度更快的BIRCH分层聚类算法，能够在类别数未知的情况下，通过对数据集的单遍扫描完成聚类，算法的时间和空间复杂度都较低(赵玉艳等，2008)。

3.1 BIRCH算法

BIRCH算法是一种基于距离的分层聚类算法(Zhang et al., 1996)，通过构建聚类特征树CF-Tree(clustering feature tree)将样本集聚类成若干个簇，CF-Tree的节点由若干个聚类特征CF组成。给定一个聚类中N个d维样本{X₁，X₂，…，X_N}，CF是对该聚类信息的汇总描述，定义为三元组：${\rm{CF}} = (N, \overrightarrow {{\rm{LS}}}, {\rm{SS}}) $其中，$ \overrightarrow {{\rm{LS}}} = \sum\limits_{i = 1}^N {} {X_i}$是样本的线性和向量，反映了聚类的质心，${\rm{SS}} = \sum\limits_{i = 1}^N {} X_i^2$是各样本的平方和，反映了聚类的直径。CF-Tree是1株高度平衡树，包括分支因子B、L和阈值T 3个构建参数，B表示非叶节点所包含的最大CF数，L表示叶节点所包含的最大CF数，T表示叶节点中簇的半径阈值。例如，B=5、L=4的CF-Tree结构如图 2所示。新样本插入时，通过度量与各CF之间的距离，将样本插入到距离最近的叶节点CF所对应的簇中，如果插入后簇的半径超过了阈值T，则对该簇进行分裂。将所有样本插入到CF-Tree后，就完成了对样本集的聚类，得到的各簇就是聚类结果。如果聚类精度达不到要求或者在构建过程中CF-Tree过于庞大，可以调整参数B、L和阈值T，对CF-Tree进行重构，从而优化聚类结果。

3.2 纹理特征聚类

将分辨率为512像素×512像素的刨花板表面图像分成1 024个16×16窗口，为每个窗口构建0°灰度共生矩阵，计算二阶矩、对比度、均值和、方差和作为窗口的特征向量，1 024个特征向量构成样本集，使用BIRCH分层聚类算法对样本集进行聚类，区分出正常窗口和含有缺陷的窗口。BIRCH分层聚类算法的关键在于CF-Tree构建参数B、L和阈值T的选取以及对聚类结果的优化策略。

3.2.1 CF-Tree构建参数的选取

CF-Tree构建参数B、L直接影响CF-Tree的规模，B、L取值减小会使CF-Tree规模增大。蒋盛益(2009)通过试验确定B、L参数值选择在[2，7]之间聚类效果较好。本研究的样本集属于中等规模大小，因此B、L参数值都设置为5。

阈值T决定了聚类结果中每个簇的半径。T越小，同一个簇中各样本的距离越近，簇的数量越多，聚类精度越高，但是CF-Tree的规模会越大，消耗的时间和内存会越多；如果T过大，则会导致不同类别的样本被划分到同一簇中，达不到聚类精度。

选取带有胶斑、杂物、大刨花和油污缺陷的4张刨花板表面图像，形成4个样本集。当T在[10，40]区间取值时，对4个样本集进行聚类，得到的簇数量如表 4所示。

对4个样本集进行手工分类，得到的簇数量均为2。从表 4可以看出，当T取30和35时，BIRCH分层聚类算法得到的簇数量与手工分类结果最相近。但当T取35时样本集4分类精度不够，当T取30时样本集3分类精度过高，这2种情况都需要对聚类结果进行优化。为了保证足够的分类精度，本研究阈值T设置为30。

3.2.2 聚类结果优化

聚类完成后，可能存在聚类精度不符合要求或CF-Tree过于庞大的情况，需要对聚类结果进行优化。通常的优化方法是对CF-Tree进行重构，增大或减小T后扫描CF-Tree中所有CF结点，进行分裂或合并；但重构过程计算量较大，T的变化量难以精确把控。为了使聚类结果更准确，同时加快计算速度，本研究提出一种优化策略，对样本集的均值和特征进行统计，绘制直方图，根据波峰数量确定样本集的理想类别数量K，当聚类产生的簇数量大于K时，将聚类结果中距离近的簇合并，从而解决分类精度过高而导致的过分割问题。

1) K的计算  均值和特征反映了窗口纹理的灰度特性。在刨花板表面图像中，缺陷部分的灰度变化较小、相对均一，灰度值与正常部分有一定差异，包含缺陷的窗口均值和比较相近，且与正常窗口相差较大。本研究将各窗口的均值和取整，统计各个值对应的窗口数量，绘制样本集的均值和统计直方图。由于缺陷窗口与正常窗口的均值和有较大差距，因此会在统计直方图中不同位置形成波峰，通过统计直方图中的波峰数量可以得出样本集的理想类别数量K。图 3a为带有缺陷的刨花板表面图像，计算其各窗口的均值和，得到统计直方图如图 3b所示，直方图有2个波峰，较大的波峰由正常窗口形成，较小的波峰由缺陷窗口形成。

2) 簇的合并  为了保证足够的分类精度，上述分析中选择了较小的阈值T，但当聚类完成后有可能出现精度过高的情况，即得到的簇数量大于K，需要进行簇的合并，合并过程如下。

① 计算各个簇的质心，计算公式为：

$ {X_o} = \frac{{\sum\limits_{i = 1}^N {} {X_i}}}{N}。$

(5)

式中：N为簇中样本数量；X_i为簇中第i个样本。

② 计算各簇之间的距离，将距离最近的2个簇合并，两簇之间的距离计算公式为：

$ D = {\left[ {{{({X_{o1}} - {X_{o2}})}^2}} \right]^{\frac{1}{2}}}。$

(6)

式中：X_o1和X_o2是2个簇的质心。

③ 如果合并后的簇数量与K相同，算法结束；若簇数量仍然大于K，则返回步骤①再次寻找距离最近的2个簇进行合并。

4 试验 4.1 试验步骤

试验使用的计算机硬件配置为Intel Core I7-7500U四核2.7 GHz CPU处理器，8 G内存，采用Java编程实现灰度共生矩阵构建、特征参数计算和聚类。采集30 cm×30 cm刨花板表面图像，得到512像素×512像素灰度图像，提取缺陷区域，提取过程和参数设置如下：

1) 使用直方图均衡法对图像进行预处理，使其纹理特征更加明显；

2) 将每幅图像划分成1 024个16×16的窗口；

3) 构建每个窗口的0°灰度共生矩阵，构造因子d=4，Ng=16；

4) 计算每个窗口灰度共生矩阵的角二阶矩、对比度、均值和、方差和，构成该窗口的四维特征向量，所有窗口的特征向量构成样本集；

5) 绘制样本集的均值和统计直方图，计算波峰数量K；

6) 使用BIRCH分层聚类算法对样本集进行聚类，得到初始聚类结果，CF-Tree的构建参数B=5、L=5、T=30；

7) 如果初始聚类结果的簇数量大于K，则通过簇的合并对结果进行优化；

8) 统计每个簇中的样本数量，数量最大的簇为正常窗口，其余为缺陷窗口，根据统计结果对原图像中各窗口进行标记，提取出缺陷区域。

4.2 结果与分析 4.2.1 缺陷提取结果

选取5张分别带有胶斑、松软、油污、杂物、大刨花胶斑混合缺陷的刨花板表面图像，使用4种方法进行缺陷区域提取，其中方法1为Otsu阈值法，方法2为吴东洋等(2010)方法，方法3为Boss等(2012)方法，方法4为本研究方法。4种方法的提取结果如图 4所示。

从图 4b可以看出，Otsu阈值法(方法1)无法将缺陷部分与正常部分区分开，其余3种基于灰度共生矩阵的纹理特征聚类方法都可以较好地提取缺陷。方法2使用5个特征参数和K-means聚类，提取出的缺陷面积最小，但存在漏检，松软缺陷提取不完全，没有将胶斑和杂物图像中面积最小的缺陷标识出来。方法3使用14个特征参数和FCM聚类，所有缺陷区域均被标识，且提取完整、边缘平滑，但松软缺陷存在较大面积的误检。本研究方法缺陷区域提取完整，虽然在提取松软缺陷时也存在误检，但较方法3大大减少。漏检往往会比误检造成更加严重的后果，因此在实际应用中更倾向于选择漏检率低的方法。

4.2.2 方法性能分析

1) 分类性能本研究使用精确度(Precision)和召回率(Recall)2个指标评价分类方法性能。精确度反映方法的查准率，召回率反映方法的查全率，指标值越高，表示分类性能越好。采集60张刨花板表面图像，其中带有胶斑缺陷15张、松软缺陷10张、油污缺陷10张、杂物缺陷15张、2种以上缺陷10张，形成5个样本集。设正常窗口为正样本，包含缺陷窗口为负样本，TP表示分类正确的正样本，FP表示分类错误的正样本，TN表示分类正确的负样本，FN表示分类错误的负样本。精确度和召回率的计算公式为：

$ {\rm{Precision = }}\frac{{{\rm{TP}}}}{{{\rm{TP + FP}}}}{\rm{;}} $

(7)

$ {\rm{Recall = }}\frac{{{\rm{TP}}}}{{{\rm{TP + FN}}}}。$

(8)

对各样本图像进行人工缺陷标识作为分类标准，分别使用方法2、方法3和本研究方法对各样本集进行分类，各方法分类结果的精确度和召回率如表 5所示。可以看出，方法2分类精确度最高，召回率最低；方法3分类精确度最低，召回率最高；本研究方法分类精度和召回率都较高。

2) 执行时间基于纹理特征的缺陷提取方法主要包括纹理特征计算和窗口聚类2个步骤。分别使用方法2、方法3和本研究方法对60张刨花板表面图像进行缺陷提取，计算平均每张图像所消耗的时间，结果如表 6所示。方法2、方法3和本研究方法使用的特征参数数量分别为5、14和4，参数数量越多，计算时间越长。在进行聚类时，本研究使用的BIRCH算法聚类时间最短，方法3使用的FCM方法时间最长。

5 结论与讨论

本研究基于灰度共生矩阵和分层聚类方法对刨花板表面图像进行缺陷区域提取, 使用Fisher准则和线性相关性对灰度共生矩阵统计特征参数进行评价，选取出4个分类能力强且相关性低的特征构成特征向量，运用BIRCH分层聚类算法对特征向量集进行聚类，并提出一种基于均值和统计直方图的聚类结果优化策略，以解决聚类精度过高而导致的过分割问题。结果表明，基于灰度共生矩阵和分层聚类的刨花板表面图像缺陷提取方法能够准确将缺陷区域提取出来，精确度达92.2%，召回率达91.8%，解决了在高精度、高分辨率情况下刨花板表面缺陷难以准确检测的问题。

下一步，需要进一步提高方法的实时性，使其能够应用到人造板企业的实际生产中。在人造板连续压机生产线上，刨花板规格为240 cm× 480 cm，为了提高检测速度，可以使用多台相机和计算机进行图像采集和并行处理以加快检测速度。此外，由于缺陷面积占板面面积比例极小，可以先将采集到的图像划分成若干个大小合适的子区域，通过计算子区域的灰度均值或方差对缺陷所在子区域进行定位，再使用本研究方法提取出缺陷具体形状，以便后续对缺陷进行分类识别。