2008数学建模论文

数码相机定位

摘要：本文主要研究的是采用数码相机进行双目定位过程中有关系统标定的问题。即对物体的某个特征点，用两部固定于不同位置的数码相机摄得物体的像，根据该特征点在两个像坐标系里的坐标和两部相机的相对位置来确定特征点的位置。

第一步，论文考虑单个固定相机摄像问题。首先，通过对图像数据处理和分析，文章对靶标和所摄得的像分别在各自的平面内建立独立的坐标系。考虑影响相机拍摄效果的各个因素，诸如相机与靶标的距离、相机的光学平面与靶标的夹角（空间上可分为三个方向的偏角）、相机透镜自身的放大率、焦距、像距、对光线折射的畸变等，再根据单个固定相机的特殊性，本文提出了一种二维图像测定系统标定的方法，给出了两个坐标系图像相应点的非线性函数对应关系。

第二步，利用上述得到的函数关系，计算若干特征点对应像平面上的位置（坐标值）。利用计算机对图像数据进行分析处理，采用网格式平行线和特征交线标定法，选取特殊点，并利用这些点的位置数据评定上述模型的精度和稳定性。分析表明，该模型对单个固定相机拍摄所得的靶标特征点在像平面上的位置预测是可行的。利用该模型对原始靶标进行在像平面的投影变换，其所得到的像基本与单个固定相机所拍摄的像相重合。

最后，考虑两个相机的定位问题。论文在二维图像测定系统标定的基础上，将上述的模型推广到三维立体空间坐标系。为了问题的求解，本文以靶标的中心为原点建立世界坐标系，将两个相机中的像平面坐标都转换为世界坐标系中的坐标。这样，两个相机的位置参数就可以用同一个坐标系（世界坐标系）来度量，从而建立了求解两个相机相对位置的数学模型和方法。

总之，本文根据实际问题情况，建立了在数码相机进行双目定位过程中的系统标定问题。给出了单个相机摄像的二维图像测定系统标定模型和两个相机摄像的空间定位模型。模型简单实用，易于理解，具有很好的测定精度和稳定性。

关键字：系统标定轴偏变换旋转变换平移变换径向畸变换非线性拟合精度评价稳定性分析相对位置

一、问题的重述：

数码相机定位在交通监管（电子警察）等方面有广泛的应用。其中最常用的定位方法就是双目定位，也就是用两部相机来定位。对物体上一个特征点，用两部固定于不同位置的相机摄得物体的像，分别获得该点在两部相机像平面上的坐标。只要知道两部相机精确的相对位置，就可用几何的方法得到该特征点在固定一部相机的坐标系中的坐标，即确定了特征点的位置。于是对双目定位，精确地确定两部相机的相对位置就是关键，这一过程称为系统标定。

本文采用的标定做法是：在一块平板上画若干个点，同时用这两部相机照相，分别得到这些点在它们像平面上的像点，利用这两组像点的几何关系就可以得到这两部相机的相对位置。然而，无论在物平面或像平面上我们都无法直接得到没有几何尺寸的“点”。实际的做法是在物平面上画若干个圆（称为靶标），比如取1个边长为100mm的正方形，分别以四个顶点（对应为A、C、D、E）为圆心，12mm为半径作圆。以AC边上距离A点30mm处的B为圆心，12mm为半径作圆，如图1所示。

图 1 靶标示意图

用一位置固定的数码相机摄得其像，如图2所示。

图2 靶标的像

根据上提示可以解决解决以下的问题：

（1）建立数学模型和算法以确定靶标上圆的圆心在该相机像平面的像坐标, 这里坐标系原点取在该相机的光学中心，x-y平面平行于像平面；

（2）对由图2、图3分别给出的靶标及其像，计算靶标上圆的圆心在像平面上的像坐标, 该相机的像距（即光学中心到像平面的距离）是1577个像素单位(1毫米约为3.78个像素单位)，相机分辨率为1024×768；

（3）设计一种方法检验你们的模型，并对方法的精度和稳定性进行讨论；

（4）建立用此靶标给出两部固定相机相对位置的数学模型和方法。

二、模型的假设及符号说明

2.1符号的说明

X：靶标上特征点的横坐标（所在坐标系）；

Y：靶标上特征点的纵坐标（所在坐标系）；

Z：靶标上特征点的竖坐标（所在坐标系）；

X：靶标上点的像的横坐标（所在坐标系）；

Y：靶标上点的像的纵坐标（所在坐标系）；

Z：靶标上点的像的竖坐标（所在坐标系）；

X：横向平移变换量；

Y：纵向平移变换量；

Z：竖向平移变换量；

M：光学系统的放大倍数；

?：主光轴旋转角（用于二维空间）；

θ：主光轴轴偏角（用于二维空间）；

α：相机透镜横向偏角（用于三维空间）；

β：相机透镜纵向偏角（用于三维空间）；

γ：相机透镜竖向偏角（用于三维空间）；

k：透镜的径向畸变系数；

r：图像中心点到相机透镜中心的距离；

s：像素单位；

.3

=。

78

C：距离单位mm对像素的当量，即mm

s

C/

2.2 基本假设

为了问题的描述和模型的建立，本文对题目做了如下假设：

1.相机透镜是均匀的凸透镜，只考虑透镜的径向畸变因素；

2.透镜所成像即数码相机得到的像，不考虑透镜光学所成像与相机内所存储的数字图

像的差异；

3.靶标的位置是固定的；

4.相机与标靶的相对距离要大于像距，但具有可比性；

5.不考虑在计算机对图像数据进行处理时造成的数据误差。

三、模型的建立与问题的分析

3.1 图形数据的预处理

本工作的主要用途是将题目所给图像存入计算机，转化为计算机可以识别的图像数据，同时对图像数据进行处理得到靶标上特征点对应的像平面的位置表示。

题目所给图像数据在计算机上的存储为是三维的BRG格式。我们为了建立原靶标上的各个点与相机所得像上坐标的对应关系，首先得寻找几个特征点，同时得到这些特征点在图像上的位置表示，这样才能利用这些初始数据建立标靶坐标与所成像坐标的对应关系。

这些特征点对应像平面上的图形中心。本文确定这些像的中心点位置采用的工作主要分为下面几个步骤：

1.加载图形，得到灰度图形数据；

采用imread()函数，得到BRG格式的数据，数据量大。为了降低之后数据处理的不便，将其转换为颜色格式只有一维的灰度图形格式，数据分布在0到255之间。

2.确定图像阈值，标定图像边界，得到标靶特征区域对应像的范围；

由于最初的图形数据中图像的颜色灰度分布是阶梯分布，从0到255渐进变化，尤其是图形边界附近的像的灰度值更是如此。要确定图形的中心，首先就要确定图形的边界，这样才能确定整个图形的形状，才能确定中心点，即标靶特征点在像平面上的位置表示。

由于本文所研究图像并不复杂，图像边界的色度变化速度比较快。所以文章在设定图像阈值时，选取了灰度值的中值，即122。简化了确定图像边界点的难度，但也一定程度上可能会对某些图像的中心发生偏离。

3.求得中心点位置。

确定图像区域之后，将像平面上所有点上对应的像的灰度值二值化。即将所有图形区域内的点的灰度值设为255（纯白色），而将区域外点对应的灰度值设为0（纯黑色）。这样做的目的：一是增加图像对比度，突出目标区域；二是利用MA TLAB图像工具箱中分析图像的图形工具regionprops（）函数（该程序针对白色图形区域）。

为了求得五个特征区域的像的中心，本文对最初的数据进行了分割，利用搜索算法得到分割点坐标，使五个特征图像完全分割为五个独立的子图像。整个过程的实现详见附件一的程序。这样就可以我们所要的位置坐标。

3.2 问题一：坐标转换模型的建立

正如前面所说的，影响相机图像拍摄效果的因素有很多。从三维空间角度考虑，有标靶的位置因素，相机的位置因素，也就是标靶和相机之间的相对距离；从相机的拍摄角度考虑，有相机的主光轴于相机透镜所在平面（与像平面平行）与标靶所处平面和中轴线所成的角度因素，这里包括横向、纵向和竖向的偏离角；从相机本身的光学特性来说，会涉及到相机的放大率（与透镜的折射率相关）、透镜对光线折射的畸变、焦距和

像距等多方面的因素。

问题一考虑单个固定的相机成像，上述众多因素对于同一次拍摄而言并不会影响到所研究的各个特征点在相机中相对位置。考虑我们得到的是一个二维平面图像，我们分别对标靶所在平面和像平面建立独立的坐标系，记为Y X -坐标系和Y X -坐标系，坐标系中心分别为标靶中心和像平面中心，单位长度均为像素s 。当相机固定位置以后，影响标靶特征点在图像坐标系位置的因素就可以进行简化。此时，像平面与标靶平面所成的角度关系，是可以简化，利用二维角度来表示，他们分别为像平面的中心主轴和标靶平面的中心轴的偏角，即轴偏角θ和像平面相对于标靶平面的旋转角?。这样我们把问题做如下处理：

1. 考虑相机光学系统的放大倍数为M ，根据成像关系可得到Y X -坐标系的一个

变换，记为S S Y X -坐标系：

MX X S = MY

Y S =

2. 轴偏变换

由于相机的像平面的水平与标靶的垂直方向不一定垂直，存在一个轴偏角，故我们可以对S S Y X -坐标系进行进一步的轴偏变换，得到的变换结果记为

θ

sin s s

w

Y X

X

-=

θ

cos s w Y Y =

即

θ

sin MY MX X

w

-=

θ

cos MY Y w =

3. 旋转变换

在轴偏变换的基础上，对图像进行平面旋转，得到新的变换结果记为

?

?sin cos w w

D

Y X

X

-=

?

?cos sin w w

D Y X

Y +=

即 ?θ?θsin cos cos )sin (MY MY MX X D --=

)sin cos cos (sin cos ?θ?θ?+-=MY MX

)sin sin cos (cos sin θ?θ??-+=MY MX Y D

得到 )sin(cos θ??++=MY MX X D

)cos(sin θ??++=MY MX Y D

4. 平移变换

由于相机主轴与标靶平面的非垂直型，会造成图像的一个平移量，记为

0)sin(cos X MY MX X t +++=θ??

0)cos(sin Y MY MX Y t +++=θ??

5. 径向畸变换

考虑摄像机透镜的径向畸变因素，将理想图像坐标(,)t t X Y 转换成实际图像坐标

(,)X Y ：

02

2

2

2

1cos sin()1111t X M M X X X Y kr kr

kr kr ??θ+==+

+

++++

22221sin cos()

1111t Y M M Y Y X Y kr

kr

kr

kr

??θ+=

=

+

+

++++

其中r 整个图像的变换过程可以简单的用示意图3表示为

图3 单个固定相机成像图

由于相机位置已定，其中的参数值不会再随特征点位置变换而变换，于是我们可以做下面的简记

?cos 1M a =；)sin(2θ?+=M a ； ?sin 1M b =； )cos(2θ?+=M b

其中角度可由11

tan b a ?=

，sin θ=

从而得到最终问题一的数学模型为

()()2

2

1

2

1k X Y X a X a Y X ++=++ （一）

()()22

1

2

1k X

Y

Y a X a Y Y ++=++

可见，靶标上特征点到像平面的转换模型是一个三次的非线性方程组形式。 下面我们对模型进行求解

为了得到模型一中未知参数的值，首先要利用初始数据进行拟合求解。在靶标所在

在数据的预处理中得到的对应特征点在像平面内的像坐标如下表2所示（单位：像

在对模型一非线性曲线进行拟合计算时，,,,X Y X Y 均为已知，此时要确定模型中未

知参数，则相当于求解一个三维的线性方程组，即一个关于未知参数的线性回归模型，其中未知参数为1a ，2a ，1b ， 2b ，0X ，0Y ，k .

(

)2

2

012X a X a Y k X

Y

X X ++-+=

(

)2

2

012Y b X b Y k X Y

Y Y ++-+=

求解该模型可以采用最小二乘法来拟合，记 E=[ 0a 1a 2a k 0b 1b 2b ]T

12A=.....n A A A ????????????

,B=??

??

??

??????n B B

B ....21,()()2

2

22

100000

1

i X Y X

Y X A X

Y

Y

X

Y ??

+??

=??

+??

?

?

,i X B Y ??=???? 利用最小二乘法得到的方程的解为 E=()B A A A T T 1

-

代入数据，利用附件中程序，可得到解为

E=[-48.637 0.77274 0.11458 -9.3202e-007 11.426 0.037826 -0.76502]T 于是，最后问题一的数学模型为

()()

2222

-48.6370.772740.11458+9.3202e-00711.4260.037826-0.76502+9.3202e-007X Y X Y X X

X Y X Y Y Y ?+++=?

??++=?

由于实际图像的分辨率为1024?768，并且图像的坐标系中心为图像中心，所以模型解空间[][]512:512384:384-?-，同时因为数据在计算机存储点均为整数，所以最终求解时可取近似整数解。

3.3 问题二：利用模型求特征点坐标

在上文中，我们通过最小二乘拟合的方法得到了模型一的具体表达式。但是当用模型一求解靶标某个特征点在对应像平面的位置时，需要解决的方程是一个非线性的方程组，我们可以借助于MATLAB 软件中的函数命令，对不同的初坐标进行非线性方程组求解。

首先在原标靶平面内所求特征点的坐标依旧是如表1所示。 求A 点为例，其求解具体程序为

S=solve('(1-9.3202e-007*(x^2+y^2))*x=-48.637-0.77274*48*3.78+0.11458*50*3.78','(1-9.3202e-007*(x^2+y^2))*y=11.426-0.037826*48*3.78-0.76502*50*3.78');

选取满足题目要求的解，并取近似整数，得到该点在图像平面上的点的坐标为

()-176,147A

同样的，最终通过模型求得靶标上所有特征点在相机像平面的像坐标的估计值如下表3

5%以下，这在工程应用和实际生活中是可容许的。同时偏差的方差为4.4307，也是一个比较合适的数值。可见该模型对预测数据的准确性。

3.4 问题三：模型的精度和稳定性分析

在上一节中，我们通过利用模型一得到的四个特征点在像平面上的坐标已经对模型的预测精度有一定的认识。本节中我们将引入网格式平行线和特征交线标定法对模型精度进行进一步的分析，同时结合差分，分析模型一的稳定性。

影响模型预测精度的因素主要来源于以下三个方面：

一是在第一节中计算机对图像数据预处理的过程，采用图像学知识确定各特征区域所成图像的边界，从而确定图像的中心。由于确定图像边界阈值的选取标准不同，就会会对图像中心位置的确定产生偏差，也就是对特征点在像平面的坐标产生偏差。

二是在模型参数确定过程中取点数量和质量对拟合效果的影响。

三是由于计算机在存储图像数据时以整点的像素为单位，不存在小数，致使在坐标计算是需要对计算结果进行取整运算，从而为模型的解带来一定偏差。

我们这里模型的精度主要分析由第二点造成的精度问题。我们要提高这方面的精度，主要是增大有效数据量。为了然增加的数据具有代表性，采用网格式平行线和特征交线标定法，它主要是利用平行线和特殊特征曲线来取点验证模型的估计精度。

网格式平行线法是充分利用在原标靶平面内各条水平或垂直的线上的特征点的坐标，通过在各条水平或垂直的线上提取额外的特征点，一起来参与模型一的拟合修正，来验证模型的精度。特殊特征曲线取点是利用特征曲线交点，比如标靶平面内各个原数据点连线所成直线的交点。通过测定这些具有鲜明特点的数据来验证模型拟合的精度。

通过取样，可以扩充特征点数据为

特征点在图像中的坐标对应关系。

反映模型拟合精度的另外一个更形象的方法就是利用模型对标靶图形进行图像重构，也就是通过模型计算原标靶特征点周围圆形区域在像平面上所成的图像。该图像与原图像进行比较，看图像的重合程度。具体情况见图4和图5。其中实体的图像是原图像，虚体的图像为通过模型重构的图像。

从图中我们可以看出来，对于原图像和由模型拟合重构的图像基本上偏差并不大。尤其是平面左上方的图像，两个图形基本是重合的。

图4 利用模型求解所得图像（虚像）和原像的重合效果图（整体）

图5 利用模型求解所得图像（虚像）和原像的重合效果图（部分）

图4和图5不仅说明了模型的估计精度，还可以说明模型的稳定性。因为当数据在变动时，整个图像并没有发生大的偏差。为了进一步说明模型的稳定性，我们从理论上对模型做差分，利用变量的变化规律来解释。

我们综合考虑模型一的求解拟合过程。在模型一中，显然，当(,)(0,0)X Y ??→时，我们是可以认为2

2

2

r X Y

=+不变

02

2

2

02

2

2cos sin()111sin cos()

111X M M X X Y kr kr kr Y M M Y X Y kr

kr

kr

ααθα

αθ+=++++++=

+

+

+++

两边取微分，可得：

2

2

2

2

cos sin()11sin cos()

11M M d X dX dY

kr kr M M d Y dX dY

kr

kr

ααθα

αθ+=++++=

+

++

即

2222cos sin()11sin cos()11M M d X dX kr kr M M dY d Y

kr

kr ααθααθ+??

????????++??

=?

?????+?

???????????++??

此线性方程组的系数矩阵的行列式为：

22

2

2

2

2

2

cos sin()11cos(2)sin cos()

(1)

11M M M

kr kr M M kr kr

kr

ααθαθα

αθ+++=

+++++

由于α和θ的角度都非常小，故cos()0αθ+≠，则此线性方程组的系数矩阵的行列式不等于0。因此，当0d X →，0d Y →时，有0dX →，0dY →。从而可以说此模型是稳定的。

问题四：双相机定位模型

在实际应用中，我们用两个相机对实物进行定位，该过程主要是要确定相机内部的几何及光学特性和相机坐标系相对于世界坐标系（以标靶为中心建立的坐标系）的三维位置和方向。我们要解决的主要问题是确定两个相机间的相对位置。为此，我们必须建立两个相机的平面坐标系相对于世界坐标系的对应关系，从而在统一的一个坐标系中得到两个相机的位置变量，建立双相机的定位模型。

和二维图像的标定模型相类似，我们也是考虑在理想的透镜成像基础上加入透镜的径向畸变因素。其示意图为图6所示：

图6 双相机成像示意图

图中(,,)W W W X Y Z 是物体表面上的某一点P 在三维世界坐标系W W W O X Y Z 中的三维坐标。l l l l o x y z 和r r r r o x y z 分别为左右相机的坐标系，坐标中心l O 和r O 分别为左右相机的光学中心，l z 和r z 分别为左右相机的光轴重合。由于相机所得像总是二维平面的，所以标记il l l O X Y 和ir r r O X Y 分别为左右相机成像平面上的图像坐标系。图像中心il O 为光轴l z 与左相机图像平面的交点，图像中ir O 为光轴r z 与右相机成像平面的交点，l X 和l Y 轴分别与右摄像机坐标系的r x 和r y 轴平行。cl l l o u v 和cr r r o u v 分别为左右相机的成像在计算机中

（帧存）的图像坐标系，采集到的图像在计算机帧存中的坐标单位以像素表示。(),u u X Y 是P 点在理想的透镜下的图像坐标。(),d d X Y 是P 点因透镜径向变形而偏离了理想图像位置的位置。有效焦距f 是光学中心到图像平面的距离。

通过上面的示意图和在二维图像标定模型的处理方法，我们可以进行下面的计算处理。

首先以左相机的成像过程为例，针对我们具体的问题，以标靶中心为原点，建立物体、相机和图像的三维世界坐标，某特征点坐标记为(,,)X Y Z 。变换到计算机图像平面坐标(,)u v 的步骤如下：

（1）从物体三维世界坐标(,,)X Y Z 到相机三维坐标(,,)l l l x y z 的三维空间刚体位置变换：

l l l l l x X y R Y T Z z ????????

=?+????????????

其中l R 为33?的正交旋转变换矩阵，l T 为31?的平衡矢量。l R 确定了相机相对于世界坐标系的方向，它的矩阵元素可由三个欧拉角γ，β和α来表示，分别代表了相机相对于坐标原点的竖向、纵向和横向的偏折角度。

11121321

222331

32

33l r r r R r r r r r r ??

??=??????

其矩阵元素表示为：

11cos cos r αβ

=，21sin cos r αβ=-，31sin r β=，12sin cos cos sin sin r αγαβγ=+，

22cos cos sin sin sin r αγαβγ

=-， 32cos sin r βγ=-，13sin sin cos sin cos r αγαβγ=-，

23cos sin sin sin cos r αγαβγ=+ ，33cos cos r βγ= 平移矢量T

l x

y

z T T T T ??=??

，则表明了摄像机相对于世界坐标系的位置关系。该映射对

应了二维图像标定中的偏折变换、旋转变换和平移变换。

(2) 相机的图像始终是平面的，则相机三维坐标(,,)l l l x y z 在理想的透镜成像模型下经投影变换到图像坐标(,)u u X Y 上为

l u l l u l

x X f z x Y f

z ==

(3) 考虑相机透镜的径向畸变因素，将理想图像坐标(,)u u X Y 转换成实际图像坐标

(,)l l X Y ：

2

1

21

(1)(1)

l u l u X X kr Y Y kr --=+=+

其中r =

i O 到实际图像坐标(,)l l X Y 的距离；k 表示透镜径

向畸变系数。

综合上述变换，我们可以建立世界坐标系中标靶特征点坐标与其在左相机中所成像的位置坐标这二者之间的关系：

2

11121331323321

22232313233(1)(1)x l z y l

z r X r Y r Z T kr X f r X r Y r Z T r X r Y r Z T kr Y f r X r Y r Z T +++?+=?+++?

?

+++?+=?+++?

（二） 参考模型一的建立过程，我们同样可以建立单个相机在世界坐标系下的关于物体点与其像对应的转换模型，也就是上面的对应关系，记为模型二。

模型二同样是一个非线性的方程组，其中模型的参数可以通过数据拟合来得到。当参数为已知后，我们就可以利用其中参数确定相机在世界坐标系中的位置T

l x

y

z T T T T ??=??

。

同理，对于右相机也有相同的结果。同样可以通过建立世界坐标系中标靶特征点坐标与右相机所成像位置关系，确定参数，从而得到右相机在世界坐标系中的位置。

对右相机而言，记从物体三维世界坐标(,,)X Y Z 到相机三维坐标(,,)r r r x y z 的三维空间刚体位置变换

r r r r r x X y R Y T Z z ????????

=?+????????????

其中同样的，r R 为33?的正交旋转变换矩阵，r T 为31?的平衡矢量。r R 确定了相机相对于世界坐标系的方向，它的矩阵元素可由三个欧拉角'γ，'β和'α来表示，分别代表了相机相对于坐标原点的竖向、纵向和横向的偏折角度。

111213212223313233'

'''

'''

'

'l r r r R r r r r r r ??

??=??????

其中'ij r 形式与左相机类似。而'''T

r x y z T T T T ?

?=??

同样的经过投影变换和径向畸变换，最终得到右相机中像与其对应物点的函数关系

2

11121331323321

22232313233''''(1)''''''''

(1)''''x r z y r

z r X r Y r Z T kr X f r X r Y r Z T r X r Y r Z T kr Y f r X r Y r Z T +++?+=?+++?

?

+++?+=?+++?

于是可以通过数据拟合得到其中各个参数的值，也就得到了右相机相对于世界坐标

系的方向'''T

r x y z T T T T ??=??。在同一个坐标系有了两个物体的绝对坐标，便很统一得

到两相机的相对位置。

综合上面叙述，最后设两相机的相对位置向量为(),,X Y Z （左相机相对右相机），整理得到利用标靶求出两部固定相机相对位置的数学模型：

2

11121331323321222323132332

1112133132332122232313233(1)(1)''''(1)''''''''(1)''''

(x l z y l z

x r z y r z r X r Y r Z T kr X f r X r Y r Z T r X r Y r Z T kr Y f r X r Y r Z T r X r Y r Z T kr X f r X r Y r Z T r X r Y r Z T kr Y f r X r Y r Z T X Y

Z ++++=+++++++=+++++++=+++++++=+++())'

'

'T

T

T

l r x

y

z x y z T T T T T T T T ????????????????=-=-??

?

（三）

其中符号如前面所示。

四、模型的改进、评价、推广

从上述对模型的研究看，该模型是一个精度很高且具有良好稳定性。利用该模型对标靶点进行图像重塑，可以基本和原图像相差无几。尤其是位于标靶左上方的特征点，该模型对这些区域重构的图像基本和原图像重合。 但是从整体的重构图像看，还是有相当一部分的点的像与原图像发生偏折，使得图像出现重影。这应该和模型一参数确定的过程是密不可分的。一是在整个过程中使用了存在误差的图像中心点坐标，令一方面用来做拟合的数据量只有五个，比起模型的六个未知参数来，数量还很是不足，这就造成了拟合出来的图像在一些地方发生了误差。

为此，我们应该建立更好的图像边界提取算法，使得图像边界归属更加明确，也就使得用于拟合的数据更加精确；再者，利用网格平行线和垂直平行线，建立更多网格特征交点，使得用于拟合的数据更充分，这样就能更大可能的保证了模型的稳定性。

用相机进行双目定位，该方法采用无接触式测量，结构简单、设备成本低、数据采集迅速，操作方便。该模型除了可以用于对物体的定位，同时也还可以用于工程和人体的测量、质量的控制等方面，具有广阔的发展空间。

参考文献：

[1] 王建民等，一种新的二维图像测量系统标定方法，哈尔滨工业大学学报，V ol.32，No.3，37-41，2006年。

[2] 罗万成，大学生数学建模案例精选，成都，西南交通大学出版社，2007年。

[3] 戴朝寿等，数学建模简明教程，北京，高等教育出版社，2007年。

[4] 江世宏，MA TLAB语言与数学实验，北京，科学出版社，2007年。

[5] 张可，基于改进遗传算法的立体视觉系统标定，计算机工程与应用，1-4，2006年。

[6] 李瑛等，基于投影变换的图像测量系统标定，传感器与微系统，77-80，2006年。

[7] 廖士中等，一种光学镜头摄像机图象几何畸变的修正方法，中国图象图形学报，593-596，2000年

附件：

程序1

N=5;

Image_data=imread('xiang.bmp');

%图像的二值化

I=Image_data;

for i=1:768

for j=1:1024

if Image_data(i,j,:)>=125

I(i,j,:)=255;

else

I(i,j,:)=0;

end

end

end

imshow(I);

%图像边缘点

n=0;

k=0;

for i=1:768

for j=1:1024

while ( I(i,j,:)==0 & I(i,j-1,:)~= 0 ) k=k+1;

circle_side.num(1,k)=n;

circle_side.x_label(1,k)=i;

circle_side.y_label(1,k)=j;

end

end

end

程序2

%切割图像平面，得到切割点

Image_datas=imread('xiang.bmp');

Image_data=rgb2gray(Image_datas);

imshow(Image_data)

[m,n]=size(Image_data);

for i=1:m

for j=1:n

if Image_data(i,j)>125

Image_data(i,j)=0;

else

Image_data(i,j)=255;

end

end

end

k=0;

for i=1:m-1

if ( max( Image_data(i,:) )==255 & max( Image_data(i+1,:) )~=255 ) k=k+1;

x=i

end

end

k=0;

for j=1:n-1

if (max( Image_data(:,j) )==255 & max( Image_data(:,j+1) )~=255) k=k+1;

y=j

end

end

程序3

%求切割平面后图像中心，针对A图像

Image_data=imread('cirl1.bmp');

[m,n,z]=size(Image_data);

for i=1:m

for j=1:n

if Image_data(i,j,:)>125

Image_data(i,j,:)=0;

else

Image_data(i,j,:)=255;

end

end

end

Image_binary=rgb2gray(Image_data);

imshow(Image_binary)

%Image_binary=im2bw(Image_binary,[level]);

Image_binary1=dither(Image_binary);

L = bwlabel(Image_binary1);

STATS = regionprops(L,'Centroid')

STATS(1).Centroid

程序4

%模型参数求解

clear

x=[-50,-20,50,-50,50];

y=[50,50,50,-50,-50];

X=[323 424 641 285 584];

Y=[189 197 213 503 504];

X=X-1024/2;

Y=Y-678/2;

x=x.*3.78;

y=y.*3.78;

A=zeros(10,7);

for i=1:5

A(i,1)=1;

A(i,2)=x(i);

A(i,3)=y(i);

% A(i,4)=-(X(i).^2+Y(i).^2)*X(i);

end

for i=6:10

% A(i,4)=-(X(i-5).^2+Y(i-5).^2)*Y(i-5);

A(i,5)=1;

A(i,6)=x(i-5);

A(i,7)=y(i-5);

end

B=[X';Y'];

E=A\B;

A

result=A*E

result=[result(1:5)+1024/2;result(6:10)+678/2];

程序5

%扩充特征点集

clear

Image_datas=imread('xiang.bmp');

Image_data=rgb2gray(Image_datas);

Image_data=Image_data';

X=[323,424,641,285,584];

Y=[189,197,213,503,504];

for i=1:5

for j=X(i):1024

if Image_data(j,Y(1))>=125

X=[X j];

Y=[Y Y(1)];

break

end

end

end

程序6

S=solve('(1-9.3202e-007*(x^2+y^2))*x=-48.637-0.77274*50*3.78+0.11458*50*3.7 8','(1-9.3202e-007*(x^2+y^2))*y=11.426-0.037826*50*3.78-0.76502*50*3.78') S.x

S.y

程序7

%利用模型，重构函数图像

Image_datas=imread('xiang.bmp');

Image_data=rgb2gray(Image_datas);

%imshow(Image_data);

%Image_data=255;

for i=-62:62

for j=-62:62

if ((((i+50)^2+(j-50)^2)-144)<=4 | (((i+50)^2+(j+50)^2)-144)<=4 | (((i-50)^2+(j-50)^2)-144)<=4 | (((i-50)^2+(j+50)^2)-144)<=4

|(((i+38)^2+(j-50)^2)-144)<=4)

m=E(2)*i*3.78+E(3)*j*3.78+E(1);

% 1/(1+E(4)*((i*3.78)^2+(j*3.78)^2))*

n=E(5)*i*3.78+E(6)*j*3.78+E(4);

m=m+1024/2;

n=n+678/2;

if (((i+38)^2+(j-50)^2)-144)<=4

Image_data(fix(n),fix(m)+50)=0;

else

Image_data(fix(n),fix(m))=0;

end

end

end

end

imshow(Image_data)

figure (2)

imshow(Image_datas)