(19)国家知识产权局 (12)发明 专利申请 (10)申请公布号 (43)申请公布日 (21)申请 号 202211238599.3 (22)申请日 2022.10.11 (71)申请人 天津工业大学 地址 300387 天津市西青区 宾水西道39 9号 (72)发明人 王红一　窦公成　宋丽梅　 (51)Int.Cl. G06T 17/00(2006.01) G06T 7/80(2017.01) (54)发明名称 一种应用于气泡三维流场重建的改进MLEM 方法 (57)摘要 本发明属于气泡三维流场重建领域， 涉及一 种应用于气泡三维流场重建的改进MLEM方法法。 本方法将多目视觉测量与改进的MLEM方法相结 合， 首先利用多台相机同步拍摄的多组标定板图 像对多相机系统进行标定， 将多相机的世界坐标 系进行统一。 然后利用多台相机对气泡流场进行 同步拍摄， 得到同一时刻多个角度的气泡图像。 最后针对气泡流场背光照射的特点， 改进MLEM方 法中迭代初值的设定方法及迭代过程中的投影 和反投影方法， 实现了基于改进MLEM方法的气泡 三维流场快速重建。 本发明所设计的应用于气泡 三维流场重建的改进MLEM方法， 在提高了计算速 度的同时， 也提升了重建精度， 可 以广泛应用于 气液两相流中气泡 三维流场重建领域。 权利要求书4页 说明书8页 附图4页 CN 115457217 A 2022.12.09 CN 115457217 A 1.本发明涉及一种应用于气泡三维流场重建的改进MLEM方法， 其特征是， 包括下列步 骤： 步骤1： 将K台相机围绕拍摄目标放置， 且处于同一水平面上， 相机的整体拍摄角度小于 180度； 步骤2： 在被拍摄的目标区域中放置棋盘格标定板， 多次改变棋盘格标定板的位置， 利 用K台相机对棋盘格标定板进行同步拍摄， 得到多组同一时刻K个角度的棋盘格标定板图 像； 步骤3： 对各相机进行 单目标定， 获得 各相机的内外参数信息； 在标定时， 计算机会获取每个相机拍摄的棋盘格图像的角点信息， 而每个相机将其拍 摄的最靠近左上角的棋盘格角点坐标定义为世界坐标的原 点， 棋盘格水平方向为X轴， 棋盘 格垂直方向为Y轴， 垂 直于棋盘格平面的方向为Z轴， 以此为世界坐标系， 计算外参时将以此 世界坐标系为基准； 由于标定过程中要改变棋盘格的位置， 利用多张图像进行标定， 因此， 不同图像所对应的世界坐标也是不同的， 计算得到的外参也不相同； 为了让K台相机的世界 坐标系具有统一的基准， 以K台相 机在某一时刻同步拍摄的标定图像为基准来计算各相 机 的外参； 相机的内部参数表示相机坐标和像素坐标之间的关系， 且只由相机本身决定， 可以 唯一确定； 分别对各台相机进行单目标定， 得到第n(1＜ ＝n＜＝K)台相机的内参矩阵Mn、 旋转矩阵 Rn和平移向量Tn， 此时各台相机的旋转矩阵和平移向量均 是相对于同一位置的棋盘格标定 板所确定的世界坐标系； 步骤4： 以相机1作为参考相机， 利用步骤3中所得到的相机内外参数， 计算其它相机与 相机1之间的位置关系参数； 对相机1和相机n进行双目标定， 计算相机n相对于相机1的位置关系参数； 假设空间中 某一点为q， 已知其世界坐标为qw， q点在相机1和相机n的相机坐标系下的坐标为 q1＝R1·qw+T1 公式(1) qn＝Rn·qw+Tn 公式(2) 由公式(1)和公式(2)， 可 得两个相机坐标q1和qn的关系为 qn＝R1n·q1+T1n 公式(3) 其中， T1n＝Tn‑R1n·T1 公式(5) R1n和T1n分别是相机1和相机n 之间的旋转矩阵和平移向量， 即为所求的相机间位置关系 参数； 步骤5： 以相机1的世界坐标系为基准， 利用步骤4中所述的相机1和相机n的位置关系矩 阵， 将相机n的世界坐标系转换到相机1的世界坐标系中， 则相机n的外参矩阵为 步骤6： 世界坐标系(Xw， Yw， Zw)和像素坐标系(u， v)之间的关系如公式(7)所示；权　利　要　求　书 1/4 页 2 CN 115457217 A 2其中， Zc是当前坐标点在相机坐标系的Z轴的坐标值； M是相机 的内参矩阵， 包括fx， fy， u0， v0， s等参数， fx和fy表示相机在像素坐标系和相机坐标系之间转化的尺度因子， (u0， v0) 是图像坐标系的原 点在像素坐标系中的对应坐标， s为倾斜因子； N为相机外参矩阵， 包括旋 转矩阵R和平 移向量T； MN 为相机投影矩阵， 矩阵维度为3 ×4； 将步骤3计算得到的各相机的内参和步骤4计算得到的各相机的外参矩阵代入公式 (7)， 可以得到世界坐标系和各相机的像素坐标系之间的投影关系表达式， 如公式(8)所示； 式中， (un， vn)是相机n的像素坐标系中的像素点 坐标； Mn是相机n的内参矩阵； 步骤7： 在完成步骤6所述的内容之后， 保持相机的相对位置不动， 在被拍摄区域中放置 水箱， 水箱中通入气泡， 利用K台相机对气泡流场进 行同步拍摄， 得到同一时刻K个角度的气 泡图像； 步骤8： 基于步骤7所采集到的K个角度气泡图像， 对图像进行去噪和二值化处理， 得到K 个角度二 值化后的气泡图像I1， I2， ...， In， ...， IK； 步骤9： 设置迭代重建方法的最大迭代次数Item和最小允许误差er rormin； 步骤10： 将步骤6所述的公式(7)中的投影矩阵MN用公式(9)来表示； 其中， Cab为投影矩阵中的矩阵元 素， a＝1， 2， 3， b＝1， 2， 3， 4； 则相机的世界坐标系到像素坐标系的投影关系， 可表示 为公式(10)所示的形式； 将三维空间区域按照世界坐标系的Zw轴进行分层， 把三维空间区域看作是连续多层的 二维图像， 对已知的Zw层， 其反投影 计算公式如公式(1 1)所示； 将二维像素坐标(u， v)及当前层的Zw坐标代入反投影计算公式， 可计算出像素坐标为 (u， v)的点的在Zw层上对应点的Xw， Yw坐标； 步骤11： 基于步骤8所述的处理后的K个角度的二值化 图像， 分别对每个角度二值化图权　利　要　求　书 2/4 页 3 CN 115457217 A 3