一种基于二维数字图像相关法的试件质量检测方法

本发明属于光测力学领域，具体涉及一种基于二维数字图像相关法的试件质量检测方法。

背景技术：

数字图像相关法(digitalimagecorrelation,dic)是一种基于数字图像处理技术的光测力学方法，通过形函数对变形前和变形后的图像子区进行相关搜索和计算来获取目标区域中每个点发生的位移，因其只依赖于简单的测试环境、能够获得全场变形信息等优点，目前dic已成功应用于多个不同的领域。

目前，二维的dic技术在计算精度、搜索速度方面具有很好的效果，然而在生产过程中仅依靠二维的dic技术是无法满足工业测量需求的，如在曲面试件表面变形的全场测量中，只能使用三维的dic技术才能实现，但是三维的dic技术在计算精度、搜索速度方面还有待提高。目前针对三维的dic测全场变形，更多的是使用三维重构技术，要求至少使用两个相机从不同的角度采集不同形变时刻的散斑图，而三维dic对实验环境提出了更高的要求，同时操作也更加的繁杂，计算量也明显的增加，如何提高三维dic测量全场变形的实时性，也是当前亟待解决的问题。然而对于一些规则的曲面试件，其表面变形的全场测量除了使用三维dic技术实现外，还可以利用坐标映射结合成熟的二维dic技术实现其表面变形的全场测量，因此亟需一种坐标映射的方法将空间曲面坐标点映射为平面坐标点。

技术实现要素：

为解决以上现有技术存在的问题，本发明提出了一种基于二维数字图像相关法的试件质量检测方法，该方法包括：

s1：向被测器件表面随机制作散斑点，将该散斑点作为位置信息载体；

s2：将试件的横截面中心点作为坐标原点，建立空间坐标系，记录所有散斑点在空间坐标系中的位置坐标；

s3：将试件表面的空间散斑点映射为平面散斑点，得到试件变形前映射平面散斑图；

s4：将图像采集设备固定在指定位置，采用张正友标定法对采集设备进行标定，确定畸变系数，根据畸变系数对图像进行畸变校正；

s5：根据畸变校正后的图像向被测试件施加轴向拉伸力，采集试件形变后的表面散斑图；

s6：将试件形变后的表面散斑图中的各个空间散斑点映射为平面散斑点，得到试件形变后的映射平面散斑图；

s7：采用反向组合高斯牛顿算法对试件形变前的映射平面散斑图和试件形变后的映射平面散斑图进行处理，得到试件表面变形测量结果；

s8：根据试件表面变形测量结果判断被测器件的质量。

优选的，将试件表面的空间散斑点映射为平面散斑点的过程包括：将原空间坐标点投影到xoy平面中；将空间点px0(r,0,0)作为映射后平面坐标的原点o(0,0)，各投影点至点px0的弧长作为映射后平面坐标的x坐标，原空间坐标系下的z坐标作为映射后平面坐标的y坐标；在原坐标系和新坐标系各个点之间一一对应的映射关系下将各个空间散斑点映射到平面坐标系中，得到平面散斑点，根据平面散斑点构建试件变形前映射平面散斑图。

优选的，将图像采集设备固定的位置为该设备采集的图像与空间坐标系的yoz平面平行，且试件表面散斑点成像于像素坐标系中。

优选的，采用张正友标定法对采集设备进行定标的过程包括：

步骤1：将图像采集设备固定在合适的位置；

步骤2：采用图像采集设备采集15张不同旋转角度、不同倾斜角的棋盘定标模板；

步骤3：对获取的15张棋盘定标模板进行定标；

步骤4：根据定标结果获取图像采集设备的内参数；

步骤5：根据图像采集设备的内参数求出畸变系数。

进一步的，将图像采集设备固定的位置为该设备采集的图像与空间坐标系的yoz平面平行，并使被测试件表面散斑点成像于像素坐标系中。

进一步的，获取的图像采集设备内参数包括：内参矩阵、外参矩阵以及畸变系数。

优选的，根据畸变系数对图像进行畸变校正的过程包括：获取理想状态下无畸变的图像坐标和畸变影响下的真实图像坐标；根据获取的坐标建立畸变校正模型；根据畸变系数计算正向畸变系数和切向畸变系数；根据正向畸变系数和切向畸变系数对畸变校正模型进行改进；采用改进的畸变校正模型对图像进行畸变校正，得到理想状态下无畸变发生时的图像坐标。

进一步的，计算正向畸变系数和切向畸变系数的公式分别为：

径向畸变系数公式：

切向畸变系数的公式为：

进一步的，改进的畸变校正模型为：

优选的，采用反向组合高斯牛顿算法对试件形变前的映射平面散斑图和试件形变后的映射平面散斑图进行处理的过程包括：

步骤1：对试件形变前的映射平面散斑图和试件形变后的映射平面散斑图进行常量预计算；预算的常量包括参考子区的图像梯度、jacobian矩阵以及hessian矩阵；

步骤2：对预计算后的常量进行迭代运算，求出各个常量的形函数；

步骤3：根据各个常量的形函数计算其增量δp；

步骤4：根据增量δp判断是否退出迭代，即当||δp||<ε时，终止迭代，否则执行步骤5；

步骤5：对形函数进行更新并返回给步骤3。

进一步的，计算其增量δp的公式为：

本发明的优点：

1)本发明只采用了一个图像采集设备就可以完成试件表面变形的全场测量，避免了标定误差所造成的精度下降影响，节约了测量成本；

2)本发明提出的空间坐标点向平面坐标点映射的方法，不仅保证原坐标系和新坐标系之间一一对应的映射关系，而且计算量小、编程易于实现。

附图说明

图1为本发明的整体实施流程图；

图2为本发明的试件映射结果示意图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

一种基于二维数字图像相关法的试件质量检测方法，如图1所示，该方法包括：

s1：向被测器件表面随机制作散斑点，将该散斑点作为位置信息载体；

s2：将试件的横截面中心点作为坐标原点，建立空间坐标系，记录所有散斑点在空间坐标系中的位置坐标；

s3：将试件表面的空间散斑点映射为平面散斑点，得到试件变形前映射平面散斑图；

s4：将图像采集设备固定在指定位置，采用张正友标定法对采集设备进行标定，确定畸变系数，根据畸变系数对图像进行畸变校正；

s5：根据畸变校正后的图像向被测试件施加轴向拉伸力，采集试件形变后的表面散斑图；

s6：将试件形变后的表面散斑图中的各个空间散斑点映射为平面散斑点，得到试件形变后的映射平面散斑图；

s7：采用反向组合高斯牛顿算法对试件形变前的映射平面散斑图和试件形变后的映射平面散斑图进行处理，得到试件表面变形测量结果；

s8：根据试件表面变形测量结果判断被测器件的质量。

优选的，被测试的试件为圆柱体试件。

如图2所示，将圆柱体试件表面的空间散斑点映射为平面散斑点的过程包括：将空间散斑点的坐标映射到xoy平面坐标系中，将空间点px0(r,0,0)作为映射后平面坐标的原点o(0,0)，各投影点至点px0的弧长作为映射后平面坐标的x坐标，原空间坐标系下的z坐标作为映射后平面坐标的y坐标；根据各个空间散斑点映射到平面坐标系的平面散斑点构建试件变形前映射平面散斑图。将空间点px0(r,0,0)作为映射后平面坐标的原点o(0,0)，各投影点至点px0的弧长作为映射后平面坐标的x坐标的公式为：

y＝z0

其中，r表示圆柱体试件的直径。

采用张友正标定法对采集设备进行定标的过程包括：

步骤1：将图像采集设备固定在合适的位置；

步骤2：采用图像采集设备采集15张不同旋转角度、不同倾斜角的棋盘定标模板；

步骤3：对获取的15张棋盘定标模板进行定标；

步骤4：根据定标结果获取图像采集设备的内参数；

步骤5：根据图像采集设备的内参数求出畸变系数。

优选的，本发明中的图像采集设备为ccd相机。

将ccd相机固定的位置为该设备采集的图像与空间坐标系的yoz平面平行，且试件表面散斑点成像于像素坐标系中。

获取的图像采集设备内参数包括：内参矩阵、外参矩阵以及畸变系数。

根据图像采集设备的内参数求出畸变系数的过程中通过matlab标定直接获得畸变系数。获得畸变系数后对图像进行校正，其具体的过程包括：获取理想状态下无畸变的图像坐标和畸变影响下的真实图像坐标；根据获取的坐标建立畸变校正模型；根据畸变系数计算正向畸变系数和切向畸变系数；根据正向畸变系数和切向畸变系数对畸变校正模型进行改进；采用改进的畸变校正模型对图像进行畸变校正，得到理想状态下无畸变发生时的图像坐标。

建立的畸变校正模型为：

p^u＝p^d+δ(p^d+ki)

其中，p^u＝(x^u,y^u)和p^d＝(x^d,y^d)分别是理想状态下无畸变产生时的像素点和考虑镜头畸变影响下的真实像素点，ki为镜头畸变系数(包括径向畸变系数和切向畸变系数)。

径向畸变的数学表达为：

δ^r(p^d,ki)＝k1r³er+k2r⁵er+k3r⁷er+…

其中，r表示畸变影响下的真实图像到图像平面中心位置的径向距离，x^d表示畸变影响下的真实图像的横坐标，y^d表示畸变影响下的真实图像的纵坐标，cx表示图像平面中心位置的横坐标，cy表示图像平面中心位置的横坐标，ki表示径向畸变系数，表示径向畸变在水平方向上的分量，表示径向畸变在竖直方向上的分量，表示切向畸变在水平方向上的分量，表示切向畸变在竖直方向上的分量，t表示切向畸变(tangentialdistortion)，pi表示切向畸变系数，er为径向单位向量。在径向畸变的数学表达中忽略高于7阶的项，将修改后的径向畸变数学表达式转换为矩阵形式可得：

切向畸变的表达式可以写成：

其中，p1和p2均为切向畸变系数。和中的r和t分别代表径向畸变(radialdistortion)、切向畸变(tangentialdistortion)。

将径向和切向畸变模型代入畸变校正模型中进行改进，得到理想状态下无畸变发生时的图像坐标为：

其中，x^u表示理想状态下无畸变产生时像素点的横坐标，y^u表示理想状态下无畸变产生时像素点的纵坐标，x^d表示畸变影响下像素点的真实横坐标，y^d表示畸变影响下像素点的真实横坐标，ki表示径向畸变系数，r表示表示像素点(x^d,y^d)到图像平面中心位置的径向距离，pi表示切向畸变系数。

采用反向组合高斯牛顿算法对试件形变前的映射平面散斑图和试件形变后的映射平面散斑图进行处理的过程包括：

步骤1：对试件形变前的映射平面散斑图和试件形变后的映射平面散斑图进行常量预计算；预算的常量包括参考子区的图像梯度、jacobian矩阵以及hessian矩阵；对各个常量预计算的公式包括：参考子区的图像梯度公式、jacobian矩阵计算公式以及hessian矩阵计算公式。

参考子区的图像梯度公式为：

其中，表示参考子区的图像梯度，x表示变形前的像素点坐标，ξ表示发生的变形量，表示变形后像素点灰度在水平方向上的偏导数，表示变形后像素点灰度在竖直方向上的偏导数。

jacobian矩阵计算公式为：

其中，δx为水平方向上的变形量，δy为竖直方向上的变形量

hessian矩阵计算公式为：

步骤2：对预计算后的常量进行迭代运算，求出各个常量的形函数。计算各个常量的形函数中，其相关函数的表达式为：

其中，cznssd(δp)表示表征变形前后图像子区相似性的相关函数，δp表示变形增量，ξ表示在子区内各个像素点相对于中心点的局部坐标，w(ξ；δp)表示形函数，p表示变形向量，和分别表示变形前图像子区和变形后图像子区内所含像素灰度的平均值，在参考子区固定的情况下，和不会随着目标子区的变化而变化，δg表示变形后图像子区的灰度误差，δf表示变形前图像子区的灰度误差。

在上述公式中，f(x)和g(x)分别表示参考图像和目标图像在坐标x＝(x,y,1)^t处的灰度值，ξ＝(δx,δy,1)^t表示在子区内各个像素点相对于中心点的局部坐标。

p＝[u,ux,uy,v,vx,vy]^t

δp＝[δu,δux,δuy,δv,δvx,δvy]^t

其中，u表示水平方向上的整像素位移，ux表示水平整像素位移在x方向上的分量，uy表示水平整像素位移在y方向上的分量，v表示竖直方向上的整像素位移，vx表示竖直整像素位移在x方向上的分量，vy表示竖直整像素位移在y方向上的分量，t表示转置，δu表示水平方向上的位移增量，δux表示水平方向上的位移增量在x方向上的分量，δuy表示水平方向上的位移增量在y方向上的分量，δv表示竖直方向上的位移增量，δvx表示竖直方向上的位移增量在x方向上的分量，δvy竖直方向上的位移增量在y方向上的分量。

采用高斯牛顿算法对相关函数进行优化；首先对相关函数做一阶泰勒展开式，得到：

其中，表示参考子区的图像梯度，可用barron梯度算子进行求解；表示形函数的jacobian矩阵。

步骤3：根据各个常量的形函数计算其增量δp；计算形函数增量δp为对一阶泰勒展开式进行求导，令求导后的公式等于0，即最终得到形函数增量δp。具体的计算公式为：

其中，表示hessian矩阵，h^-1表示hessian矩阵的逆矩阵，表示参考子区的图像梯度，表示形函数的雅可比矩阵，ξ表示子区内各个像素点相对于中心点的局部坐标(也就是说要遍历整个子区的所有像素点)，f(x+ξ)表示f(x)变形后的灰度，表示变形前图像子区内所有像素点灰度的平均值，δf表示变形前图像子区的灰度误差，δg表示变形后图像子区的灰度误差，w(ξ；p)为用于表征参考子区变形的形函数，p表示变形向量，表示变形后图像子区内所有像素点灰度的平均值。

步骤4：根据增量δp判断是否退出迭代，即当||δp||<ε时，终止迭代，否则执行步骤5；

步骤5：对形函数进行更新并返回给步骤3。对形函数进行更新的过程包括：计算出δp后，再计算增量的形函数w(ξ；p)，根据w(ξ；p)确定下一步迭代所需的形函数，完成形函数的更新。

迭代所需的形函数为：

w(ξ；p)'＝w(ξ；p)·w^-1(ξ；δp)

其中，w(ξ；p)表示用于表征参考子区变形的形函数，w^-1(ξ；δp)表示形函数增量。

在采用本发明的方法时，需要说明的是在被测圆柱试件表面人工喷射黑白光漆制作随机散斑图，散斑图的制作也影响着搜索精度和搜索速度，在喷射光漆时，不能直接将光漆喷嘴直对试件进行喷射，否则获得的散斑点会过于集中，无法实现全场测量，而应将喷嘴和试件保持一定的角度进行喷射；同时，为了避免由视角误差而引起的映射误差，一定要保证相机成像平面平行于yoz平面。

以上所举实施例，对本发明的目的、技术方案和优点进行了进一步的详细说明，所应理解的是，以上所举实施例仅为本发明的优选实施方式而已，并不用以限制本发明，凡在本发明的精神和原则之内对本发明所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。