计算机应用   2017, Vol. 37 Issue (10): 2895-2898  DOI: 10.11772/j.issn.1001-9081.2017.10.2895
0

引用本文 

申屠理锋, 奚嘉奇. 基于帧间灰度变化分析的在线光源位置计算[J]. 计算机应用, 2017, 37(10): 2895-2898.DOI: 10.11772/j.issn.1001-9081.2017.10.2895.
SHENTU Lifeng, XI Jiaqi. On-line light source position calculation based on inter-frame gray-scale variation analysis[J]. Journal of Computer Applications, 2017, 37(10): 2895-2898. DOI: 10.11772/j.issn.1001-9081.2017.10.2895.

基金项目

国家科技支撑计划项目(2015BAF22B01)

通信作者

申屠理锋, E-mail: shtulf@baosteel.com

作者简介

申屠理锋(1979—), 男, 浙江桐庐人, 高级工程师, 硕士, 主要研究方向:机器视觉、自动控制;
奚嘉奇(1989—), 男, 上海人, 工程师, 硕士, 主要研究方向:机器视觉、自动控制

文章历史

收稿日期:2017-04-25
修回日期:2017-07-07
基于帧间灰度变化分析的在线光源位置计算
申屠理锋, 奚嘉奇    
宝山钢铁股份有限公司研究院, 上海 201900
摘要: 针对机器视觉系统在实际生产中光源位置无法事先确定的问题,提出一种基于帧间特征区域灰度形态分析的在线光源位置计算方法。首先基于灰度分布来确定一特征区域作为参考点,然后运用块匹配算法确定相邻两帧中特征区域的位置变化,之后结合光照模型建立灰度和几何位置之间的关系,运用线性回归对联立方程组进行求解,最后得到光源位置。实验结果表明:光源位置的计算结果与实际所测的距离相比,误差在5%之内。所提算法已被应用于实际生产,具有较好的计算精度和实时性。
关键词: 帧间    灰度    光源    位置计算    
On-line light source position calculation based on inter-frame gray-scale variation analysis
SHENTU Lifeng, XI Jiaqi     
Research Institute, Baoshan Iron & Steel Company Limited, Shanghai 201900 China
Abstract: Aiming at the problem that the position of the light source cannot be determined in practical production, an on-line position calculation for the light source was proposed, which fully analyzed the gray-scale variation of the feature region in consecutive frames. In this approach, a featured region was first selected based on the local gray-scale distribution and used as a reference point, then its position was tracked in the following consecutive frames via block matching method. Combing with the analytic calculation from the light transport model, the relationship between the gray-scale value and geometric position was established. Finally, the set of equations were solved by linear regression method to find the position of light source. The experiment results show that the errors between the calculated results and the measured results are within 5%. The proposed approach has been successfully applied in the real-time manufacturing with good estimation accuracy.
Key words: inter-frame    gray-scale    light source    position calculation    
0 引言

在钢铁制造业中表面缺陷检测是质量控制的一个重要环节, 这一工作大都由人工操作完成, 随着图像技术的发展, 许多基于机器视觉技术的系统被集成到生产线来辅助或替代人工作业, 而表面缺陷检测算法是这类系统的最核心部分。目前主要包括三类算法[1]:基于线性及非线性滤波、基于空域-频域的变换及重构、基于机器学习的特征提取和分类。在基于滤波的算法中, 首先需根据缺陷形态确定滤波算子, 常用的算子有Laplacian和Gabor, 然后调节算子参数来实现缺陷特征的提取[2-3]。由于参数的调节往往对光照环境的变化较为敏感, 因此这类算法的鲁棒性有所欠缺。为了克服这一缺点, 文献[4]提出了一种基于双光源调节的识别方法。在光照条件给定的情况下, 通过光源的切换来消除非均匀光照的影响, 之后通过设计有限脉冲滤波器来进行缺陷定位。在基于空域-频域变化的这类方法中, 首先需依据缺陷及光照的特点选择核函数, 然后通过优化重构系数来实现前景和背景的分离[5-6]。由于增加了参数的空域维度, 所以与滤波算法相比, 其鲁棒性有一定改善。近些年也出现了许多基于机器学习的方法, 包括运用非监督分类和卷积神经网络进行缺陷检测[7-9]。其主要优势在于通过对样本的自学习, 可自适应地辨认缺陷的特点和光影的影响而无需人工干预。其难点则在于样本的选择、训练的次序及算法的收敛速度。

可以看出, 这些算法都在一定程度上实现了对图像中缺陷的自动检测和识别功能, 但对于如何克服光照环境变化带来的影响, 还缺乏有效的方法。文献[10]提出了一种基于成像过程分析的算法, 通过对坯料表面光照分布的计算, 有效去除了非均匀光照的影响, 之后根据特征分类来识别缺陷。然而在这一算法的实施过程中, 需要事先知道光源所在的位置, 在部分特殊工况下, 这一信息是很难获得的, 同时很多生产过程中, 光源相对于工件表面的距离和位置是不断变化的, 因此一定程度上限制了该算法的应用领域。为此, 本文提出了一种在线光源位置计算的方法来改进这一不足。通过比对多帧连续图像特征区域的像素信息变化, 并结合光照模型分析, 便可以准确地计算出光源相对于圆坯所在的位置。从实验结果来看, 误差控制在5%之内, 包括水平误差和垂直误差两部分。将这一算法与文献[10]算法相结合, 便可以进一步增强算法的通用性。同时由于本文算法在计算时采用了累计迭代的表示方式, 使得整体处理时间严格控制在5 s内, 因此保留了原算法的实时性。

1 在线光源位置计算算法

本文算法的基本思路是, 首先找到一帧图像中的特征区域, 即灰度发生陡变的区域。为了保证特征区域能在连续多帧图像中被稳定地辨识, 要求灰度变化大于等于50%,同时所占的大小为100个像素点。就圆坯图像而言, 典型的特征区域可以是表面缺陷处或者是氧化鳞屑覆盖的地方。在确定了特征区域之后, 接下来便要跟踪这一区域在连续几帧图像中位置的变化情况(轨迹), 考虑到速度及鲁棒性, 本文采用的是基于块匹配的算法。在得到了位置信息之后, 利用特征区域在不同位置点灰度的整体性变化趋势, 并结合光照模型的理论计算, 来推算出光源所在的位置。从计算结果来看, 它非常接近实测距离。

1.1 特征区域运动计算

对于刚体运动而言, 由于其局部和整体的运动状态具有一一对应关系, 而且可以由平移和旋转矢量完全描述, 因此特征区域跟踪这一问题等价于运动计算。对于这一类问题的求解, 通常有两类算法。第一类是基于相关性的算法[11], 另一类则是基于特征的算法[12]。对于当前的问题, 考虑到物体在帧间的运动距离较小, 而且光强分布呈连续性变化, 基于块匹配的算法较为合适。

首先根据特征区域的形态来确定块的大小, 如图 1(a)所示, 然后需要确定在下一帧图像中的搜索区域的大小, 在此不妨将其设为以参考点p0为圆心, 半径为R的圆, 其中R的值可通过对物体实际运动速度的测算而得到。在这一区域中, 本文采用菱形算法进行搜索[13],目标是在这一搜索区域中找到与参考点相匹配的块, 即对应的位置p1图 1(b)表示了这一关系, 其中p0的位置用虚线表示以区分。为了定量的反映匹配程度, 引入了匹配度函数, 表达式如下:

$ \mathit{Corr} = \sum\limits_{\left( {x, y} \right) \in M} {{{\left( {f\left( {x, y} \right) - g\left( {x, y} \right)} \right)}^2}} $ (1)
图 1 特征区域示意图 Figure 1 Schematic diagram of characteristic area

其中M为块的大小。Corr是对块内的每一个像素点作差求平方后的累加和,其值越小则匹配度越好。

1.2 光源位置计算

一般而言, 由于光强的变化, 特征区域的像素值在不同帧中会有一定差异, 而且这一差异往往具有趋势性, 主要依据光源和相机的位置而变化。接下来将利用这种差异来推算出光源的位置。在此之前, 将先简单介绍一下光照模型。在一般成像过程中, 光线到达物体表面后, 一部分被物体吸收, 另一部分则通过反射到达相机镜面, 其中反射的光线又包括漫反射和镜面反射两部分。这里主要讨论漫反射部分。根据朗伯余弦定理[14], 漫反射的入射光强与反射光强有如下关系:

$ f = {\rm{max}}\left( {k\frac{L}{{{r^2}}}\left( {{\mathit{\boldsymbol{u}}_i} \cdot {\mathit{\boldsymbol{v}}_i}} \right){\rm{, }}0} \right) $ (2)

其中:L是入射光强; f是反射光强; r是光源与被测物体间的距离; k是漫反射系数; uivi分别是物体表面pi处的法向量和入射向量。

为了便于接下来的讨论, 考虑如图 2所示的测量情形。一面光源以一定角度对着被测物体, 一CCD相机在另一侧进行图像采集, 在被测物体上拟跟踪的特征区域位于物平面S。在检测过程中, 特征区域(以圆斑表示)随着物体的运动而运动, 形成了一系列轨迹p0p1p2等。pi为其在ti时刻的所处的位置。特别地, 相邻两点pipi-1具有以下关系:

$ {\mathit{\boldsymbol{p}}_i} = {\mathit{\boldsymbol{p}}_{i - 1}} + {\mathit{\boldsymbol{d}}_i};\;i \ge 1 $ (3)
图 2 特征区域运动轨迹示意图 Figure 2 Schematic diagram of feature region trajectory

其中: di即为上一步所得到的运动矢量。本文约定t0时刻所在的位置为参考原点, 即p0=0并建立直角坐标系。我们希望得到光源相在此坐标系中的位置, 即(xs, ys, H)。xsys为光源在物平面的投影, H为到物平面的垂直距离。为了实现这一目的, 先考察光强在不同位置的分布情况。如果将特征点在ti时刻的位置pi代入式(2), 将得到:

$ {f_i} = k'\frac{{{\rm{cos }}{\theta _i}}}{{{H^2} + {{({x_s} - {x_i})}^2} + {{\left( {{y_s} - {y_i}} \right)}^2}}};k' = k \cdot L $ (4)
$ {\rm{cos }}{\theta _i} = \frac{H}{{\sqrt {{H^2} + {{({x_s} - {x_i})}^2} + {{\left( {{y_s} - {y_i}} \right)}^2}} }} $

其中:θipi处法向量ui和入射向量vi所成的夹角。

假定它们处于同一平面, 因此将pi展开成(xs, ys)的形式。对此作进一步简化, 可以得到

$ {f_i} = k'\frac{H}{{{{\left[{{H^2} + {{({x_s}-{x_i})}^2} + {{\left( {{y_s}-{y_i}} \right)}^2}} \right]}^{1.5}}}} $ (5)

为了更清晰地描述ti时刻和t0时刻光强的相对关系, 将两者的值作比可以得到:

$ \frac{{{f_i}}}{{{f_0}}} = {\left[{\frac{{{H^2} + x_s^2 + y_s^2}}{{{H^2} + {{({x_s}-{x_i})}^2} + {{\left( {{y_s}-{y_i}} \right)}^2}}}} \right]^{1.5}} $ (6)

式(6) 的左边即为相对变化量, 对此式进行化简, 并将(fi/f0)开1.5次方后记作αi, 0, 可以得到

$ \left( {1 - {\alpha _{i, 0}}} \right)\left( {{H^2} + {x_s}^2 + {y_s}^2} \right) = {x_i}^2 + {y_i}^2 - 2{x_s}{x_i} - 2{y_s}{y_i} $ (7)

对轨迹上所有的点依次作这般处理, 可以得到一组表达式, 共计n个式子:

$ \left\{ \begin{array}{l} \left( {1 - {\alpha _{1, 0}}} \right)\left( {{H^2} + {x_s}^2 + {y_s}^2} \right) = {x_1}^2 + {y_1}^2 - 2{x_s}{x_1} - 2{y_s}{y_1}\\ \left( {1 - {\alpha _{2, 0}}} \right)\left( {{H^2} + {x_s}^2 + {y_s}^2} \right) = {x_2}^2 + {y_2}^2 - 2{x_s}{x_2} - 2{y_s}{y_2}\\ \vdots \\ \left( {1 - {\alpha _{n, 0}}} \right)\left( {{H^2} + {x_s}^2 + {y_s}^2} \right) = {x_n}^2 + {y_n}^2 - 2{x_s}{x_n} - 2{y_s}{y_n} \end{array} \right. $ (8)

可以发现表达式(8) 的左边都有着类似的结构, 公共因子为H2+xs2+ys2, 因此如果两两作比, 可以得到以下n-1个式子:

$ \left\{ \begin{array}{l} \frac{{\left( {1 - {\alpha _{1, 0}}} \right)}}{{\left( {1 - {\alpha _{2, 0}}} \right)}} = \frac{{{x_1}^2 + {y_1}^2 - 2{x_s}{x_1} - 2{y_s}{y_1}}}{{{x_2}^2 + {y_2}^2 - 2{x_s}{x_2} - 2{y_s}{y_2}}}\\ \vdots \\ \frac{{\left( {1 - {\alpha _{n - 1, 0}}} \right)}}{{\left( {1 - {\alpha _{n, 0}}} \right)}} = \frac{{{x_{n - 1}}^2 + {y_{n - 1}}^2 - 2{x_s}{x_{n - 1}} - 2{y_s}{y_{n - 1}}}}{{{x_n}^2 + {y_n}^2 - 2{x_s}{x_n} - 2{y_s}{y_n}}} \end{array} \right. $ (9)

由于xiyipi的分量, 而pi又可以由di来表示, 因此均为已知量。所以此方程组为线性方程组, 利用简单的线性回归方法, 便可求得关于xsys的解。将所得的解代入式子中, 可求得H, 至此光源的位置便确定了下来。

1.3 整体表达式

对以上两个步骤作适当整理, 便可以得到完整的算法流程如图 3所示。

图 3 本文算法流程 Figure 3 Flowchart of the proposed algorithm

首先, 在第一帧中确定拟跟踪的特征区域以及它所在的位置(x0, y0)。对于此后采集到的每一帧, 将其与前一帧作比较, 通过块匹配算法来得到特征区域的位移矢量di, 并以第1帧为参考计算特征区域的相对灰度比αi, 0, 然后将其分别代入关系式(9) 的两边。以此类推, 当得到了n-1组表达式之后, 便可以运用线性回归对其求解, 最终得到光源的位置。由于在一般工业检测环境中, 光源位置一般固定, 因此算法只需在初始化中作一下自动纠正即可, 并不需要每步执行。

2 实验与结果

本文所提出的方法已成功应用于钢厂圆钢磁粉探伤表面缺陷检测系统中。磁探是钢铁制造业中一种常用的表面检测手段, 其基本原理利用的是漏磁检测技术, 被测工件被磁化后, 工件缺陷处的磁力线会发生局部畸变而产生漏磁场, 当磁悬液喷洒在工件上后, 漏磁场会吸附更多磁粉, 其将会辐射更多光能, 从而显示出缺陷的位置、形状和大小。

下面将用一组实验数据来验证本文算法的实际计算精度。实验工序如图 4所示。

图 4 实验工序示意图 Figure 4 Schematic diagram of experiment process

圆坯在一组电机的驱动下呈螺旋运动, 经过磁化以及喷洒磁悬液后, 来到检测区域。在区域的两侧置有两盏紫外灯, 以一定的角度照射着坯料表面, 在坯料的正上方安装了一台高速相机, 其采集的图像通过以太网传至机柜内的工控机, 后者在对图像进行分析处理后将最终结果显示于显示屏。为了保证清晰成像, 尽量避免因坯料运动造成的拖影, 根据实际生产现场的工况条件, 选用焦距为24 mm的工业相机, 结合坯料与相机的相对运动速度设定相机曝光时间为25 ms, 采样帧率设定为12.5 Hz。

首先需要确定待跟踪的特征区域及其位置。在此, 选择了一表面鳞屑区作为跟踪对象, 并用方框标记了其所在的位置,如图 5(a)所示。图 5(b)(c)是之后采集的第5帧和第10帧的图片, 帧间隔为80 ms, 可以看出特征区域从左上角移动到了右下角, 同样用框标记其位置。

图 5 帧间图像序列 Figure 5 Inter-frame image sequence

图 6显示的是特征区域在这3帧图像中的灰度形态, 可以发现尽管存着局部的灰度差异, 但整体形态分布非常相似, 因此利用块匹配算法, 将很容易对其进行识别, 经计算平均匹配度大于95%。有了位移矢量之后, 接下来便可利用帧间存在的局部明暗变化来计算光源所在的位置。在对连续10帧图像进行计算后, 最终光源的位置稳定在0.75 m附近, 与实际所测的距离相比, 误差在5%之内。由于运用了迭代求解的方式, 一次计算时间仅为5 s。按照这一方法又对不同规格的坯料下进行了25组实验, 具体运行结果如表 1所示, 可以看出其能较好地适应不同工况环境。

图 6 特征区域灰度形态分布 Figure 6 Gray morphological distribution of feature region
表 1 本文算法运行结果 Table 1 Running results of the proposed algorithm
3 结语

为了解决机器视觉系统在实际生产中光源位置无法确定的问题, 本文提出了一种在线光源位置计算的方法。通过比较连续多帧图像中特征区域的灰度明暗变化情况, 并结合光照模型的理论计算, 可以建立灰度与几何位置之间的关系, 之后利用线性回归的方法可高效地求解出光源位置。本文算法已被应用于钢厂圆坯磁粉探伤现场自动缺陷识别系统中, 用于计算光源位置。实验结果表明:光源位置的计算结果与实际所测的距离相比, 误差在5%之内。其速度和计算精度已得到验证, 完全满足一般工业实际生产中图像测量系统对光源位置的需求。

参考文献(References)
[1] NEOGI N, MOHANTA D K, DUTTA P K. Review of vision-based steel surface inspection systems[J]. EURASIP Journal on Image and Video Processing, 2014, 2014: 50. DOI:10.1186/1687-5281-2014-50
[2] YUN J P, CHOI S, KIM J W, et al. Automatic detection of cracks in raw steel block using Gabor filter optimized by univariate dynamic encoding algorithm for searches (uDEAS)[J]. Non Destructive Testing & Evaluation International, 2009, 42(5): 389-397.
[3] 吴秀永, 徐科, 徐金梧. 基于Gabor小波和核保局投影算法的表面缺陷自动识别方法[J]. 自动化学报, 2010, 36(3): 438-441. (WU X Y, XU K, XU J W. Automatic recognition method of surface defects based on Gabor wavelet and kernel locality preserving projections[J]. Acta Automatica Sinica, 2010, 36(3): 438-441.)
[4] JEON Y J, CHOI D C, LEE S J, et al. Steel-surface defect detection using a switching-lighting scheme[J]. Applied Optics, 2016, 55(1): 47-57. DOI:10.1364/AO.55.000047
[5] TSAI D M, HSIAO B. Automatic surface inspection using wavelet reconstruction[J]. Pattern Recognition, 2001, 34(6): 1285-305. DOI:10.1016/S0031-3203(00)00071-6
[6] GUAN S. Strip steel defect detection based on saliency map construction using Gaussian pyramid decomposition[J]. Transactions of the Iron & Steel Institute of Japan, 2015, 55(9): 1950-1955.
[7] VEITCH-MICHAELS J, TAO Y, WALTON D, et al. Crack detection in "As-Cast" steel using laser triangulation and machine learning [C]//Proceedings of the 13th Conference on Computer and Robot Vision. Piscataway, NJ: IEEE, 2016: 342-349.
[8] 吴家伟, 严京旗, 方志宏, 等. 基于AdaBoost改进算法的铸坯表面缺陷检测方法[J]. 钢铁研究学报, 2012, 24(9): 59-62. (WU J W, YAN J Q, FANG Z H, et al. Surface defect detection of slab based on the improved AdaBoost algorithm[J]. Journal of Iron and Steel Research, 2012, 24(9): 59-62.)
[9] SOUKUP D, HUBER-MÖRK R. Convolutional neural networks for steel surface defect detection from photometric stereo images [C]//Proceedings of the 10th International Symposium on Visual Computing. Berlin: Springer, 2014: 668-677.
[10] XI J Q, SHEN-TU L F, HU J K, et al. Automated surface inspection for steel products using computer vision approach[J]. Applied Optics, 2017, 56(2): 184-192. DOI:10.1364/AO.56.000184
[11] 禹晶, 苏开娜. 块运动估计的研究进展[J]. 中国图象图形学报, 2007, 12(12): 2031-2041. (YU J, SU K N. A survey of block-based motion estimation[J]. Journal of Image and Graphics, 2007, 12(12): 2031-2041. DOI:10.11834/jig.20071201)
[12] 钟平, 于前洋, 金光. 基于特征点匹配技术的运动估计及补偿方法[J]. 光电子·激光, 2004, 15(1): 73-77. (ZHONG P, YU Q Y, JIN G. Motion estimation and motion compensation based on matching technology of feature point[J]. Journal of Optoelectronics·Laser, 2004, 15(1): 73-77.)
[13] 刘海峰, 郭宝龙, 冯宗哲. 用于块匹配运动估值的正方形-菱形搜索算法[J]. 计算机学报, 2002, 25(7): 747-752. (LIU H F, GUO B L, FEN Z Z. A square-diamond search algorithm for block motion estimation[J]. Chinese Journal of Computers, 2002, 25(7): 747-752.)
[14] PHONG B T. Illumination for computer generated pictures[J]. Communications of the ACM, 1975, 18(6): 311-317. DOI:10.1145/360825.360839