数码相机定位.docx
《数码相机定位.docx》由会员分享，可在线阅读，更多相关《数码相机定位.docx（30页珍藏版）》请在冰豆网上搜索。
数码相机定位
数码相机定位研究
摘要
本文在不考虑相机镜头畸变的情况下，建立了相机成像的中心投影模型，并用搜索矩阵模型来计算靶标上圆的圆心坐标。
在模型实现上，我们进行图像分割处理，将灰度矩阵二值化。
并运用了半径变化法与椭圆拟合模型两种方法来对问题一的模型进行了检验。
最后对空间坐标变换关系进行严格推导，确定两部相机的相对位置关系。
问题一：
我们首先通过matlab对靶标的像进行处理，将其转换为数字化灰度矩阵，并对图像分割处理将灰度矩阵二值化。
通过推导，我们发现椭圆的中心与其外切矩形的中心重合，由此建立了搜索矩阵模型，通过作椭圆的外切矩形，用搜索法提取像点，求得矩形各顶点的坐标，从而得到矩形中心的坐标，即椭圆的中心坐标，通过坐标变换可得靶标的圆心在相机像平面上的坐标。
问题二：
利用问题一的搜索矩阵模型，通过matlab编程，导入相应数值和图像对模型进行求解，得到靶标上圆的圆心A,B,C,D,E在像平面的像坐标分别是（单位：
像素）：
（-189.5,-195），（-90,-186.5），（128.5,-170.5），（71.5,119.5），（-228,119）
问题三：
本文引入两种方法对模型精确度进行了检验。
方法一通过coraldraw绘图软件改变靶标上圆的半径而改变像。
由于圆心位置不变，通过模型求解将半径改变前后圆心的坐标进行对比，并从圆的半径大小角度分析了模型稳定性。
方法二通过编程提取像闭边缘的像素坐标，用matlab对离散点进行椭圆拟合得拟合方程，进而求得其圆心。
两种方法所求得的相对误差都小于2%，表明问题一中模型的精确度很高。
问题四：
本文采用空间坐标变换的方法，首先对世界坐标系、图像坐标系、相机坐标系的变换关系进行推导，然后解线性方程组分别计算两部相机的内外参数，通过坐标变换分别得到了两固定相机的光学中心在世界坐标系中的坐标，由此确定了两相机的相对位置。
关键词：
小孔成像搜索矩阵法半径变化法椭圆拟合模型坐标变换
一、问题重述
数码相机定位在交通监管等方面有广泛的应用。
所谓数码相机定位是指用数码相机摄制物体的相片确定物体表面某些特征点的位置。
最常用的定位方法是双目定位，即用两部相机来定位。
其中精确地确定两部相机的相对位置（系统标定）是关键。
标定的一种做法是：
在一块平板上画若干个点，同时用这两部相机照相，分别得到这些点在它们像平面上的像点，利用这两组像点的几何关系就可以得到这两部相机的相对位置。
然而，无论在物平面或像平面上我们都无法直接得到没有几何尺寸的“点”。
实际的做法是在物平面上画若干个圆（称为靶标），它们的圆心就是几何的点了。
而它们的像一般会变形，所以必须从靶标上的这些圆的像中把圆心的像精确地找到，标定就可实现。
现已知一靶标：
取1个边长为100mm的正方形，分别以四个顶点（对应为A、C、D、E）为圆心，12mm为半径作圆，同时以AC边上距离A点30mm处的B为圆心，12mm为半径作圆。
并用一位置固定的数码相机摄得到它的像。
本文将解决以下问题：
（1）建立数学模型和算法以确定靶标上圆的圆心在该相机像平面的像坐标,这里坐标系原点取在该相机的光学中心，x-y平面平行于像平面；
（2）对给出的靶标及其像，计算靶标上圆的圆心在像平面上的像坐标,该相机的像距（即光学中心到像平面的距离）是1577个像素单位（1毫米约为3.78个像素单位），相机分辨率为1024×768；
（3）设计一种方法对模型进行检验，并对方法的精度和稳定性进行讨论；
（4）建立用此靶标给出两部固定相机相对位置的数学模型和方法。
二、问题分析
问题一：
建立模型确定靶标的圆的圆心在该相机像平面的像坐标。
首先用matlab将靶标的像图转换为数字化灰度矩阵，并进行二值化处理。
不难证明，椭圆的中心与其外切矩形的中心重合，通过作椭圆的外切矩形，用搜索法得到矩形中心的坐标，也就得到了靶标的圆心在相机像平面上的坐标。
问题二：
确定所给靶标的圆的圆心在该相机像平面的像坐标。
将上述问题一的模型通过matlab编程，导入相应数值和图像对模型二进行求解，得到靶标上圆的圆心在像平面上的像坐标。
问题三：
对模型进行检验，并对方法的精度和稳定性进行讨论
可以改变靶标上圆的半径而改变像，同时改变圆的外切矩形各点坐标，将半径改变前后用模型所求圆心的坐标进行对比分析。
同时找寻另一种方法求得圆心坐标，并与之对比得到两种方法的相对精确度。
对于模型稳定性，可从圆的半径大小的角度进行分析。
问题四：
给出两部固定相机相对位置的数学模型和方法。
采用空间坐标变换的方法，得到了两固定相机的光学中心在世界坐标系的坐标，从而确定两相机的相对位置。
三、模型的基本假设
1.假设照相机系统严格按照小孔成像原理工作；
2.假设焦距等于像距；
3.假设透镜畸变对像的形变影响可以忽略。
四、基本变量及符号说明
表示像在相机坐标系中的坐标；
表示物在世界坐标系中的坐标；
表示图像坐标系的像素坐标；
表示光轴与图像平面的交点在图像坐标系中的像素坐标；
表示相机的焦距；
表示旋转矩阵；
表示平移矩阵。
五、模型的建立与求解
5.1成像原理
模型一:
小孔成像
数码相机成像原理与传统相机相似，实际上就是一个光学成像过程，是三维空间中的物体到像平面的投影。
理想的投影成像模型是光学中的中心投影，也称为针孔模型。
针孔模型假设物体表面的反射光都经过一个“针孔”而投影在像平面上。
物点、光心和像点在一条直线上，满足光的直线传播的条件。
这个投影中心称为光心。
在小孔成像中（如图一），焦距f等于光心到像面的距离，u为物距，v为像距,满足
。
u
f=v
图一：
小孔成像模型原理图
对于理想的透镜成像模型，如图二：
图二：
透镜成像模型原理图
其中
满足关系式：
（1）
焦距f不等于像距v,当u>>v时，可取
，此时透镜成像模型与中心投影一致。
于是我们假设
用针孔模型即中心投影作为相机成像模型。
在不考虑各种镜头的畸变的情况下，小孔成像能给实际相机一个很好的近似。
5.2问题一：
确定靶标上圆的圆心在相机像平面的像坐标
5.2.1图像预处理
（1）将题目中给的图3通过画图软件恢复其原始大小（高203.2mm,宽270.9mm），并保存为jpg格式的文件。
（2）用MATLAB对图像进行数字化处理生成灰度矩阵（768行1024列），矩阵的行号和列号分别是其像素坐标。
（3）由于相机摄像清晰度不同，图像边缘会产生模糊现象，采用阈值分割法对图像进行分割处理,即图像的二值化。
通过计算，阀值t接近50，所以用50作为背景与目标的临界，将像素大于50的作为背景，小于等于50的作为目标。
我们将目标赋值为0，对应图像中黑点区域，将背景赋值为1，对应图像中相对白色区域。
由此我们建立像素坐标系如图三所示：
图三：
靶标的像的像素坐标（坐标原点在图的左上角顶点）
5.2.2模型的建立与算法
5.2.2.1模型二：
搜索矩阵模型
由投影原理，不论在什么位置对圆进行摄影其图像不是圆就是椭圆。
考虑到相机所得图像一般会变形，我们视上图中的像均为椭圆。
通过对灰度矩阵二值化后，对此矩阵分析得，为0的元素都集中分布于一定集中的区域并且其区域连起来构成一个椭圆状，与其像闭为椭圆对应。
易知对于任意一个圆，其任意外切矩形的中心与圆的圆心是重合的。
对于椭圆由其对称性知其也有类似的结论，即其任意外切矩形的中心与椭圆的中心重合。
如下图，在椭圆上任取一点E作其外切矩形ABCD，此外切矩形是唯一的，它与椭圆切于E，F，G，H四点。
图四:
椭圆外切矩形示意图
根据椭圆的对称性，可知E点与H点，G点与F点关于椭圆中心对称的。
BE与DH等长，则四边形ADHE与四边形BCHE面积相等。
将ABCD看作质量均匀的平板，由物理中重心的概念，可知四边形ABCD的重心在EH线上，同理四边形ABFG与四边形CDGF面积相等，则ABCD的重心在FG上，所以矩形ABCD的重心在EH与GF的交线上，而EH与GF相交与椭圆的中心，即矩形的重心就是椭圆中心。
因为质量分布均匀的平板其中心即为其几何中心，所以椭圆的中心与矩形的中心重合。
由此，我们只需将图像三中五个椭圆的外切矩形的中心分别找出来，也就找到了各椭圆的中心。
5.2.2.2用搜索法求矩形中心的坐标
我们通过对0—1化后的灰度矩阵用搜索法求得矩形的各个顶点的坐标，那么椭圆的中心坐标可由矩形对角顶点的横纵坐标的平均值得到。
具体步骤如下：
首先对0—1化后的灰度矩阵按图二进行分块，然后逐行逐列搜索，对行找到第一个是0的元素记下其行号为
，对列也如此找到第一个是0的元素记下其列号
。
继续搜索下一行和列，如果有0则记下
，
，当发现这一行或列没有0时停止搜索，并停止记录（这一行或列不记）。
最后统计行向量
，
的最大元素
，
，则椭圆的中心的横、纵坐标为：
（2）
由此，我们得到椭圆的中心在像素坐标系下的坐标，然后转换到图像坐标系（图像中心为原点的坐标系）中求得其在像平面坐标系中的物理坐标。
5.2.3平移转换坐标模型（问题1、问题2）
如图五所示，设像素坐标系为
像平面的坐标系为
图五：
平移转换坐标模型示意图
设
在
中坐标为
因为相机分变率为1024×768，有
。
则可得：
（3）
将5.2.2中所计算出的像素坐标带入（3）式，可求得靶标上圆的圆心在该相机像平面的像坐标，即A、B、C、D、E五点在像平面的像坐标。
5.3问题二：
计算靶标上圆的圆心在像平面上的像坐标
5.3.1模型二的求解：
通过matlab编程，导入相应数值和图像对模型二进行求解，得到了靶标上圆的圆心在像平面上的像坐标（程序详见附录）。
用搜索法求得椭圆中心在像素坐标系中的坐标如表一：
像点
A
B
C
D
E
坐标
（322.5，189）
（422,197.5）
（640.5，213.5）
（583.5,503.5）
（284，503）
表一：
靶标像的中心在像素坐标系中的坐标（单位：
像素）
通过对所求得的像素坐标进行平移转换，得到了椭圆中心在像平面上的像坐标，如表二：
像点
A
B
C
D
E
坐标
（-189.5，-195）
（-90,-186.5）
（128.5，-170.5）
（71.5,119.5）
（-228，119）
表二：
靶标像的中心在像平面坐标系中的物理坐标（单位：
像素）
5.4问题三：
分析问题二结果的精度与稳定性
5.4.1模型二算法误差分析
在相机成像时严格按照小孔成像原理，忽略透镜畸变的影响的情况下，如果相机的分辨率足够高，使得圆周的像可以认为是连续光滑的曲线，那么可以认为模型中确定的切点即为此曲线的切点，此时模型的误差非常小。
所以说模型的误差可以认为完全是由边界点的离散特性造成的。
对模型二精确度的检验，本文引入两种方法：
（一）模型三：
改变半径法
对于模型二所用的方法，可以改变靶标上圆的半径。
改变靶标上圆的半径，像平面上所呈椭圆的像就会发生变化。
根据靶标上圆所有的点都在此圆的像平面上得知：
改变半径之后，靶标上圆的像的形状与圆心在像平面上的位置均没有发生变化，而只有像的大小改变。
改变半径后，程序确定的像的切点的坐标就会变化，由此计算出的圆心的坐标由于误差的存在也会不同。
对于题上图三给出的像，通过coraldraw绘图软件，只改变像内椭圆的轮廓大小，而不改变椭圆的位置和相片的大小，这样就相当于通过改变靶标上圆的半径而改变像。
为简化模型，我们只改变其中一个椭圆的半径（如图六），然后将求得的圆心坐标与模型二的结果相比较，算出他们的相对误差，从而确定模型二的精确度。
图六：
圆A改变半径后示意图（变化1）
对图三给出的像，把圆A进行三次轮廓变化，将其保存为png格式，然后导入MATLAB中，通过模型二的方法输出其图像和灰度矩阵，然后得出圆心的像素坐标。
分别计算半径变化前后圆心坐标的相对误差，得出下列表格:
圆A
半径未变
变化1
变化2
变化3
圆心坐标
（322.5,189）
（319.5,186.5）
（316,177）
（310,170）
相对误差
X轴
0.63%
0.71%
0.97%
Y轴
0.82%
1.24%
1.95%
表三：
改变圆的半径前后圆心坐标的比较
由上表可知，坐标轴的相对误差最大为1.95%，在误差允许的范围内，说明模型是精确的。
（二）模型四：
曲线拟合法
用Matlab将靶标的像，即题中图三读入后，编写程序选取各个椭圆上的边界点，建立椭圆线性回归模型，利用Matlab工具箱中的regress函数，得到五个椭圆的拟合方程，进而可以求得其圆心。
提取边界点的算法流程如下：
（1）逐行逐列扫描数矩阵的零元素
（2）进行判断。
若零元素相邻的上、下、左、右元素中至少有一个为非零元素，那么元素位于椭圆像图边缘，将其坐标提取；若无非零元素，则元素位于椭圆像图内部，将其忽略。
圆或椭圆的二次曲线方程的一般式为：
（4）
且
椭圆的中心
可用下式解求：
（5）
为编程方便，将一般性方程转化为：
（6）
其中
为参数。
椭圆拟合结果及与模型二相对误差见表四：
表四：
椭圆拟合与模型二结果比较
模型二
A
B
C
D
E
圆心坐标
X轴
322.5
422
640.5
583.5
284
Y轴
189
197.5
213.5
503.5
503
椭圆拟合
A
B
C
D
E
圆心坐标
X轴
322.5296
422
640.5502
583.4428
284.3894
Y轴
189.6289
197.481
214.0371
503.8498
502.225
相对误差
X轴
0.0092%
0
0.0078%
0.0098%
0.1371%
Y轴
0.3327%
0.0096%
0.2515%
0.0694%
0.1540%
由上表可知，相对误差最大仅为0.3327%，说明在误差允许范围内，模型二是非常精确的。
5.4.2模型二稳定性分析
模型二的稳定性必须依据对圆所取半径的具体长短来确定，假如所取的半径长度很短，使得圆在像平面上的像只有几个像素点，那么虽然模型中采取二值法消除误差，但是所参与的边界点很少，使得确定的切点将产生比较大的误差，从而由此半径确定的圆心像坐标的偏差比较大。
假如所取的半径长度很长，使产生的像的曲率很小，则边界点可能出现很长一段边界点都在一条直线上的情况，以此半径确定的圆心的像坐标的偏差也比较大。
根据以上分析表明，模型的适用范围为靶标上圆的半径不能太长，也不能太短，此时，按此方法的计算的结果是稳定的。
如果半径取得过长或过短，按此方法的计算的结果将偏差很大，认为此方法是不稳定的。
5.5问题四：
确定两部固定相机相对位置的数学模型
由给出的靶标，得到各圆的圆心的物点坐标和像点坐标，以此求出相机的内外参数，从而分别确定两部相机的光学中心在世界坐标系下的坐标，它们的相对位置就确定了。
此过程需要建立世界坐标系并进行空间坐标变换。
图七：
单个相机成像空间坐标关系示意图
5.5.1空间坐标的变换
（1）世界坐标系（
）到相机坐标系（
）的变换
两个坐标系都是三维直角坐标系，由世界坐标系中的点到相机坐标系的变换表示为:
（7）
其中
为正交旋转矩阵，其矩阵元素满足：
（8）
为平移变换矩阵。
（1）式用齐次坐标表示为：
（9）
式中R只含3个独立变量，再加上T中
，这六个参数决定了相机光轴在世界坐标系中空间位置，被称为相机外部参数。
（2）相机坐标系（
）与图像坐标系（
）的变换关系
相机坐标系中的物点P在图像物理坐标系中像点
坐标为：
其中f是相机的焦距，这里我们认为焦距等于像距，
。
齐次坐标表示为：
（10）
将上式进一步转化为：
（11）
齐次坐标表示为：
（12）
其中
是图像中心（光轴与图像平面的交点）坐标，
分别为一个像素在X与Y方向上的物理尺寸，
分别为X与Y方向上的采样频率，即单位长度的像素个数。
由此可得物点P与图像像素坐标系中像点P的变换关系
（13）
其中，
分别为X和Y方向的等效焦距。
四个参数只与摄像机内部结构有关，因此称为相机内部参数。
（3）世界坐标系与图像坐标系变换关系（共线方程）
世界坐标系与图像坐标系变换关系：
齐次坐标表示为：
（14）
（14）式是摄影测量学中最基本的共线方程，说明物点、光心和像点这三点必须在同一条直线上。
根据共线方程，在摄像机内部参数确定的条件下，利用若干个已知的物点和相应的像点坐标，就可以求解出摄像机的六个外部参数，即摄像机的光心坐标和光轴方位的信息。
5.5.2解线性方程组求两相机几何关系R，T
由于光轴穿过图像中心，单位长度的像素个数为3.78，（14）式可表示为：
（15）
方程（15）描述了三维世界坐标点（
，1）与相应图像点（u,v,
）之间的关系。
整理消去
后，得到如下线性方程：
因为世界坐标系的X-Y平面与物体所在平面坐标系重合，即
=0，所以它的系数对结果不影响可设为0，上式转化为：
（16）
上式对于任意一个世界坐标系中的点都成立，由题中给出的靶标，可得5个圆心的物点坐标和像点坐标，将这5组数据带入上式，可以得到10个线性方程，联立（8）式的三个方程，解方程组计算出该相机的外部参数
，从而确定第一个相机的光学中心在世界坐标系下的坐标。
按照同样的算法，确定第二个相机的光学中心在世界坐标系下的坐标，这样相对位置就出来了。
六、模型的评价与改进
（一）模型的优点：
1.照相机的成像原理由透镜成像近似为小孔成像，即中心投影。
小孔成像模型原理为光的直线传播，大大简化的做题的难度。
2.模型二对题上图3进行二值化处理，解决了边缘灰色区域属于椭圆外或内的问题，方便提取边缘点。
然后将图3进行分块，逐行逐列搜索切点进而求矩形中心，方法精确度较高。
3.用半径放大法和椭圆曲线拟合法分别检验模型二的精度，检验的结果具有说服力。
且半径放大法巧妙地将题中所给的图3中A圆的轮廓进行放大，而不改变圆心位置，避免了放大靶标上的A圆然后进行拍照获得像图的繁琐环节。
4.用空间坐标变换法确定两相机相对位置，方法的思路简单，并且理论推导严谨。
（二）模型的缺点：
本文模型没有考虑相机拍摄时的透镜畸变情况，使得结果存在一定的误差。
在对模型二稳定性进行分析时，只是提出了分析思路，并没有通过具体运算得出数据。
（三）模型的改进方向:
在第二问中算出像圆的圆心后，可以对透镜畸变引起的误差进行计算，然后得出修正后的圆心坐标，这样结果会更精确。
参考文献:
[1]姜启源，谢金星，叶俊数学模型（第三版）北京，高等教育出版社，2003
[2]郑阿奇主编，曹弋编著.MATLAB实用教程（第2版），北京：
电子工业出版社，2004,108～109。
[3]姜大志等，数码相机标定方法研究，南京航空航天大学学报，Vol.33,No.1,P55~59.2001.
[4]夏菁，椭圆拟合方法的比较研究，,2012/8/14
[5]刘文彬，周宏甫，蒋梁中基于双目立体视觉系统的精确标定算法，自动测量与控制，78-80页，2008年第27卷第6期
附录：
附录1：
模型二求解圆心主程序
clear;clc;clearall;
[x,cmap]=imread（'d:
\pic.jpg'）;
image（x）
colormap（cmap）
xlabel（'x'）;
ylabel（'y'）;
title（'thefigureofimage'）;
text（300,560,'E'）;
text（250,200,'A'）;
text（500,200,'B'）;
text（700,220,'C'）;
text（600,560,'D'）;
grid;
x1=x（1:
384,1:
1024）;
x2=x（385:
768,1:
1024）;
fori=1:
384
forj=1:
1024
ifx1（i,j）>50
x1（i,j）=0;
end
ifx2（i,j）>50
x2（i,j）=0;
end
end
end
fori=1:
384
forj=1:
370
ifx1（i,j）~=0
pa（:
i）=i;
qa（:
j）=j;
end
end
end
pa1=max（pa）;
qa1=max（qa）;
fori=1:
pa1
ifpa（:
i）~=0
pa2=pa（:
i）;
break
end
end
fori=1:
qa1
ifqa（:
i）~=0
qa2=qa（:
i）;
break
end
end
ACX=（qa1+qa2）/2
ACY=（pa1+pa2）/2
fori=1:
384
forj=371:
500
ifx1（i,j）~=0
pb（:
i）=i;
qb（:
j）=j;
end
end
end
pb1=max（pb）;
qb1=max（qb）;
fori=1:
pb1
ifpb（:
i）~=0
pb2=pb（:
i）;
break
end
end
fori=1:
qb1
ifqb（:
i）~=0
qb2=qb（:
i）;
break
end
end
BCX=（qb1+qb2）/2
BCY=（pb1+pb2）/2
fori=1:
384
forj=500:
1024
ifx1（i,j）~=0
pc（:
i）=i;
qc（:
j）=j;
end
end
end
pc1=max（pc）;
qc1=max（qc）;
fori=1:
pc1
ifpc（:
i）~=0
pc2=pc（:
i）;
break
end
end
fori=1:
qc1
ifqc（:
i）~=0
qc2=qc（:
i）;
break
end
end
CCX=（qc1+qc2）/2
CCY=（pc1+pc2）/2
fori=1:
384
forj=1:
500
ifx2（i,j）~=0
pe（:
i）=i+384;
qe（:
j）=j;
end
end
end
pe1=max（pe）;
qe1=max（qe）;
fori=1:
pe1
ifpe（:
i）~=0
pe2=pe（:
i）;
break
end
end
fori=1:
qe1
ifqe（:
i）~=0
qe2=qe（:
i）;
break
升级会员 