厦大数学建模2008 A题论文

承 诺 书

我们仔细阅读了中国大学生数学建模竞赛的竞赛规则.

我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。

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

我们郑重承诺，严格遵守竞赛规则，以保证竞赛的公正、公平性。如有违反竞赛规则的行为，我们将受到严肃处理。

我们参赛选择的题号是（从A/B/C/D中选择一项填写）： A 我们的参赛报名号为（如果赛区设置报名号的话）： 所属学校（请填写完整的全名）： 厦门大学 参赛队员 (打印并签名) ：1. 许世杰

2. 刘辰晨 3. 蔡鸿祥 指导教师或指导教师组负责人 (打印并签名)： 谭忠

日期： 2008 年 9 月 22 日

赛区评阅编号（由赛区组委会评阅前进行编号）：

2008高教社杯全国大学生数学建模竞赛

编 号 专 用 页

评 阅 人 评 分 备 注 赛区评阅编号（由赛区组委会评阅前进行编号）：

赛区评阅记录（可供赛区评阅时使用）：

全国统一编号（由赛区组委会送交全国前编号）：

全国评阅编号（由全国组委会评阅前进行编号）：

A题 数码相机定位

摘要

本文主要针对数码相机的空间定位、相机的标定、空间点在二维和三维互相转换的问题，提出了以下四个模型以及改进算法：

模型一：是通过类比于质心的求法，对图像进行二值化和分析，提取像素点的坐

xnsumx/npix标和像素的个数，根据公式ynsumy/npix

可得圆心像坐标(xn,yn)。代入第二问数据就可以求得各圆心在像平面的坐标为： 像平面坐标系（单位：像素） 以光学中心为原点，xy平面平行于像平面的坐标系 （单位：毫米） A（322.033，188.550） A(-50.255174,-54.087249,417.195767) B（422.001，195.942） B(-23.809206,-52.131667,417.195767) C（633.888，212.201） C(33.568280,-47.830370,417.195767) D（581.733，502.050） D(18.447857,28.849259,417.195767) E（283.286，500.613） E(-60.506323,28.468915,417.195767) 经过对投影前后在同一条直线A，B,C三点的交比不变以及他们是否仍然在同一条直线来进行的分析，准确性为95％，可认为模型相对准确。经过对误差的分析，发现模型相对较为稳定。 模型二：通过利用模型一得出的点的坐标在三个坐标系之间的相互转换求出投影矩阵M，M(KTK)1KTU从而得到数码相机的各种内外参数。

模型三：采用具有反馈机制、并且可以进行自我学习与调整的神经网络作为分析相机相对位置的算法。并根据其不足提出了神经网络的改进算法：动量－自适应学习率BP神经网络的权值修正算法

(k1)(k)(k)[(1mc)D(k)mcD(k1)]，使得整个模型学习速

率更低，得到的结果也更加准确可信。

模型四: 利用遗传算法在解决多参数,非线性的复杂函数的优化问题方面的优势并结合matlab编程，求出数码相机外参数，进而建立两部固定相机相对位置的数学模型，其中以残差f()为目标函数，

但数据的可依赖度很高。

1N[(XwiX'wi)2(YwiY'wi)2(ZwiZ'wi)2]Ni1 为适应度函数。该算法的收敛速度虽然较慢，

关键词 ：二值化 系统标定 神经网络 遗传算法

1

一、问题重述

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

标定的一种做法是：在一块平板上画若干个点， 同时用这两部相机照相，分别得到这些点在它们像平面上的像点，利用这两组像点的几何关系就可以得到这两部相机的相对位置。然而，无论在物平面或像平面上我们都无法直接得到没有几何尺寸的“点”。实际的做法是在物平面上画若干个圆（称为靶标），它们的圆心就是几何的点了。而它们的像一般会变形，如图1所示，所以必须从靶标上的这些圆的像中把圆心的像精

确地找到，标定就可实现。 图 1 靶标上圆的像 有人设计靶标如下，取1个边长为100mm的正方形，分别以四个顶点（对应为A、C、D、E）为圆心，12mm为半径作圆。以AC边上距离A点30mm处的B为圆心，12mm为半径作圆，如图2所示。

图 2 靶标示意图

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

2

图3 靶标的像

请你们：

（1） 建立数学模型和算法以确定靶标上圆的圆心在该相机像平面的像坐标,

这里坐标系原点取在该相机的光学中心，x-y平面平行于像平面；

（2） 对由图2、图3分别给出的靶标及其像，计算靶标上圆的圆心在像平面上

的像坐标, 该相机的像距（即焦点到像平面的距离）是1577个像素单位(1毫米约为3.78个像素单位)，相机分辨率为1024×768；

（3） 设计一种方法检验你们的模型，并对方法的精度和稳定性进行讨论； （4） 建立用此靶标给出两部固定相机相对位置的数学模型和方法。

二、问题分析

本问题涉及到的是由三维空间映射得到的二维图像经处理后还原为三维空间立体图形的问题。

这四个问题是相互关联的，第一问是建立模型来求解靶标上圆的圆心在该相机像平面的像坐标，第二问是将已知的照相机的内部参数代入第一问建立的模型中来获得圆心在该相机像平面的像坐标值。第三问是对第一问所建立模型进行精度和稳定性进行检验。第四问也是在第一问模型的基础上进行进一步的扩展。所以只有建立第一问的模型，求出圆心在该相机像平面的像坐标，第四问模型才能建立，并且第一问的模型还要兼顾到第二问的所给的数据，否则将得再建立一个新的模型才能求解。第一问可以通过建立图像坐标系、摄像机坐标系与世界坐标系及他们之间的转换关系来建立像平面系与数码相机坐标系中坐标的对应关系的模型。但是缺乏照相机相关的内，外部的参数，第二问不能直接代入数据进行求解，但是第四问模型的基础。为此，本文建立了一个新的数学模型，他建立的基本思想是仿照质心的求法，对图像进行分析，利用像素和其坐标信息进行建模，这样只要将第二问的数据代入就能得到只要圆心在该相机像平面的像坐标。对第一问建立的模型的验证，本文是通过投影前后，在同一条直线A，B,C三点的交比不变以及他们是否仍然在同一条直线来进行的。由第一问的模型所得到的靶标像平面坐标，通过模型2可以得到数码相机的外部参数，就可以建立靶标的数学

3

模型。

三、模型假设与符号说明 A 模型假设

1 假设两部相机参数完全相同；

2 假设在拍摄时无畸变发生以及无外界情况干扰。 B符号说明

世界坐标系，也成真实或现实坐标系统，是一般三维物体的参照坐标系 摄像机坐标系，以摄像机为中心制定的坐标系，以摄像机的光学轴 为Z轴

u v 摄像机像平面坐标系，其与摄像机坐标系的xy面平行。

椭圆的中心像坐标

四、模型建立及求解

模型一： 1 模型建立：

圆在成像后为椭圆，圆心对应的即是椭圆的中心，因此求圆的圆心在相机像平面的像坐标即求成像后的椭圆的坐标。

我们定义：npix：属于椭圆集合的所有像素点的总个数；

sumx：属于椭圆集合的所有像素点的横坐标累加值； sumy：属于椭圆集合的所有像素点的纵坐标累加值；

M00a0b0npix

M10a1b0asumx

M01a0b1bsumy

我们可以通过以下式子求得椭圆的中心像坐标(xn,yn)

(xn,yn)(M10/M00,M01/M00)

Mpqf(a,b)apbq

其中(a,b)：二值化后，属于椭圆集合的像素点的像坐标；

f(a,b)：像坐标为(a,b)的像素的灰度值。对于二值化后属于椭圆集合的像素，显然f(a,b)＝1；

Mpq：图像的pq阶矩阵；

4

xnsumx/npix由上可得：

ynsumy/npix由此我们知道，要求椭圆中心的像坐标(xn,yn)，我们只需知道npix、sumx、sumy即可。

2模型算法：

我们作以下两个假设： （1）.椭圆集合相互之间无重合或交错、遮断，且近邻背景区域灰度或彩色均匀，无强烈噪声； （2）.能够采用有效的边缘检测方法和严格的判别标准，准确地检测到处于相邻两行j行和j－1行的椭圆边缘上左右两个像素点，其像坐标分别为

、

,和(、(,表示，显然=j，=j-1；

我们知道椭圆有一些重要的几何性质

1. 椭圆集合是凸集合，对于任意

性质2.椭圆曲线L为连线光滑曲线，假设

和任意实数，连线

为L上两点，其平面坐标分别

为(，则当两点距离够近时：

x4x3(y4y3)1/Lxx3x

由以上性质我们可以知道，对于j行，有： npix=(x''jx'j1)

sumx=(x''jx'j)(x''jx'j1)/2 sumy=(x''jx'j1)yj

对于j－1行，则在j行像素基础上，累加像素数npix，像面横坐标值sumx和

纵坐标值sumy分别为： npix=npix+(x''jx'j1)

5

sumx=sumx+(x''jx'j)(x''jx'j1)/2 sumy=sumy+(x''jx'j1)yj

xnsumx/npix经过计算我们可以的出npix、sumx、sumy。然后根据式子

ynsumy/npix

我们就可以得到椭圆中心的像坐标，在经过转化即可得到：

坐标系原点取在该相机的光学中心，x-y平面平行于像平面的圆心的像坐标。3模型求解：经过读图，我们可以得到A点所在椭圆的点坐标见下 x起始x终止x坐标坐标 坐标 个数 316 331 16 312 336 25 308 339 32 306 341 36 304 343 40 302 344 43 302 345 44 299 347 49 297 350 54 296 351 56 295 353 59 294 353 60 292 354 63 291 355 65 291 354 64 290 357 68 289 358 70 288 358 71 287 359 73 287 358 72 286 360 75 287 359 73 286 360 75 284 362 79 283 362 80 283 363 81 282 363 82 282 364 83 281 364 84 281 364 84

x 280

281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309

y坐标和 3810 5334 5539 6894 7620 8930 9144 10123 10505 11210 11590 11970 11970 11970 12350 12730 12886 13110 13265 13265 13797 14023 14023 14402 14402 14630 15080 15309 15309 15309

y 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176

x坐标和 5176 8100 10352 11646 12940 13889 14234 15827 17469 18116 19116 19410 20349 20995 20640 21998 22645 22933 23579 23220 24225 23579 24225 25517 25800 26163 26445 26809 27090 27090

6

y起始坐标 181 177 177 174 171 167 167 165 164 161 160 159 159 159 158 157 156 156 155 155 153 153 153 152 152 152 149 149 149 149 y终止坐标 200 204 205 209 210 213 214 217 218 219 220 221 221 221 222 223 223 224 224 224 225 226 226 227 227 228 228 229 229 229 y坐标个数 20 28 29 36 40 47 48 53 55 59 61 63 63 63 65 67 68 69 70 70 73 74 74 76 76 77 80 81 81 81

281 281 281 281 280 280 280 280 280 364 364 364 364 365 365 365 365 365 84 84 84 84 86 86 86 86 86 177 178 179 180 181 182 183 184 185 27090 27090 27090 27090 27735 27735 27735 27735 27735 149 149 149 149 149 149 147 147 147 229 229 229 230 230 231 230 230 230 81 81 81 82 82 83 84 84 84 310 311 312 313 314 315 316 317 318 15309 15309 15309 15539 15539 15770 15834 15834 15834 280 365 280 365 280 365 280 365 280 365 280 365 280 365 280 365 280 365 280 365 280 365 280 365 280 365 280 363 281 363 281 363 281 363 281 361 281 361 281 361 282 360 282 360 282 360 282 360 284 357 285 355 286 354 287 353 287 353 287 353 289 353 289 352 290 350 291 349 291

348 86 186 86 187 86 188 86 189 86 190 86 191 86 192 86 193 86 194 86 195 86 196 86 197 86 198 84 199 83 200 83 201 83 202 81 203 81 204 81 205 79 206 79 207 79 208 79 209 74 210 71 211 69 212 67 213 67 214 67 215 65 216 64 217 61 218 59 219 58 220 27735 27735 27735 27735 27735 27735 27735 27735 27735 27735 27735 27735 27735 27006 26726 26726 26726 26001 26001 26001 25359 25359 25359 25359 23717 22720 22080 21440 21440 21440 20865 20512 19520 18880 18531

7

147 230 147 230 147 230 147 230 147 230 147 230 147 230 147 230 147 230 147 229 147 229 147 229 148 229 148 229 148 228 148 228 148 228 148 227 147 227 147 227 147 227 150 226 150 226 150 224 151 224 152 224 152 223 153 222 153 222 154 221 155 221 155 220 156 218 156 217 157 216 84 319 84 320 84 321 84 322 84 323 84 324 84 325 84 326 84 327 83 328 83 329 83 330 82 331 82 332 81 333 81 334 81 335 80 336 81 337 81 338 81 339 77 340 77 341 75 342 74 343 73 344 72 345 70 346 70 347 68 348 67 349 66 350 63 351 62 352 60 353 15834 15834 15834 15834 15834 15834 15834 15834 15834 15604 15604 15604 15457 15457 15228 15228 15228 15000 15147 15147 15147 14476 14476 14025 13875 13724 13500 13125 13125 12750 12596 12375 11781 11563 11190

158 215 58 354 10817 160 214 55 355 10285 161 213 53 356 9911 163 211 49 357 9163 164 210 47 358 8789 166 208 43 359 8041 166 207 42 360 7833 169 204 36 361 6714 172 202 31 362 5797 175 199 25 363 4675 179 195 17 364 3179

由此我们可以求出npix＝5754、sumx＝1852980、sumy＝1084918 sumx/npix=322.033 sumy/npix=188.55

同理，我们可以得到B，C，D，E点在像平面的像坐标。（具体过程数据见附录） 结果为： npix sumx sumy sumx/npix sumy/npix A 5754 1852980 1084918 322.033 188.550 B 5039 2126464 987353 422.001 195.942 C 4334 2768941 919680 633.888 212.201 D 3523 2049445 1768723 581.733 502.050 E 4320 1223796 2162646 283.286 500.613 即各圆心在像平面的坐标为： A（322.033，188.550） B（422.001，195.942） C（633.888，212.201） D（581.733，502.050）

E（283.286，500.613）(坐标单位为像素，原点为左上顶点) 在图中表示为：

292 294 296 297 299 300 302 305 308 313 347 347 345 343 342 341 335 336 333 327 56 54 50 47 44 42 34 32 26 15 221 222 223 224 225 226 227 228 229 230 17892 17307 16025 15040 14102 13461 10829 10256 8333 4800

8

通过转换成坐标系原点取在该相机的光学中心，x-y平面平行于像平面的圆心的像坐标为：

A (-50.255174,-54.087249,417.195767) B (-23.809206,-52.131667,417.195767) C (33.568280,-47.830370,417,195767) D (18.447857,28.849259,417.195767) E (-60.506323,28.468915,417.195767) （单位mm）

4模型检验（稳定性与精度讨论）：

对于上述模型，所求得的五个像平面上的圆心点的像坐标分别为 A 322.0334 188.5502 B 422.0012 195.9423 C 638.8881 212.2012 D 581.7329 502.0502 E 283.2861 500.6125 （单位为像素）

由上图可以得到，A、B、C三点组成的AB与BC两条线段的斜率几乎一致，可以说明，AB与BC在同一条直线上。

我们验证模型的准确性的思路是：检验在同一条直线的两条有限线段AB与AC在世界坐标系与像坐标系中的长度比值是否相同。若比值相同或者误差在一定范围内，则说明模型准确度较高，反之则说明模型在某些方面需要进一步完善。

在世界坐标系中，由已知条件可知，|AB|=30,|AC|=100,即

=3:10=0.3000。而在像平面坐标系中，通过上述求得的坐标，可知

9

|AB|＝100.2407，|AC|＝317.7362，则比例为=0.3155，可见其相对误差

为4.9128％，可见模型相对较为准确，但仍需要改进。

为了检验模型是否具有稳定型，并再次检验改模型的误差，我们在模型正确的前提下，在给定靶标的边上再取四个点（每条边取一个，使得边被分为两个具有一定比例的线段）， （10，－50）、（50，－30）、（0，－50）、（－50，40），使得所截线段比例分别为2：3、4：1、1：1和9：1。

在像坐标中，由该模型可知，其上的点与靶标上的点一一对应。则通过次此比例可求得对应的四点坐标为（0.0387，－50.3331）、（30.5442，－32.4944）、（－21.0292，28.6591）和（－59.4813，20.2133），而通过模型求得的四点坐标为（0.0361，－53.2456）、（31.5320，－34.7694）、（－19.9443，26.4431）和（－57.5445，18.9833），四点坐标的最大误差分别为6.72％、7.00％、7.73％和3.26％。可见，模型的误差在一定范围内发生变动，但不会超出某一数值，可见模型相对比较稳定。

下面我们简要分析分辨率对模型精度的影响。当分辨率降低时（例如从1024×768降低至800×600），投影点边缘的像素会变得更低，以致于使得边缘模糊化，降低了点边缘像素的可读性，使精度总体下降。

模型二

1.模型准备：图像坐标系、摄像机坐标系与世界坐标系

摄像机采集的图像以标准电视信号的形式输入计算机，经计算机中的专用数模转换板变换成数字图像。每幅数字图像在计算机内为M×N数组，M行N列的图像中的每一个元素(称为象素)的数值即是图像点的亮

(或称灰度)。如图2

所示，在图像上定义直角坐标系u，v，每一象素的坐标(u，v)分别是该象素在数组中的列数与行数。所以，(u，v)是以象素为单位的图像坐标系的坐标。由于(u,v)只表示象素位于数组中的列数与行数，并没有用物理单位表示出该象素在图像中的位置，因而，需要再建立以物理单位(本题为毫米)表示的图像坐标系。该坐标系以图像内某一点

为原点，x轴与y轴分别与u，v轴平行，如图2所

示。 (u,v)表示以象素为单位的图像坐标系的坐标，(x，y)表示以毫米为单位的图像坐标系的坐标。在xy坐标系中，原点

定义在摄像机光轴与图像平面的交

点，该点一般位于图像中心处，但由于摄像机制作的原因，也会有些偏离(本题模型中不予以考虑)，若

在u，v坐标系中的坐标为

，每一个象素在x

轴与y轴方向上的物理尺寸为dx，dy，则图像中任意一个象素在两个坐标系下的坐标有如下关系：

10

uxu0dx

yvv0dy

图2

为方便使用，我们用齐次坐标与矩阵形式将上式表示为

1dxuv01001dy0u0x （1）v0y 11

xdx0u0dxuy0dyv0dyv （2）其逆关系可写成 11100

摄像机成像几何关系如图3表示。其中O点称为摄像机光心，Xc轴和Yc轴与图

像的x轴与y轴平行，Zc轴为摄像机的光轴，它与图像平面垂直。光轴与图像平面的交点，即为图像坐标系的原点，由点O与Xc，Yc，Zc轴组成的直角坐标系称为摄像机坐标系。OO1为摄像机焦距。

由于摄像机可安放在环境中的任何位置，我们在环境中还选择一个基准坐标系来描述摄像机的位置，并用它描述环境中任何物体的位置，该坐标系称为世界坐标系．它由Xw，Yw，Zw轴组成。摄像机坐标系与世界坐标系之间的关系可以用旋转矩阵R与平移向量t来描述。因此，空间中某一点P在世界坐标系与摄像 机坐标系下的齐次坐标如果分别为(Xw，Yw，Zw，1)T与(Xc，Yc，Zc，1)T，于是存在如下关系：

11

XcYcRTZc01XwXwYwtYwM1 （3） Zw1Zw11

其中R为3×3正交单位矩阵；

r11cos*cosr21sin*cosr31sinr12sincoscossinsinr22coscossinsinsinr32cossinr13sinsincossincosr23cossinsinsincosr33coscost为三维平移向量；0＝(0,0,0)T；M1为4×4矩阵。

图3 2.模型建立

空间任何一点P在图像上的成像位置可以用针孔模型近似表示，即任何点P在图像上的投影位置p，为光心O与P点的连线OP与图像平面的交点。这种关系也称为中心射影或透视投影。由比例关系有如下关系式：

12

fXcx （4）Zc fYcy （5）Zc其中，(x，y)为p点的图象坐标；(Xc,Yc,Zc)为空间内点P在摄像机坐标系下的坐标。我们用齐次坐标与矩阵表示上述透视投影关系：

xf0Zcy010f0Xc0 0Yc0 0 （6） Zc1 01将(1)与(2)代入上式，我们得到以世界坐标系表示的P点与其投影点p的坐标(u,v)的关系：

1dxu0Zcv1001dy0u0fv00010f00 0R0 0T01 0XwtYw1Zw1

ax0u0 0R 0ayv0 0T001 00XwtYwM1M2XwMXw （7）1Zw1其中，axf/dx,ayf/dy；M为3×4矩阵； M1完全由ax,ay,u0,v0决定，由于ax,ay,u0,v0只与摄像机内部结构有关，我们称这些参数为摄像机的内部参数；M2完全由摄像机相对于世界坐标系的方位决定，称为摄像机外部参数，我们要为摄像机定标，所要做的就是确定某一摄像机的内外参数。

由（7）可见，如果已知道摄像机的内外参数，就已知道投影矩阵M，对任何空间点P，如已知它的坐标Xw=

,就可以求出它的图象点p的位

置(u,v)。

我们首先对M矩阵进行求解。将（7）改写成

uim11m12vimZci 21m22m131m32

Xwim13 m14Ywim23 m24 （10）

Zwim33 m34113

其中，标；

为空间第i个点的坐标；为第i点的图像坐

为投影矩阵M的第i行j列元素。（10）包含三个方程：

Zciuim11Xwim12Ywim13Zwim14

Zcivim21Xwim22Ywim23Zwim24

Zcim31Xwim32Ywim33Zwim34

将上式中的第一项除以第三项,第二项除以第三项分别消去两个关于

的线性方程：

后,可得如下

m11Xwim12Ywim13Zwim14(uim31Xwiuim32Ywiuim33Zwi)uim34m21Xwim22Ywim23Zwim24(vim31Xwivim32Ywivim33Zwi)vim34

上式表示，如果标定块上有n个已知点，并已知它们的空间坐标

(i=1,2,„,n)与它们的图象点坐标

(i=1,2,„,n)，则我们

有2n个关于M矩阵元素的线性方程，下面用矩阵形式写出这些方程：

14

Xw1 Yw1 Zw1 1 0 0 0 0 -u1Xw1 -u1Yw1 -u1Zw10 0 0 0 Xw1 Yw1 Zw1 1 -v1Xw1 -v1Yw1 -v1Zw1       Xwn Ywn Zwn 1 0 0 0 0 -u1Xwn -unYwn -unZwn0 0 0 0 Xwn Ywn Zwn 1 -v1Xwn -vnYwn -vnZwnm11u1m34m12v1m34m13m14m21 m22 （11）m23m24m31m32unm34m33vnm34由（7）可见，M矩阵乘以任意不为零的常数并不影响因此，在（11）中可以指定

与(u,v)的关系，

=1，从而得到关于M矩阵其他元素的2n个线性

方程，这些未知元素的个数为11个，记为11维向量m，将（11）简写成

Km＝U

其中，K为（11）左边2n×11矩阵；m为未知的11维向量；U为（11）右边的2n维向量；K，U为已知向量。当2n>11时，我们可用最小二乘法求出上述线性方程的解为

m(KTK)1KTU

m向量与

=1构成了所求解的M矩阵。由上可见，由空间6个以上已知点与它

们的图像点坐标，我们可以求出M矩阵。本题中，我们用已知的ABCDE五点再加上靶标中心点求解。

求出M矩阵后，可由(7)所表示的关系算出摄像机的全部内外参数。但是，我们所求得的M矩阵与(7)所表示的矩阵M相差一个常数因子(11)时指定

，这一点可以从解 并不影响投影关系，

中看出。虽然我们已指出，指定

15

但在分解M矩阵时必须考虑。我们将(7)中M矩阵与摄像机内外参数的关系写成

m1T m14axm34m2T m240m3T 10其中，

r1T tx0u0 0Tr2 tyayv0 0r3T tz01 0T0 1

为由(11)求得的M矩阵的第i行的前三个元素组成的行向量；为M矩阵第i行第四列元素；分别为平移向量t的三个分量。

为旋转矩阵R的第i行；

m1T m14axr1Tu0r3T axtxu0tzTTTm34m2 m24ayr2v0r3 aytyv0tz由上式可得

m3T 1 r3T tz 比较上式两边可知：

，由于是正交单位矩阵的第三行，

＝1，因

此，我们可以从求出 再由以下式子可求得：

其中×表示向量积运算符。由以上求出的参数可进一步求出以下参数：

16

综上所说，由空间6个以上已知点以及它们的图像点坐标，我们可求出M矩阵，并可按上述式子即可求出全部内外参数。

3模型求解： 由matlab编

（见附录），可得M投影矩阵为

[51.805162 -612.537515 290.648242 -89.021404

-170.317892 497.968757 60.671301 -131.014365 -0.987058 0.246764 -3.895403 1.000000 ]

在矩阵已知的情况下，可通过

XwiuiYwi ZciviMZwi11

的等价关系计算出标靶上任意一点坐标与像平面上点的对应关系。

下面两个模型用于第四问，即确定两台相机相对位置模型的建立。 模型三 神经网络

BP神经网络：BP网络采用反馈机制，采用输出信号对输入信号产生影响来实现网络的自我调整与学习。其算法的基本思想是梯度下降法，通过自适应来修改相应藏书的权值，以使得实际输出值与期望输出值的误差均方值最小。

神经网络有效的避开了非线性算法繁琐步骤的缺陷，从而可以快速有效的解决非线性映射的问题。根据透视投影原理，左右图像中同一点的像坐标与世界坐标是一一对应的，并且两者由于畸变等因素存在着非线性关系。如果用神经网络去解决这种非线性关系，那么许多不确定的因素（比如畸变因素）就可以被包含进来并得以合理的拟合与解决。得到这种隐含的映射关系后，一旦空间一点在左右摄像机中的像坐标，就可以利用训练好的网络求出其相对世界坐标系的三维坐标。

17

BP网络有以下特：

1 BP网络是一种多层网络，包括输入层、隐层和输出层； 2 层与层之间采用全互联方式，同一层神经元之间不连接； 3 权值通过学习算法进行调节；

4 神经元激发函数为S函数；

5 学习算法由正向传播和反向传播组成；

6 层与层之间的连接是单向的，信息的传播是双向的。 单隐层BP网络结构如下图所示：

网络的逼近如下图所示：

BP算法的学习过程由正向传播和反向传播组成。在正向传播过程中，输入信息从输入层经隐层逐层处理，并传向输出层，每层神经元的状态只影响下一层神经元的状态。如果有输出层不能得到期望的输出，则转至反向传播，将误差信号按连接通路反向计算，由梯度下降法调整各层神经元的权值，使误差信号减小。

（1） 正向算法：

设输出节点k与隐层节点j之间的权值为

18

，输出层作用函数为f，则输出

信号可表示为 。由逼近图可知，网络输出与理想输出之

间的误差为e(k)＝y(k)－

（2） 反向算法： 采用

，而误差性能指标函数为 E＝。

，调整各层间的权值。如上文所述，采用梯度下降法。输出层

的学习算法为：

。隐层与

与隐层之间的连接权值

yn'·(ek)?＝ekx 其中为学习效率，jkj，

jk输入层连接权值

的学习算法为：

ynEe(k) ijijij此时权值为

。为了提高训练速度并且避免

陷入最小解，需要考虑上次权值变化对本次权值变化的影响，即加入动量因子，

jk(k1)jk(k)jk(jk(k)jk(k1))

ij(k1)ij(k)ij(ij(k)ij(k1))。

对于本实例来说，网络输入为点在像平面的坐标（x,y,z），输出则为点的世界坐标（

）,两者的关系可表示为(XW,YW,ZW)f(x,y,z)，其中，f

表示神经网络精确逼近的非线性函数。利用这两点间的对应关系对神经网络进行

训练，将非线性映射关系存储于神经网络的连接权值和阈值中，从而可以得到精度较高的空间点坐标。

网络标定的大致过程为：

（1） 采集合适的点在像平面的点的像坐标作为样本输入； （2） 以点的世界坐标作为样本输出对神经网络进行训练；

（3） 在网络实现收敛后，其输出能够反映出训练样本的特征。此时，将非

标定样本输入已收敛的神经网络中，可以得到其对应的世界坐标。 采用人工神经网络进行标定，是将相机拍摄的几何图形以及两种坐标的非线性映射关系存储于神经网络的连接权值和阈值中，从而得到两台相机的相对位置关系，使得模型具有比较稳定、外延拓展型好的特点。

模型四 遗传算法

由前三问的模型的求解，我们可以得到足够多的点的世界坐标和相应的像坐标，由此就能求出照相机的内部和外部的参数，根据R=R1R2-1和t=t1-R2-1t2就能求

19

出两个数码相机之间的几何关系。因此，由靶标给出两部固定相机相对位置问题就转化为求数码相机的外部参数。数码相机参数的求解是一个多参数，非线性优化的问题。用一般的矩阵来求解会很复杂，而遗传算法在这方面具有独特的优势，并且可以用matlab方便，快速地求出解。本文就是采用遗传算法来建模。 遗传算法是借鉴生物的自然选择和遗传进化机制而开发出来的一种全局优化自适应概率搜索算法

它使用群体搜索技术，通过对当前群体施加选择、交叉、

变异等一系列遗传操作产生新一代的群体，并逐步使群体进化到包含或接近最优解状态。以下就是用遗传算法建模的分析过程和步骤： Step1：决策变量及目标函数

决策变量包括立体视觉系统数学模型中所有的摄像机内部和外部参数，总共24个。利用前三问求出的标定点三维坐标

之间的残差来建立如下的目标函数：

与由模型计算得到的三维坐标

1N'2'2'2f()[(XwiXwi)(YwiYwi)(ZwiZwi)] Ni1其中N代表标定点的数量； 代表视觉测量系统中所有的参数。

[u0l,v0l,Nxl,Nyl,fl,kl,l,l,l,Txl,Tyl,Tzl,u0r,v0r,Nxr,Nyr,fr,kr,r,r,r,Txr,Tyr,Tzr]简单记为：(1,2,24)T

因此，用遗传算法进行摄像机标定的目的，是要找到一个由所有摄像机参数构成的最优组合来定义立体视觉三维测量系统，使

的值最小。

T

Step2：编码

将函数优化问题的解空间转换成遗传算法的搜索空间的过程称为编码。 二进制编码方法具有编码、解码过程容易操作以及交叉、变异等遗传算子便于实现等优点，是遗传算法中常用的一种编码方法。将分别代表变量

的

二进制编码串连接在一起，组成一个共 204 位的二进制编码长串，它代表目标函数优化问题的染色体编码。 Step3：解码

解码是编码的逆过程，将编码所表示的数值从搜索空间转换到解空间。首先将204位长的二进制编码串分拆成24个分别表示不同变量的二进制编码串，然后把它们分别转换成相应的十进制代码。并根据下式求出实际决策变量的值。

实际决策变量的值=十进制数*搜索精度+决策变量的下界 Step4：选择

20

选择过程是利用解码后求得的各个体适应值大小,淘汰一些较差的个体而选出一些比较优良的个体。

适应度函数是遗传算法进化的指导准则，用来度量个体在优化过程中可能达到或接近于最优解的优良程度。立体视觉系统摄像机标定是一个求目标函数的全局最小值问题，因此，适应度函数

Step5：交叉与变异

选择算子在遗传算法中以个体的适应度评价为基础来对群体中的各个个体进行优胜劣汰操作。目的是为了保持基因稳定、增强全局收敛能力和计算效率。在采用回放式随机采样方式的比例选择方法中，个体被选中的概率与其适应度大小成正比。设群体的规模大小为M，第i个个体的适应度，则个体被选中的概率为：

PiFi （i1,2,M)

i由经以下转换得到：

Fi1M交叉算子在遗传算法中起着重要的作用，是产生新个体的主要方法。算法中采用了如图3所示的单点交叉方法。

变异算子相对交叉算子来说，只是产生新个体的一种辅助方法，但也不可忽视，因为它可以改善遗传算法的局部搜索能力，保持群体中个体的多样性，避免出现早熟现象。为了不破坏太多已有的较好模式，变异概率交叉操作如图4所示。

的值取得较小。

最后只要设置一下停留代数，遗传算法中群体规模的大小，交差概率和变异概率

分别调用附录中的matlab函

，就能得到一组变量值作为目标函数

与

表

的近优解，从而求出数码相机的相关的参数，如果外参数分别用示，则

表示摄像机与世界坐标系之间的相对位置，

表示摄像机与

世界坐标系的相对位置。对任意一点P，如它在世界坐标系、坐标系与坐标

21

系下的非齐次坐标分别为，则

xc1R1xwt1

xc2R2xwt2

将上式中

消去后得到

xc1R1R21xc2t1R21t2

因此，两个摄像机之间的几何关系可用以下R和t表示：

RR1R21

tt1R21t2

进而求出两部固定相机相对位置。

五 模型的评价与改进

模型一有效而简洁的求得了五个像坐标系上的圆心坐标，但由于算法是建立在有效的边缘检测方法和严格的判别标准，准确地检测到处于相邻两行的假设上面，未考虑边缘模糊情况下导致不能准确地检测相邻两行的像素，从而引起一定的误差（4.9128％），但总体上说，模型依然可行。

模型二从数学角度出发，缜密的证明并得到关于两坐标系中对应点的非线性关系，并可通过软件求解复杂矩阵运算从而得到相应点的具体坐标。但模型一和模型二有一个共同需要改进的地方。忽略相机的畸变会使得量模型均具有较大的误差。下面，我们引入了畸变系数，对两模型进行改进。 摄像头镜头畸变分为径向畸变、切向畸变以及薄棱镜畸变三种。下面我们引入畸变程度以及畸变系数两种参数对模型进行修正。 （1） 径向畸变

22

xrxd(k1r2k2r4k3r6......)

yryd(k1r2k2r4k3r6......)

其中投影值

22rxdyd，

是图像像素在像平面坐标系的坐标，根据真空模型得到规一化的

为畸变系数。

（2） 切向畸变

22xdp1(3xdyd)2p2xdyd.....

22yd2p2(xd3yd)p1xdyd.....

其中p2、p1为切向畸变系数，（3） 薄棱镜畸变

为像点的实际坐标。

xps1r2s2r2.....

yps1r2s2r2.....

其中

为薄棱镜畸变系数。

由上述可知，图像平面理想坐标点为：

xuxdxrxpxd

yuydyxryxpyxd

模型三利用神经网络有效的避开了非线性算法繁琐步骤的缺陷，从而可以快

速有效的解决非线性映射的问题。而且，神经网络的自我学习与调节能力使得其可以忽略畸变等外界因素带来的影响而得到较为精确的坐标点。但由于其收敛缓慢等缺点，故而提出动量－自适应学习率BP神经网络的权值修正算

。其采

用动量法和学习率自适应调整的策略，从而提高了学习速度并且增加了算法准确率。动量法权值调整算法为

(k1)(k)[(1mc)D(k)mcD(k1)]

其中(k)为第k步的权值，D(k)＝

为第k步时的负梯度，

，mc为动量因子。加入的动量项相当于阻尼项，它减小了学习

过程中的震荡趋势，从而改善了收敛性。另外，自适应调整学习率的权值调整算

23

法为

(k1)(k)(k)D(k)(k)2r(k1)

rsign[D(k)D(k1)]通过连续两次迭代判断其梯度方向是否相同，可以使学习率减半。

将上述两种方法结合起来，就得到动量－自适应学习率BP神经网络的权值修正算法为

(k1)(k)(k)[(1mc)D(k)mcD(k1)]

模型四利用遗传算法能方便的解决用传统优化方法很难解决的多参数。非线性的复杂函数的优化问题，但是第四问中包含24个变量，收敛速度比较慢，减少编码长度又会对结果的精度造成影响。通过改进遗传算法，使之具有变量搜寻区间的自适应调整能力，在保持染色体编码长度不变的情况下，收敛速度和编码精度的要求。

参考文献

[1]左爱秋 吴江宁等，基于几何特性的椭圆中心像坐标的快速求取，信号处理， 第十六卷第二期：108－111，2000.6

[2]马颂德 张正友，计算机视觉——计算理论与算法基础，北京：科学出版社，1998

[3]郑阿奇，MATLAB实用教程（第2版），北京：电子工业出版社，2007.8 [4]刘金琨，智能控制，北京：电子工业出版社，2005.5

[5]刘艳梅 董再励等，一种基于神经网络的视觉定位系统，机器人，第25卷第5期：461－464，2003.

[6]黄海云 威飞虎，利用遗传算法精确标定摄像机参数，上海交通大学学报， 第34卷第7期：948－951，2000.7

[7]McFall C .An Object Infrastrusture for Internet Middleware IBM on Conponent Briker J .IEEE Internet Computing,42-50,1998

[8]刘国华 包宏 李文超，用MATLAB实现遗传算法程序，计算机应用研究80－82，2001

24

附录

1 E、B、C、D四点像坐标计算数据 x起始坐x终止坐x坐标个x坐标y起始坐y终止坐y坐标个

y

标 标 数 和 标 标 数

281 293 13 466 3731 500 509 10 277 298 22 467 6325 495 513 19 273 301 29 468 8323 493 516 24 271 304 34 469 9775 490 518 29 x 245

246 247 248 y坐标和 5045 9576 12108 14616 268 305 266 307 264 308 262 311 261 312 259 313 258 314 257 315 256 315 255 315 254 316 253 317 252 317 251 318 251 318 250 319 250 319 250 319 249 319 249 319 248 319 248 320 248 320 247 321 247 321 246 321 246 321 246 321 246 321 246 321 245 321 245 321 245 321 245 321 245 321 38 470 42 471 45 472 50 473 52 474 55 475 57 476 59 477 60 478 61 479 63 480 65 481 66 482 68 483 68 484 70 485 70 486 70 487 71 488 71 489 72 490 73 491 73 492 75 493 75 494 76 495 76 496 76 497 76 498 76 499 77 500 77 501 77 502 77 503 77 504 10887 488 12033 486 12870 485 14325 484 14898 483 15730 482 16302 481 16874 480 17130 479 17385 478 17955 477 18525 476 18777 476 19346 474 19346 473 19915 473 19915 472 19915 472 20164 471 20164 471 20412 470 20732 470 20732 469 21300 469 21300 468 21546 468 21546 468 21546 468 21546 467 21546 467 21791 467 21791 467 21791 467 21791 466 21791 466 25

520 33 521 36 523 39 524 41 525 43 526 45 527 47 528 49 529 51 529 52 530 54 530 55 531 56 532 59 532 60 533 61 532 61 532 61 533 63 533 63 533 64 533 64 534 66 534 66 534 67 534 67 534 67 534 67 534 68 535 69 535 69 535 69 535 69 534 69 534 69

249 16632 250 18126 251 19656 252 20664 253 21672 254 22680 255 23688 256 24696 257 25704 258 26182 259 27189 260 27665 261 28196 262 29677 263 30150 264 30683 265 30622 266 30622 267 31626 268 31626 269 32096 270 32096 271 33099 272 33099 273 33567 274 33567 275 33567 276 33567 277 34034 278 34569 279 34569 280 34569 281 34569 282 34500 283

34500

69 69 69 69 69 69 68 68 68 245 245 245 245 245 246 246 246 246 320 320 320 319 319 319 319 319 318 76 76 76 75 75 74 74 74 73 505 506 507 508 509 510 511 512 513 21470 21470 21470 21150 21150 20905 20905 20905 20586 466 466 466 466 466 466 466 466 466 534 534 534 534 534 534 533 533 533 284 285 286 287 288 289 290 291 292 34500 34500 34500 34500 34500 34500 33966 33966 33966 247 247 247 247 247 249 249 250 251 251 252 253 254 255 256 257 259 261 262 264 267 271 278 x起始坐标 423 318 318 318 317 316 315 315 314 313 311 310 309 308 307 306 304 303 301 299 296 294 291 283 x终止坐标 427 72 72 72 71 70 67 67 65 63 61 59 57 55 53 51 48 45 41 38 33 28 21 6 x坐标个数 5 514 20340 515 20340 516 20340 517 20022 518 19705 519 18894 520 18894 521 18330 522 17766 523 17141 524 16579 525 16017 526 15455 527 14893 528 14331 529 13464 530 12645 531 11521 532 10659 533 9240 534 7854 535 5901 536 1683

y

x坐标和 156 2125 26

466 466 467 467 467 467 468 468 468 469 469 469 470 470 470 471 472 473 473 474 475 476 477 478 480 481 483 484 485 487 491 y起始坐标 197 533 532 533 532 532 532 532 532 531 531 530 529 529 528 526 526 525 524 523 522 521 520 518 517 516 515 513 511 509 505 504 y终止坐标 201 68 67 67 66 66 66 65 65 64 63 62 61 60 59 57 56 54 52 51 49 47 45 42 40

37

35

31

28

25

19

14

4320

y坐标个

数

5

293 33966 294 33433 295 33500 296 32967 297 32967 298 32967 299 32500 300 32500 301 31968 302 31500 303 30969 304 30439 305 29970 306 29441 307 28386 308 27916 309 26919 310 25922 311 25398 312 24402 313 23406 314 22410 315 20895 316 19900 317 18426 318 17430 319 15438 320 13930 321 12425 322 9424 323 6965

x

y坐标和 382 995

416 412 410 407 405 404 402 400 399 434 437 439 442 444 445 447 448 449 19 26 30 36 40 42 46 49 51 157 158 159 160 161 162 163 164 165 8075 11037 12735 15282 16980 17829 19527 20776 21624 190 186 184 181 179 177 176 175 173 208 211 214 216 217 219 220 222 223 19 26 31 36 39 43 45 48 51 383 384 385 386 387 388 389 390 391 3781 5161 6169 7146 7722 8514 8910 9528 10098 398 450 397 451 396 452 395 453 394 454 393 455 392 455 391 456 391 456 390 457 389 458 388 458 388 458 387 459 387 459 386 460 386 460 386 460 385 461 385 461 384 461 384 461 384 461 384 461 383 461 383 462 383 462 383 462 383 462 383 462 383 461 382 461 382 461 382 461 382

461 53 166 55 167 57 168 59 169 61 170 63 171 64 172 66 173 66 174 68 175 70 176 71 177 71 178 73 179 73 180 75 181 75 182 75 183 77 184 77 185 78 186 78 187 78 188 78 189 79 190 80 191 80 192 80 193 80 194 80 195 79 196 80 197 80 198 80 199 80 200 22472 172 23320 171 24168 170 25016 169 25864 168 26712 167 27104 166 27951 165 27951 164 28798 164 29645 163 30033 163 30033 162 30879 161 30879 161 31725 160 31725 160 31725 160 32571 159 32571 159 32955 158 32955 158 32955 158 32955 158 33338 157 33800 157 33800 157 33800 157 33800 157 33800 157 33338 157 33720 156 33720 156 33720 156 33720

156 27

224 53 225 55 226 57 227 59 228 61 228 62 229 64 230 66 231 68 231 68 232 70 232 70 232 71 233 73 233 73 234 75 234 75 234 75 235 77 235 77 235 78 235 78 235 78 235 78 235 79 235 79 235 79 235 79 235 79 235 79 235 79 235 80 235 80 235 80 235 80 392 10494 393 10890 394 11286 395 11682 396 12078 397 12245 398 12640 399 13035 400 13430 401 13430 402 13825 403 13825 404 13987 405 14381 406 14381 407 14775 408 14775 409 14775 410 15169 411 15169 412 15327 413 15327 414 15327 415 15327 416 15484 417 15484 418 15484 419 15484 420 15484 421 15484 422 15484 423 15640 424 15640 425 15640 426 15640

382 383 383 383 383 383 383 383 384 461 461 460 460 460 460 459 459 458 80 79 78 78 78 78 77 77 75 201 202 203 204 205 206 207 208 209 33720 33338 32877 32877 32877 32877 32417 32417 31575 156 235 80

157 235 79 157 234 78 157 234 78 157 234 78 157 234 78 157 233 77 157 233 77 158 232 75 427

428 429 430 431 432 433 434 435 15640 15484 15249 15249 15249 15249 15015 15015 14625 384 384 385 385 385 386 386 387 388 388 389 390 390 391 392 393 394 395 396 398 399 400 402 405 407 410 x起始坐标 637 633 630 628 626 624 458 458 457 457 456 455 455 454 454 453 452 451 450 449 448 447 446 445 443 442 440 439 437 434 432 427 x终止坐标 651 654 657 659 660 662 75 75 73 73 72 70 70 68 67 66 64 62 61 59 57 55 53 51 48 45 42 40 36 30 26 18

x坐标个数

15 22 28 32 35 39

210 31575 211 31575 212 30733 213 30733 214 30276 215 29435 216 29435 217 28594 218 28207 219 27753 220 26912 221 26071 222 25620 223 24780 224 23940 225 23100 226 22260 227 21420 228 20136 229 18900 230 17619 231 16780 232 15102 233 12585 234 10907 235 7533

y x坐标和 175 9660 176 14157 177 18018 178 20592 179 22505 180

25077

28

158 158 159 159 160 160 160 161 161 162 163 163 164 165 166 167 168 169 170 171 173 175 176 179 181 184 191 y起始坐标 210 206 203 201 199 197 232 232 231 231 230 229 229 228 227 227 226 225 224 223 222 221 220 219 217 216 214 213 211 208 206 202 195 y终止坐标 224 228 231 233 235 236 75 75 73 73 71 70 70 68 67 66 64 63 61 59 57 55 53 51 48 46 42 39 36 30 26 19 5 y坐标个数

15 23 29 33 37 40

436 14625 437 14625 438 14235 439 14235 440 13845 441 13615 442 13615 443 13226 444 12998 445 12837 446 12448 447 12222 448 11834 449 11446 450 11058 451 10670 452 10282 453 9894 454 9288 455 8901 456 8127 457 7566 458 6966 459 5805 460 5031 461 3667 462 965 x y坐标和 603 3255 604 4991 605 6293 606 7161 607 8029 608 8660

622 621 620 618 617 616 615 614 613 663 664 666 666 667 668 669 670 670 42 44 47 49 51 53 55 57 58 181 182 183 184 185 186 187 188 189 26985 28270 30221 31458 32742 34026 35310 36594 37207 195 193 192 191 190 188 187 186 185 238 239 240 241 242 243 243 244 245 44 47 49 51 53 56 57 59 61 609 610 611 612 613 614 615 616 617 9526 10152 10584 11016 11448 12068 12255 12685 13115 613 671 612 671 611 672 610 672 610 673 609 673 609 673 608 674 608 674 607 674 607 674 606 674 606 675 605 675 605 675 605 675 604 675 604 675 604 675 604 675 603 675 603 675 603 675 603 674 603 674 603 674 603 674 603 674 603 674 603 673 603 673 603 673 603 672 603 672 603

671 59 190 60 191 62 192 63 193 64 194 65 195 65 196 67 197 67 198 68 199 68 200 69 201 70 202 71 203 71 204 71 205 72 206 72 207 72 208 72 209 73 210 73 211 73 212 72 213 72 214 72 215 72 216 72 217 72 218 71 219 71 220 71 221 70 222 70 223 69 224 37878 184 38490 184 39773 183 40383 182 41056 181 41665 181 41665 180 42947 180 42947 179 43554 179 43554 178 44160 178 44835 177 45440 177 45440 177 45440 176 46044 176 46044 176 46044 176 46044 175 46647 175 46647 175 46647 175 45972 175 45972 175 45972 175 45972 175 45972 175 45972 175 45298 175 45298 175 45298 175 44625 175 44625 175 43953

176 29

246 63 246 63 247 65 247 66 247 67 248 68 248 69 248 69 249 71 249 71 249 72 249 72 249 73 249 73 249 73 250 75 250 75 250 75 249 74 249 75 249 75 249 75 249 75 249 75 249 75 248 74 248 74 248 74 248 74 247 73 247 73 246 72 246 72 245 71 245 70 618 13545 619 13545 620 13975 621 14157 622 14338 623 14586 624 14766 625 14766 626 15194 627 15194 628 15372 629 15372 630 15549 631 15549 632 15549 633 15975 634 15975 635 15975 636 15725 637 15900 638 15900 639 15900 640 15900 641 15900 642 15900 643 15651 644 15651 645 15651 646 15651 647 15403 648 15403 649 15156 650 15156 651 14910 652 14735

604 671 68 604 670 67 604 670 67 604 669 66 605 669 65 605 668 64 605 668 64 606 667 62 606 666 61 225

226 227 228 229 230 231 232 233 43350 42679 42679 42009 41405 40736 40736 39463 38796 176 176 177 177 177 178 178 179 180 244 244 243 243 242 241 240 240 239 69 69 67 67 66 64 63 62 60 653 654 655 656 657 658 659 660 661 14490 14490 14070 14070 13827 13408 13167 12989 12570 607 607 608 609 609 610 611 612 613 614 616 617 618 620 623 627 633 x起始坐标 580 576 574 571 569 568 566 565 563 562 561 560 559 558 557 556 665 665 664 663 662 661 660 658 657 656 654 653 650 648 646 642 635 x终止坐标 596 599 601 603 605 606 607 608 609 610 611 611 612 613 613 614 59 59 57 55 54 52 50 47 45 43 39 37 33 29 24 16 3 x坐标个数

17 24 28 33 37 39 42 44 47 49 51 52 54 56 57 59

234 37524 235 37524 236 36252 237 34980 238 34317 239 33046 240 31775 241 29845 242 28575 243 27305 244 24765 245 23495 246 20922 247 18386 248 15228 249 10152 250 1902 y x坐标和 470 9996 471 14100 472 16450 473 19371 474 21719 475 22893 476 24633 477 25806 478 27542 479 28714 480 29886 481 30446 482 31617 483 32788 484 33345 485 34515

30

180 181 182 183 183 185 186 187 188 190 192 194 197 202

y起始坐标

501 498 495 493 491 489 488 486 485 484 483 482 481 480 479 478

238 237 236 235 233 232 231 229 227 225 223 221 218 213

y终止坐标

514 518 520 522 524 525 526 527 528 529 529 530 531 531 532 532

59 57 55 53 51 48 46 43 40 36 32 28 22 12

y坐标个数

14 21 26 30 34 37 39 42 44 46 47 49 51 52 54 55

662 12331 663 11913 664 11495 665 11077 666 10608 667 10008 668 9591 669 8944 670 8300 671 7470 672 6640 673 5810 674 4565 675 2490

x y坐标和 548 7105 549 10668 550 13195 551 15225 552 17255 553 18759 554 19773 555 21273 556 22286 557 23299 558 23782 559 24794 560 25806 561 26286 562 27297 563

27775

555 555 554 553 553 552 552 551 551 614 614 615 615 615 615 615 616 616 60 60 62 63 63 64 64 66 66 486 487 488 489 490 491 492 493 494 35070 35070 36239 36792 36792 37344 37344 38511 38511 478 477 476 476 475 474 474 473 473 532 533 533 533 534 534 535 534 534 55 57 58 58 60 61 62 62 62 564 565 566 567 568 569 570 571 572 27775 28785 29261 29261 30270 30744 31279 31217 31217 550 616 550 616 550 616 549 616 549 616 549 615 548 615 548 615 548 615 548 615 548 615 548 614 548 614 548 614 548 613 548 613 548 612 548 612 548 611 548 611 549 610 549 609 549 609 549 608 550 607 550 606 551 606 551 604 552 604 553 602 553 601 554 600 555 599 556 597 558

596 67 495 67 496 67 497 68 498 68 499 67 500 68 501 68 502 68 503 68 504 68 505 67 506 67 507 67 508 66 509 66 510 65 511 65 512 64 513 64 514 62 515 61 516 61 517 60 518 58 519 57 520 56 521 54 522 53 523 50 524 49 525 47 526 45 527 42 528 39 529 39061 473 39061 472 39061 472 39610 472 39610 472 38994 472 39542 472 39542 470 39542 470 39542 470 39542 470 38927 470 38927 470 38927 470 38313 470 38313 470 37700 470 37700 470 37088 470 37088 470 35929 470 35319 470 35319 470 34710 470 33553 471 32946 471 32396 471 31185 472 30634 472 28875 473 28273 473 27119 474 25965 474 24213 475 22503

476 31

534 62 534 63 534 63 534 63 534 63 534 63 535 64 534 65 534 65 534 65 534 65 533 64 533 64 533 64 533 64 532 63 532 63 532 63 531 62 531 62 530 61 530 61 529 60 529 60 528 58 527 57 527 57 526 55 525 54 524 52 523 51 523 50 521 48 521 47 519 44 573 31217 574 31689 575 31689 576 31689 577 31689 578 31689 579 32224 580 32630 581 32630 582 32630 583 32630 584 32096 585 32096 586 32096 587 32096 588 31563 589 31563 590 31563 591 31031 592 31031 593 30500 594 30500 595 29970 596 29970 597 28971 598 28443 599 28443 600 27445 601 26919 602 25922 603 25398 604 24925 605 23880 606 23406 607 21890

559 560 562 565 568 594 592 590 587 583 36 33 29 23 16 530 531 532 533 534 20754 19008 16704 13248 9208

2 求解M矩阵matlab程序

k=[-50 50 50 1 0 0 0 0 -16101.67 -16101.67 16101.67; 0 0 0 0 -50 50 50 1 -9427.51 -9427.51 9427.51; -20 -50 50 1 0 0 0 0 -8440.024 -21100.06 21100.06; 0 0 0 0 -20 -50 50 1 -3918.846 -9797.115 9797.115; 50 -50 50 1 0 0 0 0 31944.405 -31944.405 31944.405; 0 0 0 0 50 -50 50 1 10610.06 -10610.06 10610.06; 50 50 50 1 0 0 0 0 29086.645 29086.645 29086.645; 0 0 0 0 50 50 50 1 25102.51 25102.51 25102.51;

-50 50 50 1 0 0 0 0 -14164.305 14164.305 14164.305; 0 0 0 0 -50 50 50 1 -25030.625 25030.625 25030.625; 0 0 50 1 0 0 0 0 0 0 22824.25625; 0 0 0 0 0 0 50 1 0 0 17542.67625;] u=[322.0334 188.5502 422.0012 195.9423 638.8881 212.2012 581.7329 502.0502 283.2861 500.6125

456.485125 350.853525] m=inv(k'*k)*k'*u

3 验证三点共线matlab程序

x=[322.0334 422.0012 638.8881] y=[188.5502 195.9423 212.2012] scatter(x,y) hold on plot(x,y)

4神经网络部分源程序

void NN::addNerveCellByNum(U2 n) /*n为增加多少个*/

32

477 478 479 480 482 483 485 488 493 518 517 515 514 512 510 508 505 499 42 40 37 35 31 28 24 18 7 608 609 610 611 612 613 614 615 616 20895 19900 18389 17395 15407 13902 11916 8937 3472

{

for(int i=1;i<=n;i++) {

*nc->next = (LOGICNERVECELL *)(malloc(sizeof(LOGICNERVECELL))); *nc->next->nc = (NERVECELL*)(malloc(sizeof(NERVECELL))); if(i==n) {

*nc->next->nc->nerve_cell_num = n_cell; *nc->next->nc->input = NULL; *nc->next->nc->output = NULL; *nc->next = NULL; } else

{

*nc->next->nc->nerve_cell_num = i+n_cell; *nc->next->nc->input = NULL; *nc->next->nc->output = NULL; *nc = *nc->next; } }

*nc = frist_nc; /*重新指想第一个神经元*/ n_cell = n_cell+n; return; }

void NN::delNerveCell(U2 n) /*n为神经元号*/ { if(n<=0) {

printf(\"invaild nerve cell num\\n\"); exit(ERROR_CODE); } else {

LOGICNERVECELL *tmpnc; tmpnc = findNerveCell(n); recycle[recycle_cell] = *nc; tmpnc->next = *nc->next; *nc = frist_nc; /*指向头*/ }

recycle_cell++;

return; }

void NN::connectNerveCell(U2 sou,U2 dec)

33

{

U2 init;

LOGICNERVECELL *sou_nc,dec_nc,tmpnc; tmpnc = findNerveCell(sou); sou_nc = *nc;

tmpnc = findNerveCell(dec); dec_nc = *nc;

sou_nc->nc->output = (NERVECELL*)(malloc(sizeof(NERVECELL)*(sou_nc->nc->o_n_cell))); sou_nc->nc->output[tick] = dec_nc; sou_nc->nc->o_is_cell = tick; sou_nc->nc->tick++;

dec_nc->nc->input = (NERVECELL*)(malloc(sizeof(NERVECELL)*(dec_nc->nc->i_n_cell))); dec_nc->nc->input[tick] = sou_nc; dec_nc->nc->i_is_cell = tick; dec_nc->nc->tick++; *nc = frist_nc;

5 遗传算法函数程序

function [bin_gen ,bits] =encoding(min_var ,max_var ,scale_var ,popsize) bits=ceil(log2((max_var-min_var)./scale_var) ) ; bin_gen=randint(popsize,sum(bits));

function new_gen = crossover (old_gen ,pc)

[ nouse ,mating] = sort (rand(size(old_gen ,1),1)) ; Mat_gen = old_gen(mating , :) ; pairs = size(mat_gen ,1)/2; bits = size(mat_gen ,2) ; cpairs = rand(pairs ,1)cpoints = randint (pairs ,1 ,[1 ,bits]) ; cpoints = cpairs.*cpoints; for i = 1:pairs

new_gen( [2*i-1 2*i ],:) = [mat_gen ([2*i-1 2*i ],1:cpoints(i)) mat_gen([2*i 2*i-1],cpoints(i)+1:bits)] ; end

function[var_gen ,fitness] = decoding (funname,bin_gen,bits,min_var ,max_var) num_var = length(bits); popsize = size (bin_gen,1) ;

scale_dec = (max_var-min_var)./(2.^bits-1) ; bits = cumsum(bits) ; bits = [0 bits];

34

for i = 1:num_var

bin_var{i} = bin_gen( : ,bits(i)+1:bits(i+1)) ; var{i}=sum(ones(popsize,1)*2.^(size (bin_var{i} ,2)-1 bin_var{i},2).*scale_dec(i)+min_var (i) ; end

var_gen = [var{1 , :} ] ; for i = 1 :popsize

fitness (i) = eval ( [funname ,'(var_gen(i , :) ) ']) ; end

function new_gen = mutation (old_gen ,pm) mpoints = find(rand(size(old_gen))< pm) ; new_gen = old_gen ;

new_gen(mpoints) = 1-old_gen(mpoints) ;

function [ evo_gen ,best_indiv ,max_fitness ] = selection ( old_gen ,fitness) popsize = length(fitness) ;

[max_fitness ,index1 ] = max(fitness) ; [min_fitness , index2 ] =min(fitness) ; Best_indiv = old_gen(index1 , :) ; index =[1:popsize ] ; index(index1) = 0 ; index(index2) = 0 ;

index = nonzeros(index) ; evo_gen = old_gen(index , :) ; evo_fitness = fitness(index , :) ; evo_popsize = popsize-2 ;

ps = evo_fitness/sum(evo_fitness) ; pscum= cumsum(ps) ; r = rand(1,evo_popsize) ;

selected = sum( pscum *ones (1 ,evo_popsize):-1 :0).*

35