普通数码相机数字摄影测量的数理模型研究.doc

第 22 卷 第 2 期2000 年 6 月武汉水利电力大学 (宜昌) 学报J . of Univ. of Hydr. & Elec . Eng. / YichangVol . 22 No . 2J un. 2000普通数码相机数字摄影测量的数理模型研究魏明果刘德富刘运长王仁明(基础课部)摘要给出了普通数码相机数字摄影测量软件主体部分的开发思路 ,对软件中的两大主要部分直接线性变换与影像匹配两部分分别给出了数学模型 ,为检验软件的整体性能 ,又给出了一个物理模型 1 通过数值模拟及物理模型的计算结果 ,表明软件在以上两部分的开发研究是成功的 1关键词普通数码相机 ; 数字摄影测量 ; 软件开发分类号P2321影像数据格式的转换普通数码相机拍摄的照片为彩色的 ,而影像匹配一般是用像素的灰度性质作为匹配的基础 ,故数字摄影的第一步须将数码相机的数字影像的彩色信息转化为 MATLAB1 ,2 软件平台上的灰度信息3 ,4 1 在MATLAB 软件平台上 ,可先用 ADOB E PHOTODELUXE 软件将彩色数字影像转存为 bmp 格式的灰度图像 信息 ,再用 MATLAB 语言获取其数据矩阵和强度式矩阵 ,方法为 X , ma p = imread ( filename1bmp )(1)其中 , X 为 M N 阶 8 位型整数矩阵 , ma p 为 M 3 阶的有关红 、绿 、蓝三色组合的强度式矩阵 1 两矩阵的关系为 X ( i , j) 的值为 ma p 矩阵的行号 1 由于 8 位型整数矩阵不能参加运算2般的整形或双精度 1 设 XDOUB 为普通双精度型 ,转化方法为XDOUB = DOUBL E ( X) + 1XDOUB 并不是影像匹配中需要的灰度矩阵形式 ,为此使X1 ( i , j) = ma p ( XDOUB ( i , j) , 1)则 X1 ( i , j) 即为通常意义下的灰度矩阵 1,故须先将其转化为一(2)(3)真接线性变换法与数学模型2专业型的数码相机能提供 3 个内方位元素 ,解算中只要通过设置一定的控制点由基本的成像方程即可获得其相机的外方位元素5 ,而普通型的数码相机由于不能提供其像主点的坐标 ,因而其内方位元 素也需要由设置的控制点来获得 1 文献 6 中给出了一种直接线性变换的方法 1 其计算公式为= ( B T B ) - 1 B T W( 4)V = BL WLV xV y其中 V =为误差矩阵 1ZA 10XA0YA0ZA0- xi XA1- yi XAXA YA- xi YA- yi YA- xi ZA- yi ZA( 5)B =0000 国家电力公司重大科技攻关项目 ( SDZ9802 - 2) 及霍英东基金 ( SDZ9801) 资助项目部分成果收稿日期 :2000 - 04 - 05魏明果 ,教授 ,主要从事应用数学 、数字 、摄影测量 、地质 CT 等教学与研究 1 武汉水利电力大学 ( 宜昌) 基础课部 (443002)1AxiyiTL 1 , L 2 , L 11i = 1 , 2( 6)L =W =XA 、YA 、ZA 为物方控制点坐标 1 xi 、yi , i = 1 , 2 为控制点的像平面坐标 , L i , i = 1 , 2 为待求的与内 、外方位元素有关的待定系数 1 A = L 9 X + L 10 Y + L 11 Z + 1 , 在实际计算中 ,若把物方坐标系的原点选取在被 测点附近 ,此时可令 A = 1 近似求解 1其正算公式为S = ( XA , YA , ZA ) T = ( N T N ) - 1 N T Q( 7)其中L i1L i5-xiL i9yiL i9L i2L i6- xi -yi -xiL i10yiL i10L i3L i7-xiL i11yiL i11i = 1 , 2(8)N =L i4L i8i = 1 , 2( 9)Q =数学模型如图 1 所示 1 其中 S 1 、S 2 分别为两个水平设置的左 、右站点 , P1 到 P6 为物方的 6 个控制点 1 在数学模型中已知物方坐 标 ,它在左 、右两片的像平面坐标可由下式给出( 10)X =A + B M其中xizixi -zi -XjZjxzX =, A =, B =, Yj + f - yi, i = 1 , 2 , 6 , j= 1 , 2(11)M =y - Yi j当 j = 1 时 , X 表示左片像平面坐标 ;当 j = 2 时 , X 表示右片像平面坐标 1在成像区域中 ,给出 6 个控制点坐标 ,由式 (4) 可分别求出左 、右两站点的与内外方位元素有关的系 统系数 LL 和 L R1 其中LL =- 010250010000010000- 010750010000010000- 011667010000L R =- 010250010000010000010750010000010000- 011667010000- 010250010000010000- 010250010000010000当求得 LL 和 L R 的情况下 ,在左 、右两片上任取一组同名像点 ,由式 (7) 就可以求得它们的物方坐标 1表 1 为由 6 对同名像点反算它们相应的物方坐标与实际物方坐标的结果对比 1表 1 像平面坐标与物方坐标数据表序号左片像平面坐标右片像平面坐标计算所得物方空间坐标实际物方坐标123x = - 010545 x = - 010500 x = - 011282x = - 010750y = 010273y = 010306y = - 010327y = 010000 x = 011091 x = 011167 x = 010355x = 010750y = 010273y = 010306y = - 010327y = 010000(013000 , 017000 , 014000)( - 112000 , 016000 , - 111000) (117000 , 015000 , 112000)(011943 e - 14 , - 016106 e - 14 ,012186 e - 14 ) (112000 , 013000 , - 110000)( - 110000 , 014000 , 112000)(013 , 017 , 014)( - 112 , 016 , - 111) (117 , 015 , 112)(010 , 010 , 010)456x = - 011105x = - 010536y = 012603y = - 010321x = 010474x = 011071y = 010263y = - 010321(112 , 013 , - 110)( - 110 , 014 , 112)3影像匹配与最小二乘影像匹配设 a1 与 a2 分别为左 、右两片中的像点 , 左片记为 I1 , 右片记为 I2 , 其中 ai ( i = 1 , 2) 的像平面坐标分别为 Pi = Pi ( xi , yi ) , ( i = 1 , 2) 1 ai ( i = 1 , 2) 的特性函数为 gi ( i = 1 , 2) 1 特性函数 gi 可以是以 Pi 为中心1A1A1A176武 汉 水 利 电 力 大 学 (宜 昌) 学 报2000 年 6 月的某一影像窗口的灰度矩阵 1 基于 gi ( i = 1 , 2) 定义某种测度 ( g1 , g2) 1 测度 表示像点 a1 与 a2 的匹配的程度 1匹配点的确定取决于匹配测度 的定义 1 本文中的测度 取为2= ( g1 g2) 2/g12 g22(12)当测度取最大值时 ,认为左片影像窗口中心的像点 a1 与右片影像窗口的 a2 为同名像点7 1设 a1 为目标点 ,以 a1 为中心在左片上组成一个 m n 个像素的影像窗口 ,称此影像窗口为目标窗 口 1以目标点在左片上的位置 ,可以根据先验知识 ,在右片相应位置建立一 K L ( K m , L n) 个像素作为搜索区 1 在右片的搜索区中依此取出 m n 个像素 ,并计算其与目标区的相似性测度 ij ( i = i0L n , i0 + L - n ; j = j0 - + m ,k, j0 + - m ) 1 其中 ( i0 , j0) 为搜索区的中心像素 1 由式k- +22222222(12) 当 取得最大值时 ,该搜索窗口的中心像素被认为是同名像点 1 这种方法称为二维搜索影像匹配1最小二乘影像匹配是以左片 ( 或右片) 作为匹配的目标片 ,在左片上任取一点 M ( x , y ) ,其灰度为g1 ,设 M2 ( x2 , y2) 为右片上的同名点 ,其灰度为 g2 ,右片经几何形变与幅射形变有如下关系x2 = a0 + a1 x + a2 yy2 = b0 + b1 x + b2 yg1 ( x , y) + n1 ( x , y) = h0 + h1 g2 ( x2 , y2) + n2 ( x , y)其中 n1 ( x , y) , n2 ( x , y) 为随机误差将式 (15) 线性化后可有最小二乘影像匹配误差方程v = c1 d h0 + c2 d h1 + c3 d a0 + c4 d a1 + c5 d a2 + c6 d b0 + c7 d b1 + c8 d b2 - g( 13)( 14) (15)( 16)其中 , c1 、c2 、c3 、c4 、c5 、c6 、c7 、c8 为二阶 TAYLOR 展开式的系数 ,g 是相应像素的灰度差 ,d h0 、d h1 、d a0 、d a1 、d a2 、d b0 、d b1 、d b2 为待定参数的改正值 1 记X =h0 h1 a0 a1 a2 b0 b1 b2 T则有最小二乘影像匹配的误差方程V = CX - L其中 , L = g1 ( x , y) - g2 ( x2 , y2) ,由误差方程可得法方程式D X = ( C T C) - 1 ( C T L )(19) 式是一个严重病态的方程组 ,可用阻尼最小二乘方法求解 1考虑数学模型如下 :其中左片灰度函数按如下方法生成g1 ( x , y) = x012 + y013( 17)(18)( 19)(20)(21)0105113T取X = 015201618112010120175按 (13) 、(14) 、(15) 式方法可生成右片灰度矩阵 1在左片上任取一点 ,在右片上寻找同名像点 ,取初值 X = 0 1 0 1 0 0 0 1T 当迭代次数 时 ,以 5 5 个像素为搜索窗口 ,其解收敛于K 大于 30017887010526113671TX = 0151197016664112940010129比较 (21) 、(22) 两式 ,可见 (22) 式基本收敛于 (21) 式 1(22)4软件流程与步骤由两张普通数码相机的平面彩色影像求取物方三维空间坐标计算过程如下 : (1) 由彩色数据影像格式转化为灰度数据影像格式 ; (2) 输入控制点的三维地面坐标及控制点的二维像平面坐标 ; (3) 由控制点 的同名像平面坐标应用直接线性变换法 (4) 式计算系数 L ; (4) 进行影像匹配寻找左右两片的同名像点 ;(5) 读取左右两片的同名像点坐标 ,由公式 (7) 求取物方的三维空间坐标向量 S15物理模型图 1 、图 2 为一物理模型的左 、右两张灰度照片 ,该模型是在一张为 75 105 ( cm2 ) 坐标纸上 ,设置了若干个直径为 4mm 、高度不等 、顶端颜色各异的小木棍 1 设置的地面坐标的原心在坐标纸的中心点处 1 采用 KODAK265 数码相机在距模型 118m 处拍摄了两张彩色照片 ,经数据格式处理后转化为如图 1 ,图 2 的 750 900 个像素的灰度照片 1 其中各点的坐标与实际高程如表 2 所示 1 各点经本软件处理后的物方坐标均见表 21图 2物理模型左片灰度图图 3 物理模型右片灰度图表 2 同名点三维坐标左片同名点右片同名点横坐标 纵坐标横坐标 纵坐标计算的物方空间三维坐标实际的物方三维空间坐标- 15010990- 1491423601412510011775- 15012129- 012863- 01413115014352- 150131499912950- 14911502- 4613161913674151693219129014416688511439656105651234562165854055312053984012065942696004661238910273656469603270470406209597270607468- 1501500100- 15000150- 150100- 150- 50101520455055 7 395 211 12 465 212 013314 15013658 6412512 0 150 65 6结语数字摄影测量系统是一个非常复杂的系统 ,由于从二维平面影像恢复三维空间坐标在理论上是一个病态的问题 ,因而还存在着许多需要研究的问题 ,目前数字摄影测量正处于发展的初期7 ,加之对不 同的工程问题其处理方法上的差异 ,因而本文只就普通数码相机在 MATLAB 软件平台上的一些基本的 方法给出了相应的数学或物理模型 1在数字摄影测量系统中 ,影响最后计算结果精度的因素很多 ,但成像区域中控制点的设置 、控制点像平面坐标的测量及最小二乘影像匹配中的插值方法的选取等也是其中重要的因素 1 由于式 (19) 是一 个严重病态的方程组 ,故插值结果的好坏对影像匹配的结果有着决定性的影响 ,因而如何有效地提高插值的精度是最小二乘法的关键 1 另外在控制点的选取上应尽量使得所需匹配点位于控制点所成区域的 内部 ,实践表明 ,所要匹配点越靠近控制点反算结果精度越高 1 控制点像平面坐标是数字摄影测量的基 础 ,各个要素的求取都是以它的结果为计算的根本 ,在实践中控制点坐标的测量除可先用光标定位法外 还当结合特征点定位等多种方法通盘考虑 1178武 汉 水 利 电 力 大 学 (宜 昌) 学 报2000 年 6 月参考文献1 张宜华 1 精通 MATLAB 51 北京 :清华大学出版社 ,19992 高俊斌 1MATLAB 510 语言与程序设计 1 武汉 :华中理工大学出版社 ,19983 史文革 1 微机图像格式大全 1 北京 :海洋出版社 ,19924 王润生 1 图像理解 1 长沙 :国防科技大学出版社 ,19955 金为铣 ,扬先宏 1 摄影测量学 1 武汉 :武汉测绘科技大学出版社 ,19966 丁 辋 ,刘友光 ,张松波 1 工程与工业摄影测量 1 武汉 :中国地质大学出版社 ,19957 张祖勋 ,张剑清 1 数字摄影测量学 1 武汉 :武汉测绘科技大学出版社 ,1996Study of Mathematic and Phys