第2章 三维数散斑相关法

第2章三维数字散斑相关法

三维数字图像相关方法（简称是基于双目立体视觉原理和数字图像相

关方法，测量物体表面三维形貌以及三维变形的方法。本章将讨论的

原理及方法。本章首先介绍数字散斑相关法，然后再介绍双目立体视

觉技术。

三维数字散斑相关方法（3D-DSCM）是一种光学测|量方法，通过采集「但注融:基于机翼的日强

目标变形前后的四幅散斑图像，利用双目立体视觉技术进行空间点

的重构、二维数字散斑相关方法（2D-DSCM）进行变形前后的空间点

的对应，在此基础上完成三维坐标及三维变形的测量。3D-DSCM克服

了

2D-DSCM只能测量平面物体二维形变的局限，可以获得任意被测表面

的空间位移及形变，而且具有实时性、对测量环境要求低、试样准备

简单、适用范围广等优点。

2.1二维数字散斑相关法

二维数字散斑相关方法（2D-DSCM）又称为数字图像相关方法（

DIC）,是基于物体表面散斑图像的灰度特征来进行测量的，根据灰度

特征的相关性完成被测物体位移和变形信息测量。下面是相关搜索的

原理，如图2T所示。

数字散斑相关法是一种对试件（受载荷作用下）发生形变前后的散斑

场进行相关运算并以此来获得位移全场信息的测量方法。数字散斑相

关法起源于机器视觉的发展，它具有机器视觉的优点一一非接触式、

全场在线测量等。数字散斑相关方法是由计算机技术、图像处理技术

以及光学技术结合而成的。相比于前文提及的传统光测法，它的光路

相对简单，对测量环境要求低，故其应用面更加广泛。随着数字化技

术迅速发展，其在生物力学、微观结构、材料力学等诸多领域都得到

了相对广泛的应用，同时也促进了其他学科的发展。

在基于数字图像相关法的测量实验中，先采集试件变形前后的散

斑图像，分别将这两幅图像表示为图像A和图（像,B。）如图2.1所

示，在参考图像（2A（+即1变）X形前2的+图1像）中随机选择一

种子点，以户点为中心选择一个（）像素大小的样本子区。然后在图

像（即变形后的图像）中通过搜索算法寻找目标子区。该子区以与

样本子区的互相关系数符合要求的点为

00

中心。从而进一步确定计算点P在x和y方向的位移量u和V。之

所以为P点为中心选择一个样本子区作为搜索点是由于样本子区比

单独的计算像素点包含更多的灰度值信息，所以更易于识别。

Y

（样本子区变形位移）

而在实际计算中，常将样本子区以一定的宽和高划分成若干个虚拟网格，再通过

计算每个虚拟网格节点的位移来近似得到试件表面在变形中的全场位移信息。

数字散斑相关的图像采集系统

图像采集的突出优点之一是所使用的实验设备装置比较简单且

对样本的图像获取比较容易实现。图像采集系统主要由光源（自然光、

白色光源等）、CCD摄像机、高速图像采集卡和计算机组成，如图2.1

所示。在图像采集过程中，首先使CCD摄像机保持不动，然后给试件

加一定的载荷并使之产生形变，从而使其与CCD摄像机产生相对位

移。由于CCD摄像机在图像采集过程中始终保持不动，所以计算机通

过高速图像采集卡所获得的相对位移便是试件在该载荷作用下的实

际位移。为防止计算机在图像采集过程中获得较大的测量误差，首先

CCD摄像机一定要安装在固定的三脚架上不移动，并且从采集第一幅

图像到最后一幅，这个过程中，三脚架的位置始终固定不动。CCD摄

像机与三脚架这两者之间的位置也不能发生相对移动，即使是CCD摄

像机本身产生微小抖动也会使所测量的结果受到很大的影响。其次

CCD摄像机的光轴要与试件表面成近似90°,使试件表面与摄像机镜

头垂直，以确保获得精准的位移值。

以上几个条件更加限制了二维数字散斑相关方法的应用范围。在工程

应用中，物体的变形情况往往超出以上要求，为了更好地满足实际应

用中对空间形变的测量需求，将2D-DSCM拓展到3D-DSCM是非常必要

的，因此在2D-DSCM基础上引入双目立体视觉技术以实现离面位移的

测量。

2.1.3变形数学模型

数字散斑相关方法在对材料的力学分析，归根到底是在求取散斑图像

上像素点的位移。如上一节所述，相关搜索是以窗口为单位进行搜索,

那么就得有合适的变量来表征物体变形前后散斑图上所选定窗口位

移和变形量。由于在实际工业中，大多变形都是非线性的变形，非线

性变形不仅直接导致窗口位置发生位移，而且窗口的形状也会产生变

化。在参考图片中选取变形前的窗口中心点，设为P(x,y),变形后

的目标图片对应窗口的中心点P*(x*,y*),则变形前后相关点的形变

关系为：

“*

II-X-X

<*

iv-y-y“

其中U,v分别表征变形前后的所选点P的x,y方向形变分量。

对于大多数试件加载时候，试件形变包含两部分：刚性位移和非线性

位移，则在此处引入位移梯度来表征非线性为部分。如图2.4所示：

△x

与

图2.4试件形变数学模型示意图

Fig2.4Thematliematicalmodelschematicofspecimendeformation

在计算窗口选取中心点PPy,xP和点P邻域中点QQy,xQ,根据Q与P

坐标之间关系Q点又可表示为yy,xxQPP。根据变形前后点之间一一

对应关系，参考泰勒公式可将形变数学模型表示为：

*.dii.du.

x=x+AX+“H--Ax+-Ay

°nrdx

*A5VAdVA

dxay

diidiidvdv

—>—和一，—

dxdydx,方

式中v,u代表选取窗口的中心点的位移,

分别表示Q点y,x方向的位移梯度，y,x表示计算窗口中任意点

QQy,xQ与窗口中心点PPy,xP之间y,x方的距离。根据式（2.3）可

知，准确的量化试件在形变过程具体数据，就得计算出式中六个参数

diidiidvdv

dx办dx力o

2.1.4亚像素相关搜索算法（放在最后）

在基于数字图像相关法的相关测量中，计算点识别的精度尤为重

要，甚至将影响整个测量的精度，所以搜索算法的选择特别重要。

采集得到的散斑图像中记录的是离散的灰度信息，数字图像相关

法所处理的是经过数字化、灰度化的图像（最小单位是1个像素点）。

在进行样本子区的搜索，样本子区的平移只能以像素为单位进行，

因此最后得到的位移量肯定是像素点大小的整数倍。目前，数字图

像相关法常用的整像素搜索方法包括粗-细法，交叉搜索和遗传算

法。而在实践中，种子点的位移量可能不是整数个像素点，并且由

于有限的像素的CCD照相机，在整数像素精度的测量的位移的定位

精度是不够的。为了提高测量精度，需要在整个像素进一步的亚像

素位移的计算求解的结果的基础上。

近年来，在二维数字图像相关领域，许多研究人员视图通过利用

软件处理的方法来解决数字图像中目标的高精度定位问题。例如，

当算法的精度为0.1个像素时，就相当于测量系统的硬件分辨率提

高了十倍。因此，对图像中目标进行高精度的亚像素位移定位就成

为提高光学测量系统测量精度的重要环节之一。这种亚像素定位技

术具有十分重要的理论意义和实践意义，亚像素位移测量算法对提

高测量精度就非常重要，而且已被认为是数字图像相关法中的关键

技术之一。

从数字图像相关法提出至今，从众多参考文献中可以查阅到的亚像

素位移测量算法主要有以下几种：

1）亚像素灰度插值法

2）曲面拟合、插值法

3）坐标轮换法（十字搜索法）

4）Newton-Rapshon法（简称N-R）

5）拟牛顿法

6）灰度梯度法（简称梯度法）

7）频域相关法

8）后验概率法

9）遗传算法

10）神经网络算法等等

这些算法所能达到的位移测量精度在0.005-0.1像素之间不等。在

以上所列的可供选择的算法中，亚像素灰度插值算法因计算量大、

精度较低，现已少见于应用；坐标轮换法和N-R方法的基本假设相

同，两者区别在于具体的计算方法不同，坐标轮换法通过测试来寻

求事相关函数满足机织条件的待定参数，N-R方法出现后，该算法

也少见于应用；拟牛顿方法与N-R方法的却别在于算法实现方法不

同，但并不能提高亚像素位移测量精度；频域相关法的计算效率较

高，但对变形和转动较为敏感，其定位精度认为是最低的；而后验

概率法、遗传算法、神经网络算法等只限于理论研究。因此，在上

述亚像素位移测量算法中，实际上最常见于应用的三种算法是曲面

拟合法、梯度法、N-R法。

文献对这三种亚像素位移测量方法进行了比较，结果显示在计算子

区大小相同的情况下，N-R方法的计算精度最高，但是其计算时间

约为其它两种方法的30倍左右，梯度法的计算效率与曲面拟合的方

法相当，但却有稍高的计算精度。由于本文的主要工作是将二维的

数字图像相关法拓展到三维数字图像相关，并有效、精确的计算出

被测物体内部的三维位移场，考虑到三维位移场计算量大的问题，

在此只介绍下相对快速的梯度法和曲面拟合法。

2.3.1曲面拟合法

曲面拟合法假设整像素位移相关搜索结果及其相邻8点（见图

2.1,拟合窗口设为3X3像素)的相关系数矩阵可以拟合为连续曲

面，然后以该曲面的极值点位置作为变形后图像子区的中心位置。

通常假设连续曲面可用如下二元二次函数表示：

-22

。(七,为)=。0+。|%+a4xiyJ+asyJ

函数C(x,y)在拟合曲面的极值点即为变形后图像子区中心位置

(x0‘，yO')。因此,就可由u=xO-xO,,v=yO-yO'求出变形前图像子

区的位移，其中(xO,yO)为变形前图像子区的中心位置，u,v即为通

过上述方法计算得到的x和y方向的位移。

图2.1相关系数矩阵中的相邻8点

函数C(x,y)在拟合曲面的极值点应满足以下方程组:

ac（x,y）

=q+2ax+ay=0

~dx-35

ac（x,y）

=%+2a5«y+%x=O

于是，由上式就可求出拟合曲面的极值点位置:

flj—4^3%

2a2a}-ata4

即变性后图像子区的中心位置，进而可由u=xO-xO',v=yO-yO'求出

最终的位移。

2.3.2梯度法

基于梯度的亚像素位移算法是由Freeman等最初从光流计算中引入

并应用到图像定位(Imageregistration)中。文献将该方法引入到

数字图像相关方法中。张军等进一步完善了该方法，给出了基于微小

区域统计特性的亚像素位移梯度算法的四种不同模式。

基于梯度的亚像素位移算法的基本思想为：令f(x,y),g(x',y')分

别表示变形前和变形后的图像子区。根据变形及数字图像的基本假

设，当图像子区足够小且物体做微小位移时，则该子区可看成做近似

的刚体运动，此时有：

f(x,y)=g(x',y')(])

此处：x'=x+"+Ax,y'=y+v+Ay

其中：u,v分别为已经求得的当前参考图像子区中心点的整像素位

移，AX,Ay分别为x,y方向整像素位移对应的亚像素位移。将式

(1)对Ax,Ay进行一阶泰勒展开式并舍去高阶小量，可得:

g(x+〃+Ar,y+u+Ay)=g(x+〃,y+u)+Arxg,(x+〃,y+u)+4yxg、(x+〃,y+v)(2.15)

其中：gx,gy为灰度的一阶梯度，其计算方法稍后讨论。

对于真实的微小位移Ax、Ay,应使下面的最小平方距离相关函数取

驻值：

MM

C(Ax,Ay)=ZZ[f(x，y)r(x+〃+Ax,y+'，+△')『(2.16)

x-My—A/

即:

dC(Ax，4y)(

(2.17)

d(Ax)’5(Ay)

”推导可得(2.18)

其中:

ii等“"ds

sZu-g仪

x—M>—.W隐Lx

A=B=,△=(2.19)

Z"Z"(7-蟾ds△y_

X-M—(<y

由公式(2.19)可以看出，利用基于梯度的亚像素位移算法求解亚像

素位移的关键是灰度梯度gx,gy的计算。不同学者提出了许多灰度

梯度算法。下面列举其中比较常用的几种：

(1)Horn算子。Davis和Freeman在文献中提到的灰度梯度是由与

所求点相邻的4个像素点灰度的一阶差分求得。在此之前这是由

Horn和Schunck最先在文献中提出的。文献认为这是一种比较粗略

的方法，会引起较大的误差：

g.=T［g（x+i'y）-g（X,y）+g（x+i,y+i）-g（x,y+i）］

g,=；［g（x，y+1）-g（xy）+g（x+Ly+1）-g（x+l,y）］

（2）边缘检测算子。在数字图像处理中，常常利用小区域模版来近

似计算图像在某点g（x,y）的灰度梯度来进行边缘检测。对gx,gy

各使用一次模版，这种模版可以看作是二维梯度算子。经常使用的模

版有Prewitt算子，Sobel算子以及各向通性算子。其模版分别如图

所示

图2.2Prewitt算子

图2.3Sobel算子

ZZ二二E

A二二二二

二二3EA

图2.4各向通性算子

(3)Barron算子。文献在讨论光流计算时引入了该算子:

1c881

-+—.?(x+l,y)--i?(x+2,y)

[88]

,?=—^(^y-2)-—^(x,y-l)+—^(x.y+l)-—^(%,>-+2)

'v12121212

该算子是在非节点边界约束条件下对感兴趣点g(x,y)邻域5点的离

散灰度值进行三次样条插值的结果。

(4)自然样条算子。该算子是在自然边界约束条件下对感兴趣点

g(x,y)邻域5点的离散灰度值进行三次样条插值可以得到的另一种

梯度算子：

纵观以上几种梯度算子，在文献中作者对比了各算子的优劣，其中

Prewitt算子，Sobel算子以及各向同性算子的亚像素位移计算精度

在同一数量级上，其误差要比Barron算子的计算误差大，且误差波

动也很大。Horn算子误差的波动范围较大，计算结果稳定性也不好。

从整体效果上来看，Barron算子的计算结果是最为精确的，且很稳

定，因此在我们的实验程序中，计算亚像素位移采用梯度算法中的B

arron算子，来确保算法的精确度和稳定性，在以后的三维亚像素位

移测量当中，我们在此基础上拓展出适合三维亚像素计算的算子。

2.3.3N-R法

Newton-Raphson偏微分方法基本理论就是牛顿迭代法，在赋予合适

的初值后，经过不断迭代改善待求的参数，以期函数达到最大或者最

小值。在2.1.2节中提到数字散斑相关方法要求六个参数

dudiidvdv

’&'为'ax'®将Pi与式(2.3)代入式(2.1),可得到S是一个

关于P的函数S(Pi)。在给定合适的初值P0后，N-R

算法通过不断的迭代生成新值，达到收敛条件即找到最相关的点。由

于N-R算法是在牛顿迭代基础上演变而来，那么迭代公式为：

+川=々_玄飞片）*▽仍）

其中iP为雅各比向量，每个参数的一阶导数构成雅各比向量的元素;

iHP为海森矩阵，其元素是由相关公式（2.1）对各个形变参数分量

的二阶偏导数所组成。

dv

dS

_du、

况

dx

dS

V(甲

_du、

比一)

力

dS

d(-T)

dx

dS

心）

力.

在这个迭代过程中不断生成新的位移参数后，根据式（2.3）生成新

的相关计算点的坐标，但由于计算点坐标并不是处于整像素位置，那

么就需要利用双三次多项式插值求出新计算点的灰度值。（1）若在迭

代过程中Pi+1满足|Pi+「Pi|〈时，迭代终止，为允许的误差。（2）

设定最大的迭代次数maxN,若迭代的次数达到maxN则停止计算。

2.2双目立体视觉技术

人类是通过两眼分别同时获取外部场景的二维图像，然后经过大脑的

处理，从而得到外部场景的三维信息。双目立体视觉的基本原理与此

类似，即利用两个摄像头记录下空间同一场景的图像，然后寻找这两

幅二维图像中的对应点，根据已知的两个摄像头的内部参数，计算得

到其相对于空间中某个坐标系的三维坐标。

图2.5双目立体视觉示意图

图2.5是关于双目立体视觉原理的简单示意图，其中符号c表示两个

摄像机的光心。从图中可以看出，空间中的点P分别成像于点P1

（位于摄像机1的像平面上）与点P2（位于摄像机2的像平面上），

双目立体视觉的目标，就是要从点P1和点P2确定点P在预先设定的

世界坐标系中的坐标。从图2.5中可以很直观的看到，如果只有一个

摄像机，那么只能得到空间中的一条直线，至少要两个摄像机，才可

以唯一确定点P的坐标。

双目立体视觉测量中，有两个关键性的步骤：标定和匹配。标定是确

定摄像机内外参数的过程，匹配是寻找测量系统中两个摄像机分别所

记录图像中的对应点。如图2.5所示，如果已经得到了摄像机内外参

数，利用每一对如P1和点P2这样的对应点的图像坐标，就可以确定

空间中的一个点。

2.3三维数字散斑相关测量技术

2.3.1三维数字散斑相关法的测量原理

如图3.2所示，三维数字图像相关法使用两个固定的相机从不同

角度拍摄被测物表面变形前后的数字图像，再通过图像匹配算法计算

得到变形前后数字图像中待测点的图像坐标，再结合事先标定好的两

个相机的参数和相对位置关系，计算得到变形前后待测点的三维空间

坐标，最后以此解算得到全场三维应变。从上述描述中可以看出，三

维数字图像相关法中的关键技术主要包括以下四项：

a.图像采集一一同步获取两个相机被测物表面的数字图像；

b.图像匹配一一在两幅图像中寻找对应点；

c.相机模型参数标定及三维重建一一标定两个相机成像模型参

数，借助这些参数结合数字图像中对应点的图像坐标计算该点的三

维空间坐标；

d.空间位移计算一一由变形前后待测点的三维空间坐标求解空

间位移。

2.3.2摄像机成像模型和标定

［批注［A2］：数字散斑相关法的研究与应用钟健

2.3.3图像三维空间坐标的建立

摄像机成像模型，是用数学语言描述空间中的点成像于摄像机靶面

上，并经过AD转换，形成数字化图像上的像素点的简化过程。根据

是否考虑摄像镜头畸变系数的影响，可将摄像机成像模型分为线性模

型和非线性模型。

一般的摄像机成像模型涉及到如图2.6所示的五个坐标之间的转换

关系：

1.空间世界坐标wwwwzyxO到摄像机坐标cccczyxO的转换，表示如

下：

X

二R+T

ycyw

ZcZw」

其中，R是旋转矩阵,具有正交性，即T=,T是平移向量

R=，4，5，6

「7%「9

2.摄像机坐标系cccczyxO到理想摄像机像平面坐标系uuu

yxO的转换:

4

4X,F

1z

其中,

f00

0/0

001

这里，f是摄像机焦距。

注意，摄像机坐标系cccczyxO中坐标原点是摄像机光心，cz

轴与光轴重合。理想摄像机像平面坐标系uuuyxO中的坐标平面与摄

像机坐标系cccczyxO中的ccyx平面平行，这两个平面间的距离即是

摄像机焦距fo理想摄像机像平面坐标系uuuyxO坐标原点是

CZ轴与uuyx平面的交点，且UX轴与uy轴分别与ex和c轴平行。

2.4亚像素搜索算法

数字图像是以离散的灰度信息存储在计算机中并进行位移计算的且

其灰度单位是像素。因此所计算出的位移量只能是像素的整数倍。为

了提高位移测量的精准度，可以应用更先进的图像采集系统，如提高

分辨率或采用高放大倍数的光学成像系统。但随着分辨率升高相应的

图像采集成本也会增高且分辨率在技术上不可能无限制地提高，而采

用高放大倍数的光学成像系统虽然能提高测量精度但相应的测量面

积却会缩小很多。经过对对亚像素位移测量方法的多年研究，证明了

该方法对提高测量精度非常有效。经过学者们对亚像素位移测量算法

的不断研究与创新，发展出了多种亚像素位移测量方法，如：亚像素

灰度插值法、迭代法、曲面拟合法、梯度法、拟牛顿方法、坐标轮换

法、频域相关法、遗传算法、神经网络算法等，而其中尤以曲面拟合

法、N-R法和梯度法应用最为广泛。

2.4.1曲面拟合法

2.4.2N-R法

2.4.3梯度法

2.5本章小结

3.3摄像机的成像模型

为了适应三维计算机视觉的需要，我们首先采用能精确反映成像物理

过程的摄像机模型一一针孔模型(Pin-holemodel),如图3-1所示。

Ze

图4-1针孔模型

Fig.4-1Pin-holeimagingmodel

该模型为针孔模型。针孔模型是一种最常用的理想透视投影模型，它

忽略了成像光路中各种误差的影响，成像关系为线性的，也就是相当

于物理意义上薄膜透镜成像。

模型中各参考坐标系定义如下：（1）世界坐标系，根据自然环境所选

用的坐标系，坐标值用O,,ZYX表示；（2）摄像机坐标系，以摄像机

的光心cO点为坐标原点，ex轴、cy

轴分别平行于CCD平面的垂直边，摄像机的光轴为cz轴，坐标值用

（）ccc,,zyx表示，（3）图像坐标系，坐标原点在CCD图像平面的中心

O,x轴、y轴分别为平行于CCD平面两条垂直边，坐标值用（yx,）表

示；（4）像素坐标系，坐标原点在CCD图像平面的左上角，u轴、v轴

分别平行于图像坐标系的x轴、y轴，坐标值用

（VU,）来表示。

根据透视投影模型，可得摄像机理想成像数学模型：

XW

匕

5v=K[R

Zw

1

1

(1)

式中s——比例因子；

K——内参数矩阵；

：RT]——外参数矩阵。

内参数矩阵定义为：

a0

K=0%vo

001

(2)

式中,xyaa---分别表示在图像u,v轴的比例因子，且xxa=sf,

yya.=sf；f是主星巨，，xyss

是图像平面单位距离上的像素数；00u,v——主点坐标，是图像坐标

系原点在像素坐标系中的坐标。

外参数中旋转矩阵定义为：

JKK

R=r21r22r23

J31「32r33_

另外，旋转矩阵R还可以表示为绕xc、yc、zc轴的旋转角（即欧拉角，，

Ke3）的形式，如图4-2所示。

图4・2欧拉角表示法

Fig.4-2Euclidangleexpression

其表达式为:

R=RKA

%r22

%

COS。COSK-cossinA:

二cosd?sinKT+sin<i?sincosxrCOSGCOSK-sinGsin夕sin/c-sin^cos。

sin6ysin^-cos6WsincosKTsinCOSAT+costysin^sinA：COSGCOS，

因此，由式(4-4)可以确定欧拉角a、6和K,

分别为：

m-arctan--

k「33，

<O=arcsin(r31)

(〃，、

K-arctan--

、IrH7

外参数中平移向量定义为：

旋转矩阵和平移向量统称为摄像机的外参数。

对公式(4-1)消去比例因子s,即得理想成像模型:

，iX”,+6五+不2力+T、

11=%

+，32工0+乃3Z1V+T二

+尸”匕+r,3zw+T

v=ay---：：