数学建模的格式

承诺书我们仔细阅读了中国大学生数学建模竞赛的竞赛规则.我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。

我们知道，抄袭别人的成果是违反竞赛规则的, 如果引用别人的成果或其他公开的资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。

我们郑重承诺，严格遵守竞赛规则，以保证竞赛的公正、公平性。

如有违反竞赛规则的行为，我们将受到严肃处理。

我们参赛选择的题号是（从A/B/C/D中选择一项填写）：我们的参赛报名号为（如果赛区设置报名号的话）：所属学校（请填写完整的全名）：参赛队员(打印并签名) ：1.2.3.指导教师或指导教师组负责人(打印并签名)：日期：年月日赛区评阅编号（由赛区组委会评阅前进行编号）：编号专用页赛区评阅编号（由赛区组委会评阅前进行编号）：全国统一编号（由赛区组委会送交全国前编号）：全国评阅编号（由全国组委会评阅前进行编号）：双目高清定位摘要:本文是研究数码相机定位问题，先后建立了圆心定位的模型与算法和数码相机标定的模型．对于问题一，考虑只有一个圆的圆心在像平面的定位情况．首先，建立相机坐标系,利用边缘抽取法，得到圆的边缘像素坐标集．然后基于广义霍夫变换，建立圆心定位算法，求得峰值点坐标，作为圆心在像平面的像坐标．为求解问题二，先用sobel 算子边缘抽取法获取五个圆的边缘像素坐标集．采取聚类分割方法对其加以分离，获得五个圆各自的边缘像素坐标集．沿用问题一的模型，利用Matlab 软件分别求得五个圆的圆心在像平面的像坐标为 (-49.47,-51.06,417.20)，( -22.75,-49.47,417.20)，( 34.39,-45.24,417.20)，( 19.05,31.22,417.20)，(-59.52,31.22,417.20)． 问题三要求设计一种模型检验方法,由于圆A,B,C 的圆心在同一直线上,所以可检验这三个圆心的像坐标是否保持在同一直线上,以此检验模型求解结果的精度. 基于不同的检测算子边缘抽取法,分别求解出相应的结果,并利用这种检验方法求得各个结果的精度.根据精度的偏差,可以分析该检验方法的稳定性. 问题四即是相机标定的问题，先以靶标所在平面建立世界坐标系，求得靶标五个圆的圆心以及AD 与CE 之间的连线交点F 的世界坐标．通过问题1的模型，求得五个圆心的像坐标，另外，A ’D’与C’E’之间的连线交点的坐标即是点F 在像平面的像坐标．然后利用世界坐标系下的空间点坐标以其相应的像坐标之间的关系：11121314212223243132333411w i i w i ci i w iX m m m m u Y z v m m m m Z m m m m ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦， 求得相机的内外参数．然后通过空间任一点P 分别在两部相机中的坐标关系：1111221122C C X R R X t R R t --=+-，得到这两部相机的相对位置关系．关键词:检测算子 边缘抽取法 霍夫变换 聚类分割法 相机标定一、问题重述数码相机定位是指用数码相机摄制物体的相片确定物体表面某些特征点的位置．最常用的定位方法是双目定位．对物体上一个特征点，用两部固定于不同位置的相机摄得物体的像，分别获得该点在两部相机像平面上的坐标．只要知道两部相机精确的相对位置，就可用几何的方法得到该特征点在固定一部相机的坐标系中的坐标，即确定了特征点的位置．于是对双目定位，精确地确定两部相机的相对位置就是关键，这一过程称为系统标定．标定的一种做法是：在一块平板上画若干个点，同时用这两部相机照相，分别得到这些点在它们像平面上的像点，利用这两组像点的几何关系就可以得到这两部相机的相对位置．然而，无论在物平面或像平面上我们都无法直接得到没有几何尺寸的“点”．实际的做法是在物平面上画若干个圆（称为靶标），它们的圆心就是几何的点了．而它们的像一般会变形，如图1所示，所以必须从靶标上的这些圆的像中把圆心的像精确地找到，标定就可实现．图 1 靶标上圆的像现在，有人设计靶标如下，取1个边长为100mm的正方形，分别以四个顶点（对应为A、C、D、E）为圆心，12mm为半径作圆．以AC边上距离A点30mm处的B为圆心，12mm为半径作圆，如图2所示．用一位置固定的数码相机摄得其像，如图3所示．图3 靶标的像下面我们需解决一下4个问题：（1）建立数学模型和算法以确定靶标上圆的圆心在该相机像平面的像坐标,这里坐标系原点取在该相机的光学中心，x-y平面平行于像平面；（2）对由图2、图3分别给出的靶标及其像，计算靶标上圆的圆心在像平面上的像坐标, 该相机的像距（即光学中心到像平面的距离）是1577个像素单位(1毫米约为3.78个像素单位)，相机分辨率为1024×768；（3）设计一种方法检验你们的模型，并对方法的精度和稳定性进行讨论；（4）建立用此靶标给出两部固定相机相对位置的数学模型和方法．二、符号说明和基本假设2.1基本假设(1)用相机拍摄相片不会出现叠影或曝光残缺.(2)圆圈上各像素点的灰度值相差不大,相邻像素点的灰度值相差不大.(3)不考虑边缘抽取的像素点残缺较多且分散.2.2符号说明三、模型准备为了便于建立模型和求解问题,先建立像平面的坐标系.3.1坐标系的建立[1]：如图4所示，在像平面上定义直角坐标系0O uv ，每一像素点的坐标（u ，v ）分别是该像素点在数组中的列数与行数，而（u ，v ）是以像素为单位的图像坐标称为像素坐标；像素坐标的集合称为像素坐标集．由于(u,v)只表示像素位于数组中的列数与行数,并没有用长度单位表示出该像素在图像中的位置，因而，为了解题的需要，我们建立一个以物理单位(毫米)表示的像平面坐标系．如图4，坐标原点O 1为相机光轴与像平面的交点，建立与u,v 轴平行且同向的x ，y 轴，并以毫米为单位的1o xy 坐标系．若（u ，v ）表示像素坐标，（x ，y ）表示以毫米为单位的像坐标系的坐标． x ，y 轴方向上的物理尺寸为dx ，dy ，则图像中任意一个像素在两个坐标系下的坐标有如下关系：x u u d x y v v yd=+=+． 用齐次坐标与矩阵形式表示上式为：00101011001u dx u x v v y dy ⎡⎤⎢⎥⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎢⎥⎢⎥⎣⎦， 逆关系可写成：0000.1101dxu dx x uy dy v dy v ⎡⎤-⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥=-⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦图4 像平面坐标系四、问题1的求解4.1 问题1的分析：首先,为了便于求解,根据题目要求,建立一个相机坐标系.接着,要确定圆心的像坐标,需要建立模型,确定圆心的像坐标与圆面的其他像素点的关系. 边缘是图像最基本的特征 [2] ,圆的边缘点与圆心的关系最为紧密.所以需要建立利用边缘点定位圆心坐标的模型.然后,要确定圆的边缘点位置,可用检测算子的边缘抽取法,获得边缘点的像素坐标集. 最后,根据模型的算法,利用边缘点的像素坐标集,计算出圆心的像坐标.4.2 问题1的模型建立: 4.2.1 建立三维坐标系首先，我们根据题目要求建立相机坐标系，如图5．也就是说，以光学中心为原点，过原点垂直于像平面并以指向像平面为正方向的有向直线作为Z c 轴．可知，Z c 轴经过像平面坐标系的原点O 1 ,X c -Y c 平面平行于像平面，其中X c 轴，Y c 轴分别平行于x 轴，y 轴且同向．4.2.2 图像空间中圆心的GHT(广义霍夫变换)求解基本原理[3]要从一幅图像中检测出半径已知的圆形来.取和图像平面一样的参数平面，以图像上每一个前景点为圆心，以已知的半径在参数平面上画圆面，并把结果进行累加.最后找出参数平面上的峰值点，这个位置就对应了图像上的圆心，其中图像平面上的每一点对应到参数平面上的一个圆面.若图像空间中圆的方程为：2200()()x x y y r -+-=.对应于图像空间，参数空间中以(x ，y)为圆心，半径为R 的圆方程为：22()()u x v y r -+-=.在圆的半径r 为常量的前提下，GHT 是将图像空间中圆边界上的一个点(x ，y)映射参数空间(u ，v)一个圆上.GHT 求取圆心坐标的算法步骤简述如下： (1)建立起一个离散的参数空间(u ，v)，参数空间的量化间隔取决于所求圆心坐标参数的精度. (2)建立一个二维累加器N[f,g]，并置每个元素初始值为0.(3)对参数空间的每个点(u,v)，考虑图像空间中每个边界像素点(x,y)，若222()()u x v y r -+-≤图5 相机坐标系成立，则相应的累加器N[f,g]加1.(4)找出N[f,g]的极大值，所求的圆心参数坐标为(f,g).4.2.3 基于GHT的圆心定位模型对于靶标上的圆,用相机摄得其像,在像平面上不再是一个规则的圆.所以,无法简单的用广义霍夫变换求解圆心.而且,对于圆的像,一般会变形,其边界像素点也不能精确的确定.边缘是图像最基本的特征[2],基于边的分割方法,找出灰度不连续的点作为边缘.为了在图像中确定边缘点,已经有很多人设计了各种边缘检测算子,然后根据测出的边缘点作为边缘.其中Matlab软件自带的算子检测函数有sobel算子, 梯度算子, prewitt算子, log算子, canny 算子.各种检测算子的边缘抽取法,在图像上所确定的边缘点会出现重叠或缺损.如果利用所有边缘点的几何关系,来确定图像的质心,误差较大,且抗干扰性差.所以,利用边缘点的几何关系确定圆心的方法不可取.如果利用图像的所有像素点来确定图像的质心,计算量非常大,不易实现.无论相机从哪个角度摄得圆的像, 靶标上圆内的每一小块(大小取决于摄像时的像距和相机的分辨率)都会在图像的像素点上获得一个灰度值.[4]虽然圆的像在像平面上不是一个规则的圆,但是其像素点上的灰度值并不相等.而且,对于一个平面图形的像,相邻像素点的灰度值相差不大[2].所以,圆心沿着半径到边缘上各像素点的灰度值是连续(在这里,连续指灰度值相差不大)变化的.而且, 边缘上各像素点的灰度值两两彼此相差不大.这样,圆心在图像上的定位受到所有边缘点共同的限制.基于广义霍夫变换,利用边缘点定位圆心的算法如下:(1)给像平面的每个像素点建立一个计数器,初始值为0.(2)以检测算子抽取出来的边缘像素点(x,y)为中心,在像平面上画圆.(3)对于像平面的每个像素点(u,v),如果满足222-+--≤.u x v y rε()()ε>为预先设置的误差值,r为半径.则像素点(u,v)的计数器N[u,v]的值增加１,其中0(4)以每个边缘像素点为中心画圆,重复(2)(3)步.(5)找到N[u,v]的极大值,其对应的像素点(u,v)则为圆心的定位.在该模型的求解中, ε取值近似于0；r为靶标上圆的半径在图像中的长度,其单位为像素单位.而且,边缘点的重叠和缺损,对模型的求解影响较小,同样对r的取值精度要求也不需要太高.五、问题2的求解5.1问题2的分析：问题2的求解可沿用问题１的模型.第一，要求解圆心的像坐标,则需要先抽取出圆的边缘坐标集.我们可用Matlab软件读取图像的数据,并用sobel检测算子抽取图像的边缘坐标集.第二，要利用边缘坐标集给圆心定位,则需先把五个圆的图像的边缘像素坐标集分割出来,分别放在五个边缘像像素坐标矩阵中. 由于五个圆在像平面所成的像之间的欧式像素距离足够大，所以可用聚类的方法加以分割，获得五个圆在像平面上各自的边缘像素坐标集．第三，由于边缘抽取法得到的图形不一定是个规则的椭圆,而且有重叠和残缺.所以，利用边缘点的几何关系来确定圆心在像平面的像坐标的方法不可取.因此，我们采用基于广义霍夫变换的算法,即问题1的圆心定位模型．5.2 问题2的模型求解：5.2.1 坐标系的建立：我们先建立坐标系如图6,与问题一的相机坐标系和像平面坐标系一致.Z c 轴垂直于像平面，且光学中心到像平面的距离固定.因此,点在Z c 轴方向的坐标固定为417.20毫米．对于点的坐标，只需要确定X c 轴，Y c 轴方向的坐标，也就是要找出点在像平面1o xy 坐标系的坐标．5.2.2 边缘像坐标集的抽取和分割：用sobel 检测算子抽取边缘，得到五个圆在像平面上图形的边缘像坐标集．由于五个圆的像之间的欧式像素距离足够大，并且同一个圆的边缘像素点非常接近.因此,可以用如下方法对边缘像坐标集进行分割.对于任意两个边缘像素点11(,)u v ,22(,)u v ,如果满足1212u u v v αα⎧-≤⎪⎨-≤⎪⎩， 则将这两个边缘像素点,放进同一个边缘像坐标矩阵中,其中α一般取值为1.如果像素点出现残缺,可能一部分像素点被遗漏,但对圆心定位的结果几乎没什么影响.根据假设3, 不考虑边缘抽取的像素点残缺较多且分散的情况.用上述方法将所有抽取出的边缘像素点类聚，获得五个圆在像平面上各自的边缘像素坐标集．5.2.3 圆心的定位：根据5.1的分析,可利用基于广义霍夫变换定位圆心的算法(如4.2.3),用Matlab 软件求得结果如下：图6 相机坐标系由图7观察可知，圆心的定位比较准确．而且，表1的结果比较符合实际．六、问题３的求解6.1 问题3的分析：靶标上圆A,B,C 的圆心在同一直线上,用相机摄得其像后,其像也依然要保持直线.所以,可以检验这三个圆心的像坐标是否保持在同一直线上.并且,可以以此衡量模型求解结果的精度.基于不同的检测算子边缘抽取法,分别求解出相应的结果,并利用这种检验方法求得各个结果的精度.根据精度的偏差,可以分析该检验方法的稳定性.6.2 问题3的求解:6.2.1模型检验及误差衡量:一条直线,用相机拍摄后,其像依然是一条直线.在靶标上,圆A,B,C 的圆心在同一条直线上.因此, 这三个圆心的像素点坐标(,)a a A u v ,(,)b b B u v ,(,)c c C u v 应该保持在同一直线上,即AB 的斜率A B k 和AC 的斜率A C k 应该相等.但事实上,圆心的定位有误差,而且像素点坐标系是一个离散的整数点坐标系,因此A Bk 图7 sobel 算子边缘抽取及圆心定位结果和A C k 计算结果必有偏差.设AB 与AC 的夹角为β,满足关系式tan 1AB AC AB ACk k k k β-=+,则可利用A B k 和A C k 的计算结果求得夹角β为arctan1A B A C A B A Ck k k k β-=+.β的值越小,则说明了模型求解的误差越小,精度越高.用Matlab 软件求解得, AB 与AC 的夹角为00.84β=,误差非常小.6.2.2稳定性的分析:分别用Matlab 软件自带的算子检测函数sobel 算子, 梯度(Roberts)算子, prewitt 算子, log 算子, canny 算子抽取边缘像素点.基于不同的边缘抽取法,用基于广义霍夫变换的圆心定位和误差衡量,得到不同的结果,以此进行稳定性的分析.计算结果如下干扰性较好.其检验结果的误差较小,而且衡量误差的β值基本保持不变,说明该检验方法的稳定性较好.七、问题4的求解：7.1 问题分析：精确地确定两部相机的相对位置，这一过程称为系统标定．标定的一种做法是：在一块平板上画若干个点， 同时用这两部相机照相，分别得到这些点在它们像平面上的像点，利用这两组像点的几何关系就可以得到这两部相机的相对位置． 7.2 模型假设：在标定的过程中，采用小孔相机模型，忽略相机的非线性畸变．7.3 模型建立：首先，用两部相机同时观测靶标，在标定中，我们可以用单相机标定的方法分别得到，两部相机各自的内外参数，如果外参数分别用1R 、1t 与2R 、2t 表示，则1R ，1t 表示1C 相机与世界坐标系之间的相对位置，2R 、2t 表示2C 相机与世界坐标系之间的相对位置．对任意一点P ，在世界坐标系、1C 坐标系与2C 坐标系下的非齐次坐标分别为12,,w C C X X X ，则111222C w C w X R X t X R X t =+=+ （1）将上式中w X 消去后得到1111221122C C X R R X t R R t --=+- （2）因此，两部相机之间的几何关系可用以下R 和t 表示：11211122R R R t t R R t --==- （3）7.4 模型求解：第一步：建立世界坐标系以靶标的中心为原点，x,y 平面与靶标所在平面重合，背向相机的那一端为z 轴的正方向．如图所示：依据题目所给数据，求出圆心点A ，B ，C ，D ，E 的世界坐标值，以及AD 与CE 之间的连线交点F 的世界坐标值，(),,1,2,6wi wi wi x y z i = ．第二步：求相机内外参数[１]首先，根据问题1中建立的模型，求解靶标上五个圆的圆心在该相机像平面的像素坐标(),1,2,5i i u v i = ．因为点F 是AD 与CE 之间的连线交点，所以点F 在图像坐标中所对应的点'F 也是''A D 与''C E 连线的交点．因此，利用两直线交点的计算公式可以求得点'F 的图像坐标()66,u v ．因为，世界坐标系表示的空间P 点坐标(),,w w w x y z 与其投影点p 的图像平面坐标(),u v 之间存在着如下关系：图8 世界坐标系图9 相机坐标系和世界坐标系00001210000100000110010100100000100101W xW c T y W W x Wy W WT W u X d f u R t Yz v v f d Z X a u R t Y a v M M X M XZ⎡⎤⎢⎥⎡⎤⎡⎤⎢⎥⎡⎤⎢⎥⎢⎥⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎢⎥⎢⎥⎢⎥⎣⎦⎢⎥⎢⎥⎣⎦⎣⎦⎢⎥⎣⎦⎡⎤⎡⎤⎢⎥⎢⎥⎡⎤⎢⎥⎢⎥⎢⎥===⎢⎥⎢⎥⎢⎥⎣⎦⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦（4）其中，c z 是空间点P 在相机坐标系下的ｚ值，f 是相机的焦距；/,x x a f d =/y y a f d =（,x y d d 是每一个像素分别在,x y 轴方向上的物理尺寸，单位为毫米）；M 为3×4矩阵，称为投影矩阵；1M 为相机的内部参数，2M 为相机的外部参数． 所以，每个靶标上的点及其相应的投影点之间坐标关系可写成：11121314212223243132333411w i i w i ci i w iX m m m m u Y z v m m m m Z m m m m ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦（5） 其中，(),,,1w i w i w i x y z 为空间第i 个点的坐标；(),,1i i u v 为第i 个点的图像坐标；ij m 为投影矩阵M 的第i 行j 列元素．式（5）包含三个方程:111213142122232431323334ci i w i w i w i ci i w i w i w i ci w i w i w i z u m X m Y m Z m z v m X m Y m Z m z m X m Y m Z m =+++=+++=+++ （6） 将式(6)中的第一式除第三式，第二式除第三式分别消去ci z 后，可得如下两个关于ij m 的线性方程：11121314313233342122232431323334w i w i w i i w i i w i i w i i w i w i w i i w i i w i i w i i m X m Y m Z m u m X u m Y u m Z u m m X m Y m Z m v m X v m Y v m Z v m +++---=+++---= （７）将6 个圆心及其投影点之间的关系可以写成12个关于M 矩阵的线性方程，下面用矩阵形式写出这些方程：111211111111113141111111112122236661616162431111111111323310000000011000000001w w w w w w w w w w w w w w w w w w w w w w w w m mX Y Z u X u Y u Z m m X Y Z v X v Y v Z m m m X Y Z u X u Y u Z m m X Y Z v X v Y v Z m m ⎡⎢⎢---⎡⎤⎢⎢⎥⎢⎢⎥---⎢⎢⎥⎢⎢⎥⎢⎥⨯⎢⎥⎢⎥---⎢⎥⎢⎥⎢⎥---⎣⎦⎣ 134134634634i u m v m u m u m ⎡⎤⎤⎢⎥⎥⎢⎥⎥⎢⎥⎥⎢⎥⎥⎢⎥⎥⎢⎥⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎦⎢⎥⎣⎦（8） 由式（4）可见, M 矩阵乘以任意不为零的常数,不影响(),,w w w x y z 与(),u v 的关系，因此，在式（8）中可以指定34m ＝1，从而得到关于M 矩阵其他元素的12个线性方程，这些未知元素的个数为11个 ，记为11维向量ｍ，将式（8）简写成K ｍ＝U （9）其中，K 为式（8）左边12×11矩阵；ｍ未知的11维向量；U 为式（8）右边的12维向量；K 、U 为已知向量．因为12＞11，所以可以用最小二乘法求出上述线性方程的解为：T-1Tm =(K K )K U ． (10)ｍ向量与34m ＝1构成了所求解的M 矩阵．但是现在所求解的M 矩阵与式（4）所表示的矩阵M 相差一个常数因子34m ．这一点可以从解方程（8）时指定34m ＝1中可以看出．虽然指定34m ＝1不影响投影关系，但在分解M 矩阵时必须考虑．我们将式（4）中的M 矩阵与相机内外参数的关系写成1114023422403300001001001T x Tx T y T y T zT Tr t m m a u r t m m m a v r t m ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎢⎥⎣⎦⎣⎦（11） 其中，(1~3)Ti m i =为由式（8）求得的M 矩阵的第i 行的前三个元素组成的行向量；4(1~3)i m i =为M 矩阵的第i 行第4列元素；(1~3)Ti r i =为旋转矩阵R 的第i 行；,,x y zt t t 分别为平移向量t 的三个向量．比较上式两边可得，3433m r m =，由于3r 是正交单位矩阵的第三行，3r 1=，因此，由343m 1m =求得3431m m =．再由以下式子可以求得12300r ,r ,r ,,,,,,,x y x y z u v a a t t t ．33432T 034132T 034232341323423r m m m m m m m m m x y m u m v m a m a m====⨯=⨯ ， 341103342203343414034240r (m m )r (m m )()()x yx y x z ym u a m v a t m m t m u a m t m v a =-=-==-=-．至此，第一部相机的外部参数11,R t 已经求解出来了．同理，可求得第二部相机的外部参数22,R t ．第三步：求解两部相机的相对位置 将1122,,,R t R t 代入式（2）1111221122C C X R R X t R R t --=+-这样，两部相机之间的几何关系可用以下R 和t 表示：11211122R R R t t R R t --==-八、模型的分析与推广(1)问题一和问题二的圆心定位模型和算法克服了图像变形及边缘缺损的影响，抗干扰性强．(2)该模型的圆心定位算法，基于抽取的边缘点，减少了计算量，有利用计算机实现． (3)圆心定位模型的算法精度较高，且计算时间较短，可推广用于快速报靶或无线传感器的图像传送．(4)在问题4的求解中，用式子(10)所表示的线性方法求解这些参数时，并没有考虑这些变量间的约束条件，因此，当数据有误差的情况下，计算结果是有误差的，而且误差在个参数间的分配也没按它们之间的约束关系考虑．如果对M 矩阵进行适当的约束，可提高计算精度．参考文献[1] 马颂德、张正友，计算机视觉，北京：科学出版社，1998．[2] 管伟光，体视化技术及其应用，北京：电子工业出版社，1998．[3] 张红民、何健鹰，用改进的广义Hough 变换获取靶纸图像子像素圆心坐标，计算机与现代化，2003,10:43-50，2003．[4] 曹建海、路长厚，基于灰度图像和矢量的圆心定位，光电子•激光，2004,15(6):714-718,2004．九、附录：附录9.1 问题3根据不同边缘抽取法的圆心定位结果图10 sobel算子边缘抽取及圆心定位结果图11 Roberts算子边缘抽取及圆心定位结果图12 prewitt算子边缘抽取及圆心定位结果图13 log算子边缘抽取及圆心定位结果图14 canny算子边缘抽取及圆心定位结果附录9.2求解问题1和问题2的M atlab程序建立M文件bianyuan.m，程序如下：。