第３４卷第２期 测绘与空间地理信息 Ｖｏ１．３４，Ｎｏ．２ ２０１１年４月 ＧＥＯＭＡＴＩＣＳ＆ｓＰＡＴｌＡＬ ｌＮＦｏＲＭＡ，ｌ１ＩｏＮ ＴＥＣＨＮＯＬＯＧＹ Ａｐｒ．，２０１１ （辽宁省基础测绘院，辽宁锦州１２１００３） 摘要：针对影像匹配问题，对最小二乘法影像匹配的算法进行了研究，理解了该算法的匹配原理，并应用测试 程序对该算法进行测试。 关键词：最小二乘法；影像灰度重采样；影像匹配 中图分类号：Ｐ２３１．２ 文献标识码：Ｂ 文章编号：１６７２—５８６７（２０１１）０２—０１９０—０３ Ｄｅｓｉｇｎ ａｎｄ Ｉｍｐｌｅｍｅｎｔａｔｉｏｎ ｏｆ ｔｈｅ Ａｌｇｏｒｉｔｈｍ ｏｆ Ｌｅａｓｔ Ｓｑｕａｒｅｓ Ｉｍａｇｅ Ｍａｔｃｈｉｎｇ Ｕ Ｆｕ—ｊｉｎ，ＱＩＮ Ｚｈｉ—ｗｅｉ （Ｌｉａｏｎｉｎｇ Ｐｒｏｖｉｎｃｉａｌ Ｂａｓｉｃ Ｓｕｒｖｅｙｉｎｇ ａｎｄ Ｍａｐｐｉｎｇ Ｉｎｓｔｉｔｕｔｅ，Ｊｉｎｚｈｏｕ １２１００３，Ｃｈｉｎａ） Ａｂｓｔｒａｃｔ：Ａｉｍｉｎｇ ａｔ ｔｈｅ ｐｒｏｂｌｅｍ ｏｆ ｉｍａｇｅ ｍａｔｃｈｉｎｇ，ｔｈｉｓ ｐａｐｅｒ ｓｔｕｄｉｅｄ ｔｈｅ ａｌｇｏｒｉｈｔｍ ｏｆ ｌｅａｓｔ ｓｑｕａｒｅｓ ｉｍａｇｅ ｍａｔｃｈｉｎｇ ａｎｄ ｄｉｓｃｕｓｓｅｄ ｔｈｅ ｍａｔｃｈｉｎｇ ｐｒｉｎｃｉｐｌｅ ｏｆ ｔｈｉｓ ａｌｇｏｒｉｔｈｍ．Ｉｔ ａｌｓｏ ｔｅｓｔｅｄ ｔｈｅ ａｌｇｏｒｉｔｈｍ ｕｓｉｎｇ ｔｅｓｔｉｎｇ ｐｒｏｇｒａｍ． Ｋｅｙ ｗｏｒｄｓ：ｌｅａｓｔ ｓｑｕａｒｅｓ；ｉｍａｇｅ ｇｒａｙ ｌｅｖｅｌ ｒｅｓａｍｐｌｉｎｇ；ｉｍａｇｅ ｍａｔｃｈｉｎｇ ０ 引 言 ｍａｔｃｈｉｎｇ）或“多片”影像匹配（ｍｕｌｔｉ—ｐｈｏｔｏ ｍａｔｃｈｉｎｇ）。另 外，在最小二乘影像匹配系统中，可以很方便地引入“粗 随着科学技术的迅猛发展，影像匹配已成为现代信 差检测”，从而大大地提高了影像匹配的可靠性。它甚至 息处理，特别是图像信息处理领域中的一项非常重要的 还可以用于解决影像遮蔽问题（ｏｃｃｌｕｓｉｏｎ）。 技术，并在计算机视觉、数字摄影测量、遥感图像处理、资 源分析、医学图像配准、肖像和指纹识别、光学和雷达跟 １最小二乘法求法及影像灰度重采样 踪、飞行器巡航制导、导弹地形匹配及投射系统的导航目 １．１ 最小二乘法的求法 标制导等领域得到了广泛的应用。 最小二乘法是处理各种观测数据进行测量平差的一 最小二乘法在影像匹配中的应用是２Ｏ世纪８０年代 种基本方法。如果以不同精度多次观测一个或多个未知 发展起来的。德国Ａｃｋｅｒｍａｎｎ教授提出了一种新的影像 量，为了求定各未知量的最可靠值，各观测量必须加改正 匹配方法——最小二乘影像匹配（１ｅａｓｔ ｓｑｕａｒｅｓ ｉｍａｇｅ 数，使其改正数的平方乘以观测值的权数的总和为最小。 ｍａｔｃｈｉｎｇ）。由于该方法充分利用了影像窗口内的信息进 因此称作最小二乘法。所谓“权”就是表示观测结果质量 行平差计算，使影像匹配可以达到１／１０甚至１／１００像素 相对可靠程度的一种权衡值。下面以多项式拟合为例简 的高精度，即影像匹配精度可达到子像素（ｓｕｂｐｉｘｅ１）等级。 要介绍最小二乘法的求法。 为此，最小二乘影像匹配被称为“高精度影像匹配”，但也 设已知点（ ，Ｙ ）（ｉ＝１，２，…ｎ），求ｍ次多项式Ｐ 有人习惯于称其为“高精度影像相关”。它不仅仅可以被 （ ）来拟合函数，（ ）。设Ｐ （ ）＝ａ０＋ａｌ ＋…＋ＵｍＸ ，求 用于一般的产生数字地面模型，生产正射影像图，而且可 系数口０，０ｌ，…，０ ，使Ｆ（Ⅱ０，Ⅱ１，…，口 ）＝乏 ＝ｌ∑ 以用于控制点的加密（空中三角测量）及工业上的高精度 【Ｙ 一Ｐ （ ）］ 达到最，Ｊ、。 量测。由于在最小二乘影像匹配中可以非常灵活地引入 当拟合函数是一元函数时，所对应的函数图形是平 各种已知参数和条件（如共线方程等几何条件、已知的控 面曲线。这时，数据拟合问题的几何背景是寻求一条近 制点坐标等），从而可以进行整体平差。它不仅可以解决 似通过给定离散点的曲线，故称为曲线拟合问题。 “单点”的影像匹配问题，以求其“视差”；也可以直接解求 由于函数Ｆ（ａ。，ａ 一ａ ）达到最小，所以由高等数学 其空间坐标；且可以同时解球待定点的坐标与影像的方 知识有： 位元素；还可以同时解决“多点”影像匹配（ｍｕｌｔｉ—ｐｏｉｎｔ 收稿日期：２０１０—１２—２２ 作者简介：李福金（１９５５一），男，辽宁锦州人，高级工程师，主要从事遥感与ＧＩＳ的应用与研究工作。 第２期 李福金等：最小二乘影像匹配算法设计与实现 ｌ９１ 芸＝２砉 —Ｐｍ（ ］ ＝２砉 一笔。 ］ ＝。 左影像坐标间的关系为： （１） 其中ｎ。（ ，Ｙ ），／７，：（ ，Ｙ２）为偶然误差，而右影像坐标与 ２＝Ⅱ０＋口ｉｘｌ＋０２Ｙｌ Ｙ２ ｂ０＋ｂｌ １＋ｂ２ｙｌ ｆ６） 即：∑ｏ （Ｅ ｋ ）＝∑ｙ （ ＝０，１，２，…，ｍ）（２） 于是得线性方程组： 口１∑ ｌ＋…＋０ ∑ ＝∑Ｙ ｉ＝ｌ ‘＝ｌ ｆ＋ｏｌ∑ ＋… ＋％∑ ＝∑Ｙｌ ＋口ｌ∑ “＋ 一ｉ：Ｉ ＋ｏ ∑ ２； ＝∑ｙ （３） 以上方程组称为法方程组。 可以证明，该方程组有唯 一解： ａｏ＝ （ｋ＝０，１，…，ｍ） （４） 并且相应的函数Ｐ （ ）＝ ＋口。 ＋…＋０ 就是满足方 程组的最小二乘解。 １．２影像灰度重采样 影像重采样在数字影像处理中是不可缺少的过程。 当欲知不位于矩阵（采样）点上的原始函数ｇ（ｘ，Ｙ）的数值 时就需要进行内插，此时称为重采样，意即在原采样基础 上再一次采样。 常用的灰度重采样方法有３种： １）最近邻法 这种方法是距离投影点最近元灰度值代替输出像元 灰度值。此法优点是简单、省时，但精度较差，图斑边缘 会出现“锯齿”现象。 ２）双线性插值法 此法考虑投影点周围４个相邻像元灰度值，并根据各 自权重计算输出像元灰度值。 ３）双三次卷积法 该方法是取与投影点邻近的１６个像元灰度值计算输 出像元灰度值。 ２最小二乘影像匹配算法设计与实现 ２．１ 最小二乘影像算法的设计 首先在左影像中选择一个适当大小的区域，作为匹 配的目标区。目标区的大小按照实际情况而定，不宜过 大也不宜过小。过大就会使匹配准确度降低，也会增加 运算量，使匹配效率降低；过小会导致参与平差计算的像 素点数目不足，从而会影响最小二乘法计算出的参数不 能满足要求。因此，合适的目标区域是非常重要的。 选择完目标区后即确定了以左影像为基准，利用变 换参数对右影像进行几何变换和辐射变换，变换模型 如下： ｇ１（ １，Ｙ１）＋ｎ１（ ｌ，Ｙ１）＝ｈｏ＋ ｌｇ２（ ２，Ｙ２）＋ （ １，Ｙ１） （５） 以上变换参数的初值为ｈ０：０，ｈＩ＝１，‰＝０，ａ。＝１，口２＝ ０，ｂ０＝０，ｂｌ＝０，ｂ２＝０。 最ｄｘ－－乘影像匹配算法可分为以下几个步骤来 理解： 第一步：先将目标区各个像素点的坐标按（６）式进行 变换，如图１所示。 一一 （ａ）左影像 （ｂ）右影像 图１选取的目标区域与经参数变换后的影像对比图 Ｆｉｇ．１ Ｔｈｅ ｃｏｍｐａｒｉｓｏｎ ｂｅｔｗｅｅｎ ｔｈｅ ｓｅｌｅｃｔｅｄ ｏｂｊｅｃｔ ａｒｅ８ ａｎｄ ｔｈｅ ｉｍａｇｅ ｏｆｆｅｒ ｐａｒａｍｅｔｒｉｃ ｔｒａｎｓｆｏｒｍａｔｉｏｎ 目标区经过变换后，在右影像中的形状和位置的情 况（上图中的变换参数取的是初始值，因此右影像中的区 域形状与目标区相比，没有发生变形）。 第二步：完成几何畸变改正后，图１中Ｂ区域内每个 像素点的坐标值一般都不会为整数，因此需要灰度重采 样（第一次几何畸变改正后，只需读取该区域中的各像素 点在右影像中的灰度值即可），可以采用间接法重采样， 即对Ｂ区域中的每个像素点用双线性内插的方法计算其 灰度值，将计算后得到的灰度值赋给Ｂ区域中的对应像 素点，像素点的坐标值与目标区域中像素点的坐标值是 对应的。 第三步：对于重采样后的像素点用（５）式进行辐射畸 变的改正，按照假设此时Ｂ区域的各像素点的灰度值估 值应该与目标区域对应像素点的灰度值估值相等。可以 用下式表示： ｇ１（ １，Ｙ１）：ｈｏ＋＾】ｇ２（ １，Ｙ１） （７） 第四步：按照最小二乘影像匹配算法计算出畸变参 数的改正值，从而计算出畸变参数。接着计算目标区域 与Ｂ区域的相关系数，由于匹配过程是一个迭代过程，因 此此次计算的相关系数若大于上一次的相关系数，就可 以认为匹配结束。否则用计算得出的畸变参数重新进行 计算，得到新的Ｂ区域灰度值。再次计算相关系数并判 断迭代是否结束。 以上过程可以认为Ｂ区域中的小区域影像是每次用 畸变参数计算后的更加接近匹配区域的小区域影像代替 １９２ 测绘与空间地理信息 ２０１１卑 后得到的。随着畸变参数的不断变化，Ｂ区域中的小区域 影像就无限接近于匹配区域的影像，通过计算最佳匹配 点坐标，就可以得出目标区与右影像的最佳同名点的 坐标。 配，最终求出目标区与匹配区的最佳匹配点的坐标值。 ３实验分析 为了对匹配算法进行测试，选取三对同名点进行测 ２．２最小二乘影像算法的实现 最小二乘法影像匹配虽然匹配精度高，但是需要目 标区每个像元参加计算，运算时间过长，加之需要多次迭 试，测试程序共分为６个主要部分，即粗匹配部分、几何畸 变处理部分、灰度重采样部分、辐射畸变处理部分、最小 二乘法计算参数改正数部分和计算相关系数部分等。 从实验的效果来看，目标区域选择３０ Ｘ ３０的区域效 代后才可以得出最佳参数值。因此，可以先在右影像中 人工选出匹配区域，然后用目标区与匹配区域中的像素 果会比较好，选择区域较大或较小时效果都不是很好。 进行最小二乘匹配，得出精确匹配点的坐标值。本文中 所编写的实现算法的程序采用了这种思路，即首先在左、 ４结束语 右影像中选出几个同名区域，然后在每个区域中选出一 本文针对影像匹配问题，对最小二乘法影像匹配的 个最有意义的点，选择的条件可以是对每个区域的像素 算法进行了深入研究，理解了该算法的匹配原理，并应用 点坐标求加权平均值，权值可以取每个像素点的灰度值 测试程序对该算法进行测试。最小二乘影像匹配算法匹 除以区域中所有像素点灰度值的和。可用下式简要 配精度高，能够很方便地引入变形参数进行整体平差，从 ＿ ：—茛 ：．Ｊ一一 而抵消左、右影像之间的几何差异和辐射差异，既考虑了 ＸｌｌＸＰｌ＋ Ｉ２， ２丰…＋ ｌ＾Ａ ｎ 辐射误差对影像匹配的影响，同时也考虑了几何畸变对 ＝———————●—————————————一 印１＋印２＋…＋印ｎ ， 影像匹配的影响。避免了只考虑灰度的匹配方法出现的 ｕ， ＹＩ１ＹＰｌ＋Ｙｌ２ＹＰ２＋…＋ｙｌ ，，ｐ 误匹配，如目标点与最终计算得出的匹配点虽然灰度相 Ｖ＝————————————————————＿————————一 卯ｌ：＋ ２＋…＋卯 似性很大，但它们可能并不是同名点，因为它们的几何位 （ ，ｙ）为最有意义的点的坐标，Ｐ。，Ｐ ，…，Ｐ 为点（ 置可能存在很大的差异，而最小二乘影像匹配的方法由 Ｙ ），（ Ｉ２，Ｙ ２），…，（ Ｙ】 ）韵权值，即ｔ 于考虑了目标区与匹配区的几何畸变，所以能有效地避 ｐｌ＝ 免此类问题的产生。 ｇｌ（ ｌｉ，Ｙｌｉ） ｇｌ（ ｌ１， １）＋ （ ｌ２，， ）＋…＋ｇｌ（ １ｉ，Ｙｌ１）＋一・＋ （ ｌｎ，ｙ１ ） 参考文献： （９） ［１］张祖勋，张剑清．数字摄影测量学［Ｍ］．武汉：武汉大 用类似的方法可以求出右影像中所选区域内最有意义的 学出版社，２００２． 点的坐标。这样就可以粗略求出左、右影像同名点的坐 ［２］李德仁．摄影测量与遥感的现状与发展趋势［Ｊ］．武汉 标平移值 ，ｄｙ可以用下式表示： 测绘科技大学学报，２０００，２５（１）：１—５． Ｘ２＝ Ｉ＋ｄｘ ［３］ 王岩松，阮秋琦．基于最小二乘互相关算法的图像定位 （１Ｏ） Ｙ２＝Ｙ１＋ 匹配研究［Ｊ］．通信学报，２００３，２４（９）：８０—８６． 此处只考虑左、右影像的平移变换。 ［责任编辑：王丽欣］ 接着就可以用最小二乘影像匹配算法进行精确匹 （上接第１８９页） ２．３碎部测量数据采集软件的实现体系结构 一个动态链接库，可自动地加载到ＡｕｔｏＣＡＤ绘图环境中， 在ＡｕｔｏＣＡＤ中加载ＯｂｊｅｃｔＡＲＸ程序有多种方法，如 使其成为一个大比例尺数据测图系统，该系统实现了外 修改注册表法，命令加载法，自动加载法。最方便的方法 业控制、碎部测量和数据成图一体化作业。 是修改ＡｕｔｏＣＡＤ安装目录中的ａｃａｄ．ＩＸ的文本文件，在其 中加入碎部测量数据采集程序模块ｓｕｒｖｅｙ．ａｒｘ，这样Ａｕｔｏ． 参考文献： ＣＡＤ在运行时就会自动加载ｓｕｒｖｅｙ．ａｒｘ模块，该模块自动 ［１］张正禄，郭际明，黄全义，等．超站仪定位系统概述 注册前述的各项碎部测量功能命令。 ［Ｊ］．测绘信息与工程，２００１（２）：４Ｏ一４２． ［２］郭际明，罗年学，张正禄，等．ＨＤ—ＳＴＧＰＳ超站仪的实 ３结束语 时坐标转换［Ｊ］．测绘信息与工程，２００１（４）：４１—４３． 基于ＡｕｔｏＣＡＤ的面向对象的二次开发工具一Ｏｂｊｅｃｔ． ［３］郑家健．利用Ａｕｔｏ ＣＡＤ制作台阶程序［Ｊ］．辽宁测绘， ＡＲＸ，采用面向对象的ＶｉｓｕａｌＣ＋＋程序设计技术，实现了 ２００９（１１：５０—５２． ＨＤ－－ＳＴＧＰＳ超站仪的碎部测量数据采集系统。系统作为 ［编辑：胡雪］
