1 引言

核电站燃料组件工作环境恶劣,核燃料棒作为燃料的载体是核电站安全运行的重要环节。

按照核燃料棒的安全运行标准,直径为9.67mm的燃料棒表面的裂缝的安全宽度不应超过0.3mm,加之受水深15m和水下杂质等因素的干扰,准确测量裂缝尺寸的难度较大,但是对核燃料棒表面裂缝尺寸的准确把控对核电站的安全运行具有重要意义。在这种特殊环境下采用非接触测量最为合适,机器视觉检测是非接触测量的典型代表^[1],相比于超声波检测法^[2]、声发射检测法^[3]、光纤传感检测法^[4]等,视觉检测具有更加便捷、直观、效率高的特点,考虑到水下的复杂环境以及核电站的环境问题,采用线结构光视觉检测方法最为合适。

目前广大学者利用视觉检测对物体表面缺陷的处理已经取得了诸多研究成果,但是对于核燃料棒表面缺陷的研究文献较少。卢荣胜等^[5]以智能制造业表面缺陷在线自动检测为应用背景,系统地综述了自动光学(视觉)检测技术。刘丹等^[6]对于复杂背景下微小圆柱体芯块的端面缺陷检测,采用了基于线结构光的三角剖分技术进行了缺陷检测,有很好的检测效果。高军强等^[7]为实现钢轨表面缺陷的实时全方位检测,提出基于多摄像机线结构光的三维检测方法,实现了钢轨表面缺陷的准确检测。Tsai等^[8]对规则纹理表面(机械加工表面、天然木材、纺织面料)的缺陷进行复原,利用高频的傅里叶分量对应表面纹理线型区域,低频的傅里叶分量对应表面缺陷区域,实现了缺陷的检测。李书华等^[9]将投影法、区域生长法、分块法应用到表面缺陷检测中,能够精确检测PVC管材的细小缺陷。Alippi等^[10]利用图像对铁轨的表面质量进行自动检测,利用车载检测系统对铁轨的表面质量进行实时检测和分类。王颖等^[11]通过判断相邻点法线夹角的变化进行管道内表面缺陷的识别,简化了三角剖分过程。刘婷婷等^[12]基于Zernike矩亚像素的高反光金属工件缺陷检测算法,有效地解决了放光问题的干扰。王耀东等^[13]提出了一种基于连通区域的多级滤波算法,结合全局与局部的预处理算法,较好地实现地铁隧道表面裂缝的智能检测。苏俊宏等^[14]针对圆柱型高精密零件高曲率表面缺陷检测的问题,设计了一套基于机器视觉的检测系统,可有效解决圆柱型高精密零件表面缺陷在线检测的问题。李克斌等^[15]采用高、低角度光源组合打光方式提取感兴趣区域,构建划痕形态学的中值滤波核,获取了准确的背景图像,再经背景差分后提取出划痕缺陷。

燃料棒表面裂缝细小,并且相机被置于深水下距离燃料棒65cm处进行拍摄,受水下的杂质和辐射等影响,测量难度大。针对以上问题,搭建了图 1 所示的水下燃料棒模拟测量系统,在距离65cm处利用伺服电机带动测量系统进行旋转扫描拍摄,提出了基于线结构光扫描的视觉测量方法,可以满足水下核燃料棒测量要求,实验可得使用所提方法能简单快捷地测量出燃料棒裂缝尺寸,并且误差保持在0.03mm以内,有效地满足了测量的便捷性和精确性。

图 1. 水下核燃料棒测量系统

Fig. 1. Underwater nuclear fuel rod detection system

下载图片 查看所有图片

2 测量系统模型及标定

基于线结构光的水下视觉检测系统,主要由CCD相机、激光器以及旋转驱动组成,测量系统的数学模型如图2所示,利用激光发射器发出的线结构光平面与待测物体表面相交形成激光线条,利用相机拍摄带有结构光线条的被测物体图像,利用系统与待测物体之间的相对运动完成对被测物体的扫描。测量系统的标定包括相机标定、水下折射校正、线结构光平面标定、旋转中心标定。基于相机拍摄的待测物不同位置的线结构光二维图像,利用相机与激光器之间的固定参数关系,通过旋转平移矩阵将二维的线结构光像素坐标转化为空间三维坐标,以反映被测物体的三维形态。

图 2. 系统测量数学模型

Fig. 2. Mathematical model of mearsured system

下载图片 查看所有图片

系统标定时需要建立从图像二维坐标(u,v)到空间三维坐标(X_w,Y_w,Z_w)之间的转换关系,转换公式为

Zcuv1=fx000fy0u0v01000RT01XwYwZw1,(1)

式中:Z_c为转换矩阵;f_x= fdx、f_y= fdy分别为x轴和y轴上的归一化焦距;图像平面坐标系x O y的原点 u0,v0表示摄像机光轴与图像平面的交点;R= r11r12r13r21r22r23r31r32r33为3×3的单位正交旋转矩阵;T= txtytzT表示平移矩阵。

利用立体靶标法标定激光平面方程和旋转中心^[16],得到激光平面方程为

AXW+BYW+CZW+D=0。(2)

联立 (1) 式和(2) 式实现了二维坐标到空间三维坐标之间的转换,获得了被测物体的(X_w,Y_w,Z_w)空间坐标。

检测设备是利用伺服电机带动摄像机和激光器的平台对待检测物体进行扫描拍摄,所以电机每转一个角度,对应的空间坐标系也随之发生改变,如图3所示。故需要将旋转后的坐标点转换到统一的坐标系下,即以旋转中心点为新的空间坐标系原点,将原坐标系平移,X轴对应的坐标值不变,转换后的坐标为

X'w,Y'w,Z'w=(XW,Yw-Yrcos⁡cosθ+Zw-Zrsin⁡sinθ,Zw-Zrcos⁡cosθ-Yw-Yrsin⁡sinθ,(3)

其中,(X_r,Y_r,Z_r)为旋转中心在原始空间坐标系下的坐标,θ为相机与激光器平台转过的角度。

图 3. 旋转中心坐标转换关系

Fig. 3. Conversion relationship of rotation center coordinate

下载图片 查看所有图片

3 核燃料棒表面裂缝测量算法

3.1 感兴趣区域的定位框选

由于拍摄距离远,图像尺寸为2464pixel×2056pixel,结构光区域大小为40pixel×130pixel。为了观察方便,将燃料棒表面图像裁剪放大为如图4(a)所示的原始扫描图像,根据结构光的特点设定框选的矩形框大小,将待处理结构光框选出来,以避免噪声的干扰和减少计算量,结构光待处理区域的定位过程步骤如下:

1) 根据测量的图像灰度值分布特点,设定一个灰度阈值K,作为灰度值的对比参考值;

2) 对图像进行逐行灰度值遍历,搜索出该行灰度值最大值G_max对应的像素点,当第一次出现G_max≥K时,若该点周围像素点的灰度值均不小于K,则把该点作为区域的起始定位点P_s;

3) 继续对图像进行逐行搜索,搜索出每行灰度值最大值G_max的点,如果G_max<K,则将该点作为待测量区域定位的结束点P_e;

4) 根据所搜索出的P_s和P_e点,计算出待处理区域的中心坐标为

Co=Ps+Pe/2。(4)

5) 根据结构光的像素尺寸大小特点,设定框选矩形的长H、宽W,以C_o为矩形的中心坐标,框选出待处理区域。

该方法能够准确快速地搜索出需要处理的结构光区域,能够有效地将部分噪点隔离在矩形框之外,为之后的图像处理缩小了范围,提高了处理效率,框选效果如图4(b)所示。

图 4. 感兴趣区域的定位框选。(a)原始图像;(b)框选效果图

Fig. 4. Position and frame selection for region of interest. (a) Original image; (b) diagram of frame selection effect

下载图片 查看所有图片

3.2 灰度值水平投影与微分处理

通过分析裂缝的灰度特点,对已经框选定位的图像区域进行灰度值水平投影和一阶微分差值求取。然后将每行值统计到统一的坐标系下,对灰度值变化情况进行分析,最终确定出裂缝区域。

第i行的图像灰度的水平投影值为

Ph=∑j=Co-W2Co+W2xi,j,(5)

式中:x(i,j)为图像中第i行第j列像素点的灰度值。对燃料棒灰度水平投影值中相邻两行的灰度值求取一阶微分差值,可得

Di=Phi-Phi-1,(6)

式中:C_o- W2是矩形区域的左端点对应的像素列值,C_o+ W2是矩形区域的右端点对应的像素列值。因为核燃料棒表面的裂缝区域和正常区域相比,在灰度值上会有很大差异,如图5所示,其中NR为图像中像素行数,HPJ为每行所有点像素灰度值之和,FOD是一阶微分差值,裂缝位置定位步骤如下。

图 5. 燃料棒水平投影与一阶微分结果

Fig. 5. Horizontal projection and first-order differential results for fuel rods

下载图片 查看所有图片

1) 利用水平投影可以确定结构光水平投影值最大时所处的像素行i_h,再利用一阶微分差值可以确定最大值所处的像素行 i_f。

2) 比较i_h与i_f的大小,确定裂缝处在水平投影灰度值最大处的左右区域。

ih≥if,leftih<if,right。

3) 判断裂缝位置。对于右区域:以微分最大值点为基础点,在其左右寻找微分大于等于0的点,直到值小于0为止,如该区域最小行为i_s,最大行为i_m,则裂缝的结束像素行为i_m,因为裂缝近似对称,可以确定裂缝的开始行为(2i_s-i_m)。对于左区域:以微分最小值点为基础点,在其左右寻找微分小于等于0的点,直到值大于0为止,该区域最小行为i_s,最大行为i_m,则裂缝的起始像素行为i_s,可以确定裂缝的结束行为(2i_m-i_s)。

4) 为了更完整地定位裂缝区域,将框选区域适当地扩大。

(1) |i_m-i_s|≥6,则有:右区域中起始行为(2.5i_s-1.5i_m),结束行为(1.5i_m-0.5i_s);左区域中起始行为(1.5i_s-0.5i_m),结束行为(2.5i_m-1.5i_s),如图6(a)所示。

(2) |i_m-i_s|<6,则有:右区域中起始行为(3i_s-2i_m),结束行为(2i_m-i_s);左区域中起始行为(2i_s-i_m),结束行为(3i_m-2i_s),如图6(b)所示。

图 6. 框选区域的扩大。(a) |i_m-i_s|≥6; (b)|i_m-i_s|<6

Fig. 6. Expansion of the frame selection area. (a) |i_m-i_s|≥6; (b)i_m-i_s|<6

下载图片 查看所有图片

3.3 裂缝边界点判定

对已经确定的裂缝区域进行进一步处理,需要确定出裂缝的具体边界点,首先利用阈值分割分割出裂缝区域,进一步缩小裂缝范围;然后,设定种子点类别,根据种子点的特征缩小裂缝区域,判定边界点,并将其分为上下两个边界点集合,具体步骤如下:

1) 对待处理区域进行阈值分类,并对裂缝区域进行灰度统计,如图7所示,可以看出裂缝区域出现了明显的灰度两极化现象。采用最大类间方差法^[17]对区域进行阈值分割,对于灰度值大于限定阈值的点,将合格点用1代替,将不合格点用0代替,其中Grd为灰度区间划分段,Sum为灰度区间内像素点数量,对处理后的裂缝区域可以进一步地缩小区域范围,以裂缝区域第一行的第一个标记为1的点为起点,最后一个标记为1的点为结束点,缩小选定区域的宽度,如图8所示。

图 7. 灰度统计图

Fig. 7. Gray statistics image

下载图片 查看所有图片

图 8. 裂缝局部放大图

Fig. 8. Enlarged view of cracks

下载图片 查看所有图片

2) 对裂缝区域进行像素灰度值遍历,遍历分为上下边界两种边界类型,利用种子点进行判定,如图9所示的种子点。裂缝上边界判断方法为:对裂缝定位区域按照从左往右、从下往上的顺序进行遍历,以当前遍历像素位置为C(i,j)点,且对应标记为1,若(i+1,j)位置的标记为0,则C为上边界点,如图9(a)所示,此时,设定C点标记为2。对于图9(b)~(d)所示的位置C(i,j),且对应标记为1,判断(i+1,j)位置标记是否为2,若是则C点也属于边界点,将上边界点存于集合S_up。与上边界的判断类似,裂缝下边界判断方法为:对裂缝定位区域按照从左往右、从上往下的顺序进行遍历,边界点的判断方法是将图9所示的图上下翻转,将下边界点存于集合S_down。

图 9. 上边界判断方法。(a)边界像素点类型一;(b)边界像素点类型二;(c)边界像素点类型三;(d)边界像素点类型四

Fig. 9. Upper boundary judgment method. (a) Boundary pixel type one; (b) boundary pixel type two; (c) boundary pixel type three; (d) boundary pixel type four

下载图片 查看所有图片

通过以上步骤可以得到裂缝的上下边界点,实现边界点划分和提取,为之后裂缝边界的准确定位提供基础。

3.4 结构光中心提取及裂缝边界点定位

线结构光视觉测量技术的关键在于线结构光中心的提取,这意味着单张图像中只需在裂缝缺陷的上下边界上各获取一个像素点即可。针对该问题,采用了中心线曲线拟合,利用中心线方程与裂缝边界求交集,确定裂缝的准确位置,具体步骤如下:

1) 自适应阈值处理

计算两像素点间灰度梯度值和每行灰度梯度的平均值,表达式分别为

Tij=Gi-Gij-1,(7)

Ai=STWd-1(8)

式中:G_ij是指图像第i行第j列的灰度值;S_T是第i行灰度梯度值T_ij之和;W_d表示线结构光宽度。

阈值G_T和线结构光宽度W_d之间的关系为

GT=1-Ai×WdGim(9)

其中G_im表示每行灰度值的最大值。

通过该方法求取每行阈值,可以保证不同区域的灰度值不存在相互干扰。依赖与像素所在行的灰度特点设定阈值,如图10所示。

图 10. 自适应阈值分割

Fig. 10. Adaptive threshold segmentation

下载图片 查看所有图片

2) 结构光中心提取及曲线拟合

对非裂缝区域进行光条中心提取,记录每行灰度值大于0的起始点S_p和每行灰度值不为0的最后一个点E_p。经过以上的自适应行阈值处理,可以确定梯度值对应的梯度位置P_ij为

Pij=j+0.5。(10)

根据梯度重心法,计算第i行的结构光条中心,表达式为

Ci=∑j=SpEpTij×Pij∑j=SpEpTij。(11)

对待处理区域采用梯度中心法处理,即可得到结构光条纹的中心亚像素坐标, 根据非裂缝区域的结构光条纹中心坐标,将亚像素坐标点用曲线拟合,拟合公式为

y=a0xn+a1xn-1+…+an-1x+b。(12)

(12) 式是关于a₀,a₁,…,a_n-1和b的线性方程。利用最小二乘法建立曲线拟合的数学模型,根据结构光条纹点集合,建立需求目标函数:

Q(a0,a1,…,an-1,b)=∑i=1Na0xn+a1xn-1+…+an-1x+b-y2。(13)

拟合裂缝区域的正常结构光条中心曲线方程的关键在于曲线最高阶次的确定,拟合过程需防止过拟合和拟合不准确的问题。对于同一个实验,结构光能量强度和分布特点是相似的,所以带有裂缝的结构光[图11(a)]与其附近不带有缺陷的结构光[图11(b)]特点类似。以图11(c)给出的裂缝区域的去除位置为例,对图11(b)中同位置进行去除光条处理,可得图11(d)所示。

图 11. 裂缝处结构光中心线的获取。(a)带有裂缝;(b)不带有裂缝;(c) (a)图去除裂缝区域;(d) (b)图去除裂缝区域;(e)正常光条中心;(f)残缺光条中心

Fig. 11. Acquisition of structured light bar center at the crack. (a) With a crack; (b) without cracks;(c) removing cracked area for Fig. 11(a); (d) removing cracked area for Fig. 11(b); (e) normal light bar center; (f) incomplete light bar center

下载图片 查看所有图片

利用梯度重心法提取图11(b)和(d)的结构光线条中心,如图11(e)和(f)所示,对图11(d)提取的光条中心线进行最小二乘拟合,计算曲线拟合误差的平方和,可得

σs=∑i=Co-H2Co+H2y-yi2。(14)

接着利用去除段的拟合曲线值(x,y(x))与图11(e)裂缝段真实结构光中心坐标值(x,y) ,计算误差平方和,得

σ0=∑i=ibiey'i-yi2,(15)

其中,i_b为裂缝段起始行,i_e为裂缝段结束行。采用不同阶次的曲线进行最小二乘拟合,计算结果如表1所示。

观察发现,σ_s,σ_o在采用8阶拟合函数拟合时总体误差趋于平缓,当拟合阶次达到12时,出现了过拟合现象,拟合误差较大,且误差值上下波动。所以考虑到计算的精度和复杂度,最终确定拟合曲线的最高拟合阶次为8,拟合的曲线如图12所示。

3) 利用拟合曲线方程在裂缝区域的坐标值,计算裂缝区间内曲线在每行处的坐标值并取整,即

y_i= floory(i),ko≥0.5ceilingy(i),ko<0.5,

其中,y_i为曲线在第i行处的取整值,k_o为y(i)的小数部分,floor[·]为向上取整,ceiling[·]为向下取整。

若计算出的上边界点多于一点,如图13(a)所示的点p₁,p₂,然而p₂点与曲线点的误差值小于p₁,则p₂处曲线上的点(i,y(i))为裂缝的上边界点,下边界处理方法相同,中心线与边界相交示意图如图13(b)所示。上下边界点的判别公式为

(i,yi)∈Sup⋃Sdown,boundarypoint(i,yi)∉Sup⋃Sdown,noboundary。

图 12. 结构光中心曲线拟合

Fig. 12. Structured light center curve fitting

下载图片 查看所有图片

图 13. 裂缝边界处理。(a)边界判定示意图;(b)中心线与边界相交示意图

Fig. 13. Crack boundary treatment. (a) Diagram of boundary decision; (b) intersection diagram of center line and boundary

下载图片 查看所有图片

3.5 裂缝尺寸计算方法

利用标定好的旋转平移矩阵将裂缝边界像素的二维坐标转化为三维空间坐标,再将得到的三维坐标点集分为上边界裂缝点集P_u以及下边界裂缝点集P_d,由于裂缝倾斜,测量出来的裂缝宽度其实是倾斜裂缝的斜线长度,如图14所示,点A与点B是同一条光线上的上下边界点,但是点A与点B的距离D_AB大于点A与点C的距离D_AC,直接按照同一条光线上的点计算裂缝宽度会加大裂缝尺寸,导致计算结果偏大,故可将问题转化为计算两点集P_u与P_d间的最短距离。计算两点集间的最短距离最常用的方法就是穷举求解,但是穷举法耗时长,计算空间复杂度较大。本文提出基于距离采样法减小计算量,根据点集P_u与P_d之间裂缝点的一一对应关系,计算出每条线条上裂缝上下边界的距离,表达式为

Di=(xui-xdi)2+(yui-ydi)2+(zui-zdi)2,(16)

其中,(x_ui,y_ui,z_ui)为上边界裂缝点集P_u的第i个点P_ui,以P_ui为球心,以D_i为半径划分出一个采样球形范围如图14所示,区域内点的坐标取值范围为

x∈ xui-Di,xui+Di,

y∈ yui-Di,yui+Di,

x∈ zui-Di,zui+Di。

下边界点集在此范围以内的点为点P_ui的局内点N_ui,计算P_ui与局内点N_ui的距离,再取距离最小值代表裂缝在上边界点P_ui处的实际宽度D_ui,计算出的基于上边界的平均宽度A_u为

Au=∑i=0nDuin。(17)

同理可以计算出下边界在每个边界点P_di处的实际宽度D_di及下边界裂缝平均宽度A_d,计算燃料棒裂缝平均宽度A_ud,即

Aud=Au+Ad2。(18)

图 14. 距离采样图

Fig. 14. Diagram of distance sampling

下载图片 查看所有图片

4 实验部分

为了验证本文算法的适用性,采用win7 64位系统,使用可跨平台的QT5.8.0和VS2013作为开发工具进行程序编写,开发了核燃料组件的测量平台,完成了对核燃料棒的在线扫描测量。该平台的运行流程如图15所示。

图 15. 平台运行流程图

Fig. 15. Platform operation flowchart

下载图片 查看所有图片

1) 运行程序需要先连接好相机和电机,并设置好相关参数。电机旋转带动相机及激光器旋转,通过旋转拍摄实时获取图像。设定每20ms获取一次当前位置的扫描图像。

2) 基于相机实时拍摄的图像,使用文中设计的裂缝测量算法对图像进行处理,采用多线程运行,使图像处理和点云生成同时进行,互不干扰,计算燃料棒表面裂缝的尺寸,并反映出裂缝所处的位置和形态。

图16给出核燃料棒的水下实验平台。本文基于此平台进行了水下测量实验。在正常光照条件下,设置相机曝光度为200、300、400、500,依次拍摄同一处裂缝,所得图像如图17(a)~(d)所示,可以发现随着曝光度的提高,裂缝区域与正常区域的区分度越来越低,故需要根据环境光照条件选择适合的曝光度。图17(e)~(g)为在曝光度为200时拍摄不同尺寸裂缝所得的图像。

图 16. 水下实验装置

Fig. 16. Underwater experimental device

下载图片 查看所有图片

图 17. 相机拍摄图片。(a)曝光度为200;(b)曝光度为300;(c)曝光度为400;(d)曝光度为500;(e)裂缝尺寸为190.0μm;(f)裂缝尺寸为181.6μm;(g)裂缝尺寸为 128.6μm

Fig. 17. Photos taken. (a) Exposure of 200; (b) exposure of 300; (c) exposure of 400; (d) exposure of 500; (e) crack size of 190.0μm; (f) crack size of 181.6μm; (g) crack size of 128.6μm

下载图片 查看所有图片

图18(a)给出核燃料棒表面裂缝。对燃料棒表面三处裂缝进行测量,每处裂缝转换三种不同角度测量三次,得到如图18(b)所示的带有裂缝的燃料棒扫描点云。表2给出对核燃料棒表面的三处裂缝的测量结果,从中可以看出本文算法的测量误差在0.03mm以内,测量精度高,符合精度要求。由此可知,文中提出的基于线结构光燃料棒表面裂缝测量算法能够很好地解决裂缝测量的难题,同时实验发现该方法能够有效地避免在线测量的水下干扰问题,测量系统具有稳健性。

图 18. 核燃料棒表面裂缝。(a)燃料棒表面裂缝;(b)燃料棒表面裂缝点云图

Fig. 18. Cracks on the surface of nuclear fuel rods. (a) Crack on the surface of fuel rod; (b) point cloud of crack on fuel rod surface

下载图片 查看所有图片

5 结论

针对核燃料棒表面细小裂缝测量问题,提出了一种基于结构光扫描的视觉测量方法。该方法通过灰度二维搜索、水平投影和一阶微分,实现了裂缝区域的准确定位,提高了后续处理的计算效率和准确率。然后,利用裂缝区域阈值分割法、种子点判断法,将边界像素点准确地提取和分类。最后,使用最小二乘法拟合结构光中心线,用中心线方程与裂缝边界点集求交点,确定裂缝具体位置,再将裂缝边界点转为三维坐标进行距离采样,计算实际尺寸。实验证明,使用该方法实现了水下65cm处燃料棒表面裂缝尺寸的测量,测量误差保持在0.03mm以内,有效地满足了测量的便捷性和精确性,对核电站的安全运行具有重要意义。