一种基于数字散斑的高速三维应变测量方法

摘要

本发明公开了一种基于数字散斑的高速三维应变测量方法。采用高速图像采集系统进行高速形变物体的采集，利用图像的散斑颗粒均值或者灰度梯度评价散斑图像质量，确保测量过程中，图像具有高相关性，高质量，有效提高测量的准确性；在此基础上，根据快速获取的变形物体图像序列对变形过程中的被测物体进行实时连续的三维重建和应变测量，从而最终实现对高速形变物体的三维应变测量。

说明书

技术领域

本发明属于应变测量技术领域，更具体地，涉及一种基于数字散斑的 高速三维应变测量方法。

背景技术

随着材料科学和测量技术的飞速发展，高速应变测量技术在许多领域 有着越来越广泛的应用前景。在实验力学领域，可利用该技术对板料进行 高速拉伸测试，进而得到板料在高速拉伸状态下的应变信息及力学参数， 这对于材料的疲劳实验有着重要意义；在车辆工程领域，可利用该技术对 汽车碰撞实验过程中，车体的高速变形信息进行测量和监控，帮助相关技 术人员改善结构设计并提高车辆的安全性能；在航空航天领域，可利用该 技术测量风洞实验中飞机外壳、机翼等的变形数据，也可以用来分析航空 发动机叶片在高速运行时的形变和振动信息，帮助技术人员快速找到问题， 改善设计，缩短研发周期并提高稳定性。因此，高速应变测量技术对于许 多行业都有着巨大的帮助与促进作用。

现有的高速应变测量技术主要分为接触式和非接触式测量。接触式测 量中，电阻应变计是最为常用的应变测量方法，其测量精度较高，设备成 熟。然而，应变计只能测量单向的应变数据，若要获得密集、多向的应变 数据便需要粘贴非常多的应变计，操作繁琐，而且在高速运动状态下，粘 贴的应变计容易脱落造成测量失败。非接触式测量中，光学测量方法最具 有代表性，主要分为干涉法和非干涉法两大类。干涉法又包括：云纹相干 法、全息相干法等。干涉法一般需要建立一个复杂精密的干涉光路以及干 涉光源，对于振动信息非常敏感，因此，其并不适用于高速三维应变测量。 非干涉法具有代表性的是网格法，是通过采集变形过程中的图像，并分析 图像的特征信息进而获得变形数据，无需复杂光路且对于测量环境要求低。 然而，该方法需要在被测物表面刻画精细标准的网格或者栅线，在实际应 用中具有一定的局限性，此外，对于高速形变的物体，普通的图像传感器 无法获取清晰同步的散斑图像，成为限制高速三维应变测量的壁垒。

发明内容

针对现有技术的以上缺陷或改进需求，本发明提供了一种基于数字散 斑的高速三维应变测量方法，能测量高速应变、振动等信息，并获得高精 度的结果。

为实现上述目的，本发明提供了一种高速三维应变测量方法，其特征 在于，包括如下步骤：

(1)在被测物表面制备随机的散斑图案；

(2)利用第一和第二图像传感器从不同角度拍摄标定物体，提取得到 的两个标定物体图像上的特征点，根据特征点与其在世界坐标系中对应点 的关系，计算第一和第二图像传感器的内参数以及第一和第二图像传感器 间的位置转换关系；

其中，世界坐标系的原点为标定物体上的任一固定点，X轴和Y轴在 原点所在的标定物体表面且相互垂直，Z轴垂直于X轴和Y轴；图像传感 器的内参数包括焦距、主点位置和畸变参数；

(3)利用第一和第二图像传感器从不同角度拍摄被测物，分别计算得 到的两个被测物图像的散斑质量评价参数Q₁和Q₂，将其与阈值Q进行比 较，若Q₁≥Q且Q₂≥Q，则判定散斑质量合格，顺序执行步骤(4)；否则返 回步骤(1)，并调整步骤(1)中的散斑制作工艺，同时调整步骤(2)中 第一和第二图像传感器的位置关系；

(4)利用第一和第二图像传感器从不同角度同步采集被测物图像；

(5)将第一和第二图像传感器采集的被测物图像分别记为左图像和右 图像，获取同一时刻采集的左图像和右图像对应点的匹配关系，根据第一 和第二图像传感器间的位置转换关系，得到被测物形变过程中的三维数据；

(6)根据多个时刻被测物形变过程中的三维数据，比较各个时刻的三 维形貌变化，获得全场三维应变分区图。

优选地，所述步骤(3)中，标记得到的两个被测物图像中的任意一个 为R，按如下方法计算R的散斑质量评价参数Q_R：

(A1)计算R中散斑颗粒的平均尺寸A＝S/N，其中，S为R中散斑 颗粒的总面积，N为R中散斑颗粒的总数；

(A2)计算R的灰度梯度均值其中，W和H分别为 以像素为单位的R的宽度和高度，x_ij代表R中第i行第j列的像素点， 为x_ij的灰度梯度矢量的模，f_u(x_ij)和f_v(x_ij)分别为 x_ij在R的像素点的行所在方向u和列所在方向v的灰度的导数；

(A3)根据R中散斑颗粒的平均尺寸和R的灰度梯度均值，计算R的 散斑质量评价参数Q_R＝α·A+β·B，其中，α为尺寸系数，β为灰度系数。

优选地，所述步骤(5)中，按如下方法获取同一时刻采集的左图像和 右图像对应点的匹配关系：

(B1)选取左图像上的像素点P(x₀,y₀)，以P为中心选取子集，子集 的边长为2M+1，子集中的像素点坐标为(x,y)，计算得到左图像子集的归 一化最小平方距离其中，I(x,y)表示左图像 上坐标为(x,y)的像素点的灰度值，I_m是左图像子集中像素点的灰度均值， $I_m={[1/(2M+1)]}^2\underset{x=-M}{\overset{M}{\Sigma }}\underset{y=-M}{\overset{M}{\Sigma }}I(x,y);$

(B2)令P(x₀,y₀)在右图像上对应的像素点为P'(x₀+ξ₀,y₀+η₀)，以P' 为中心，以与左图像的子集同样大小的窗口在右图像上选取子集，得到右 图像子集的归一化最小平方距离其中，ξ₀和η₀分别为P'点相对于P点的横纵坐标偏移，I'(x',y')表示右图像上像素坐 标为(x',y')的像素点的灰度值，(x',y')为右图像上与左图像上坐标为(x,y) 的像素点对应的像素点的坐标，I'_m是右图像子集中像素点的灰度均值，

${I^{\prime }}_m={[1/(2M+1)]}^2\underset{x=-M}{\overset{M}{\Sigma }}\underset{y=-M}{\overset{M}{\Sigma }}{I}^{\prime }(x^{\prime },y^{\prime }),$

$x^{\prime }=x_0+\mathrm{\Delta x}+{\xi }_0+{\xi }_x\mathrm{\Delta x}+{\xi }_y\mathrm{\Delta y}+\frac{1}{2}{\xi }_{\mathrm{xx}}{\mathrm{\Delta x}}^2+{\xi }_{\mathrm{xy}}\mathrm{\Delta x\Delta y}+\frac{1}{2}{\xi }_{\mathrm{yy}}{\mathrm{\Delta y}}^2,$

$y^{\prime }=y_0+\mathrm{\Delta y}+{\eta }_0+{\eta }_x\mathrm{\Delta x}+{\eta }_y\mathrm{\Delta y}+\frac{1}{2}{\eta }_{\mathrm{xx}}{\mathrm{\Delta x}}^2+{\eta }_{\mathrm{xy}}\mathrm{\Delta x\Delta y}+\frac{1}{2}{\eta }_{\mathrm{yy}}{\mathrm{\Delta y}}^2,$

其中，Δx＝x-x₀，Δy＝y-y₀，ξ_x、ξ_y、η_x和η_y为左图像子集和右图像子 集内像素点的灰度差异的一阶梯度，ξ_xx、ξ_xy、ξ_yy、η_xx、η_xy和η_yy为左图 像子集和右图像子集内像素点的灰度差异的二阶梯度；

(B3)根据左图像子集的归一化最小平方距离和右图像子集的归一化 最小平方距离，计算优化目标函数其中， $\overrightarrow{p}=({\xi }_0,{\eta }_0,{\xi }_x,{\xi }_y,{\eta }_x{,\eta }_y,{\xi }_{\mathrm{xx}}{,\xi }_{\mathrm{xy}}{,\xi }_{\mathrm{yy}}{,\eta }_{\mathrm{xx}}{,\eta }_{\mathrm{xy}}{,\eta }_{\mathrm{yy}}),$计算得到为最小值时ξ₀和 η₀的值，确定P'(x₀+ξ₀,y₀+η₀)的位置；

(B4)改变像素点P的位置，重复执行上述步骤(B1)至(B3)，得 到左图像与右图像对应点的匹配关系。

总体而言，通过本发明所构思的以上技术方案与现有技术相比，具有 以下有益效果：采用高速图像采集系统进行高速形变物体的采集，利用图 像的散斑颗粒均值或者灰度梯度评价散斑图像质量，确保测量过程中，图 像具有高相关性，高质量，有效提高测量的准确性；在此基础上，根据快 速获取的变形物体图像序列对变形过程中的被测物体进行实时连续的三维 重建和应变测量，从而最终实现对高速形变物体的三维应变测量。

附图说明

图1是高速三维应变测量系统示意图；

图2是本发明实施例的基于数字散斑的高速三维应变测量方法流程图；

图3是电化学腐蚀形成的高质量散斑图像；

图4是左图像和右图像的立体匹配示意图。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图 及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体 实施例仅仅用以解释本发明，并不用于限定本发明。此外，下面所描述的 本发明各个实施方式中所涉及到的技术特征只要彼此之间未构成冲突就可 以相互组合。

如图1所示，典型的高速三维应变测量系统包括工作站1、时序同步系 统2、第一高速图像传感器3和第二高速图像传感器4。第一高速图像传感 器3和第二高速图像传感器4连接时序同步系统2，用于从不同角度拍摄被 测物5的图像。通过配置合适的延时使第一高速图像传感器3和第二高速 图像传感器4严格精确同步地采集图像，同时通过调整图像的曝光度、传 感器方位和环境光照设置，获得视场合适且清晰的图像。

如图2所示，本发明实施例的基于数字散斑的高速三维应变测量方法 包括如下步骤：

(1)根据被测物的表面材质和变形速率，在被测物表面制备随机的散 斑图案，并确保散斑图案在外界冲击下不会脱落。

由于物体会进行高速的形变，因此普通的散斑制作工艺，例如喷涂制 备方案将不再适用。一般采用电化学腐蚀方法和激光打标方法进行被测物 表面的散斑制备，电化学腐蚀形成的高质量散斑图像如图3所示。

(2)利用第一和第二图像传感器从不同角度拍摄标定物体，提取得到 的两个标定物体图像上的特征点，根据特征点与其在世界坐标系中对应点 的关系，计算第一和第二图像传感器的内参数以及第一和第二图像传感器 间的位置转换关系。

其中，世界坐标系的原点为标定物体上的任一固定点，X轴和Y轴在 原点所在的标定物体表面且相互垂直，Z轴垂直于X轴和Y轴。图像传感 器的内参数具体包括：焦距、主点位置和畸变参数(包括径向畸变、切向 畸变)。

(3)利用第一和第二图像传感器从不同角度拍摄被测物，分别计算得 到的两个被测物图像的散斑质量评价参数Q₁和Q₂，将其与阈值Q(通常取 25)进行比较，若Q₁≥Q且Q₂≥Q，则判定散斑质量合格，顺序执行步骤(4)； 否则返回步骤(1)，并调整步骤(1)中的散斑制作工艺，同时调整步骤(2) 中第一和第二图像传感器的位置关系。

具体地，标记得到的两个被测物图像中的任意一个为R，按如下方法 计算R的散斑质量评价参数Q_R：

(A1)计算R中散斑颗粒的平均尺寸A＝S/N，其中，S为R中散斑 颗粒的总面积，N为R中散斑颗粒的总数；

(A2)计算R的灰度梯度均值其中，W和H分别为 以像素为单位的R的宽度和高度，x_ij代表R中第i行第j列的像素点， 为x_ij的灰度梯度矢量的模，f_u(x_ij)和f_v(x_ij)分别为 x_ij在R的像素点的行所在方向u和列所在方向v的灰度的导数，利用常规梯 度算子确定；

(A3)根据R中散斑颗粒的平均尺寸和R的灰度梯度均值，计算R的 散斑质量评价参数Q_R＝α·A+β·B，其中，α为尺寸系数，一般取0.2，β 为灰度系数，一般取0.8。

(4)利用第一和第二图像传感器从不同角度同步采集被测物图像。

在被测物开始变形时，根据被测物的材料和测量场景，选择开始标识 信号并发送开始标识信号到两个图像传感器，例如，开始标识信号为声信 号、光信号或电信号等；两个图像传感器收到开始采集标识后，进行高速 同步采集，将数据通过线缆传递至工作站。

(5)将第一和第二图像传感器采集的被测物图像分别记为左图像和右 图像，获取同一时刻采集的左图像和右图像对应点的匹配关系，根据第一 和第二图像传感器间的位置转换关系，通过立体视觉算法得到被测物形变 过程中的三维数据。

具体地，按如下方法获取同一时刻采集的左图像和右图像对应点的匹 配关系，其中，像素点的横坐标为该像素点在其所在图像中的行数，像素 点的纵坐标为该像素点在其所在图像中的列数。

(B1)如图4所示，选取左图像上的像素点P(x₀,y₀)，以P为中心选 取子集，子集的边长为2M+1，子集中的像素点坐标为(x,y)，计算得到左 图像子集的归一化最小平方距离其中， I(x,y)表示左图像上坐标为(x,y)的像素点的灰度值，I_m是左图像子集中像 素点的灰度均值，$I_m={[1/(2M+1)]}^2\underset{x=-M}{\overset{M}{\Sigma }}\underset{y=-M}{\overset{M}{\Sigma }}I(x,y);$

(B2)令P(x₀,y₀)在右图像上对应的像素点为P'(x₀+ξ₀,y₀+η₀)，以P' 为中心，以与左图像的子集同样大小的窗口在右图像上选取子集，得到右 图像子集的归一化最小平方距离其中，ξ₀和η₀分别为P'点相对于P点的横纵坐标偏移，I'(x',y')表示右图像上像素坐 标为(x',y')的像素点的灰度值，(x',y')为右图像上与左图像上坐标为(x,y) 的像素点对应的像素点的坐标，I'_m是右图像子集中像素点的灰度均值，

${I^{\prime }}_m={[1/(2M+1)]}^2\underset{x=-M}{\overset{M}{\Sigma }}\underset{y=-M}{\overset{M}{\Sigma }}{I}^{\prime }(x^{\prime },y^{\prime }),$

$x^{\prime }=x_0+\mathrm{\Delta x}+{\xi }_0+{\xi }_x\mathrm{\Delta x}+{\xi }_y\mathrm{\Delta y}+\frac{1}{2}{\xi }_{\mathrm{xx}}{\mathrm{\Delta x}}^2+{\xi }_{\mathrm{xy}}\mathrm{\Delta x\Delta y}+\frac{1}{2}{\xi }_{\mathrm{yy}}{\mathrm{\Delta y}}^2,$

$y^{\prime }=y_0+\mathrm{\Delta y}+{\eta }_0+{\eta }_x\mathrm{\Delta x}+{\eta }_y\mathrm{\Delta y}+\frac{1}{2}{\eta }_{\mathrm{xx}}{\mathrm{\Delta x}}^2+{\eta }_{\mathrm{xy}}\mathrm{\Delta x\Delta y}+\frac{1}{2}{\eta }_{\mathrm{yy}}{\mathrm{\Delta y}}^2,$

其中，Δx＝x-x₀，Δy＝y-y₀，ξ_x、ξ_y、η_x和η_y为左图像子集和右图像子 集内像素点的灰度差异的一阶梯度，ξ_xx、ξ_xy、ξ_yy、η_xx、η_xy和η_yy为左图 像子集和右图像子集内像素点的灰度差异的二阶梯度；

(B3)根据左图像子集的归一化最小平方距离和右图像子集的归一化 最小平方距离，计算优化目标函数其中， $\overrightarrow{p}=({\xi }_0,{\eta }_0,{\xi }_x,{\xi }_y,{\eta }_x{,\eta }_y,{\xi }_{\mathrm{xx}}{,\xi }_{\mathrm{xy}}{,\xi }_{\mathrm{yy}}{,\eta }_{\mathrm{xx}}{,\eta }_{\mathrm{xy}}{,\eta }_{\mathrm{yy}}),$通过牛顿-拉夫逊 (Newton-Raphson,NR)算法迭代找到的最小值，根据为最小值 时ξ₀和η₀的取值，确定P'(x₀+ξ₀,y₀+η₀)的位置，此时P'(x₀+ξ₀,y₀+η₀)精 确匹配P(x₀,y₀)；

(B4)改变像素点P的位置，重复执行上述步骤(B1)至(B3)，得 到左图像与右图像对应点的匹配关系。

(6)根据多个时刻被测物形变过程中的三维数据，比较各个时刻的三 维形貌变化，获得全场三维应变分区图。

本领域的技术人员容易理解，以上所述仅为本发明的较佳实施例而已， 并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等 同替换和改进等，均应包含在本发明的保护范围之内。