胡得林 胡 勇* 于 洋 何朋朋 俞 峰
(武汉理工大学船海与能源动力工程学院1) 武汉 430063) (江南造船(集团)有限责任公司2) 上海 201913)
船体曲板的制造是船舶建造中的重要环节，目前船体曲板的加工主要依赖水火弯板工艺和机械冷冲压工艺[1-2].船体曲板成型后，需要对曲板进行测量，用以计算回弹补偿量和确定制造精度.传统的接触式船体曲板测量依靠人工样箱，测量精度低，且测量结果不能量化.为此，需要对船体曲板测量方法进行改进，以满足大尺寸船体曲板高精度测量要求，实现大尺寸船体曲板的自动化测量.
随着计算机视觉和摄影技术的发展，基于计算机视觉的摄影测量技术在船舶制造领域得到了广泛的应用.Goldan等[3]利用近景摄影测量技术与CAD/CAM系统结合，测量在建船体曲板，提高了测量效率和精度，并且适合大尺寸船体曲板的测量，但该方法需要辅助靶标，对拍摄要求较高.Heo等[4]将多条结构光投影到船体曲板表面，实现了多结构光同步测量，该方法测量精度和效率都很高，但只适合小尺寸船体曲板的测量.Paoli等[5]将基于面结构光的双目立体视觉传感器安装在机械臂上，实现了对船体表面的全方位测量，该方法测量范围大、自动化程度高，但测量设备昂贵.王振兴[6]利用基于结构光的双目立体视觉技术对船体曲板进行在位测量，实现了加工过程中船体曲板的自反馈检测，但该方法测量场景偏小.Deng等[7]用三维激光扫描仪扫描船体曲板，得到高精度点云，直接获取船体曲板三维数据，但测量效率较低.考虑到船体曲板具有线性复杂、尺寸较大的特点，因此需要研究出一种能适应大尺寸测量且测量效率高的测量方法.
文中提出了基于计算机视觉的大尺寸船体曲板摄影测量方法，通过拍摄目标物体，将获取的数字图像通过算法进行特征提取与匹配，还原特征的三维坐标，得到目标物体的三维点云数据.相比于传统摄影测量方法，投入成本低，测量视场大，可完成大尺寸船体曲板的摄影测量，且测量精度小于1mm.文中整体算法框架见图1.
图1 本文算法框架
摄影测量是基于立体像对的前方交会视觉成像模型，利用立体像对上的同名像点，得到两条同名射线在空间的相交点，从而计算出目标点的三维坐标.图2为双目摄影测量的基本原理图，Ol、Or为左右摄像机的光学中心.设世界坐系下上任意一点P，经过摄像机成像，在左右摄像机的成像平面Il、Ir分别对应有Pl(xl,yl)、Pr(xr,yr)两个同名像点.左右摄像机的光轴Olzl、Orzr互相平行，f为摄像机焦距，PlPr满足极线约束，d=Δl-Δr为视差.
图2 双目摄影测量原理图
根据视差原理，P点在摄像机坐标系下可以表示为
(1)
由式(1)可知，在满足极线约束的情况下，只要知道P点在图像中的位置，并知两台摄像机的内参数和外参数，就可以得到P点在摄像机坐标系下的位置，进而通过坐标转换还原出P点的实际位置.摄像机内外参数的获取必须通过实验和计算求解，这个过程就是摄像机标定.P点在图像坐标系下的位置，则需要通过特征匹配算法得到.
1.1 摄像机标定
摄像机成像的基本数学模型描述了一个物理过程：物理世界中的一点P，经过相机成像呈现在图像上的像素坐标为(u,v)，见图3.点P经过坐标转换在图像坐标系Oi-xy下的坐标为(x,y)；
经过透视变换在摄像机坐标系Oc-xcyczc下的坐标为(xc,yc,zc)；
通过投影变换在世界坐标系Ow-xwywzw下的坐标为(xw,yw,zw).
图3 摄像机成像模型
在不考虑摄像机光学畸变的情况下，设Oi在像素坐标系Op-uv下的坐标为(ui,vi)，由透视变换理论和投影变换理论可以推导出物理世界任意一点P(xw,yw,zw)经过摄像机成像最终投影为摄像机成像平面上一点Pi(x,y)(其像素坐标为(u,v))，其数学模型为
(2)
式中：dx、dy分别为一个像素点在x、y方向上的单位物理尺寸；
s为一比例因子，f为摄像机焦距；
I为摄像机外参矩阵；
T和R分别为平移向量和旋转矩阵，构成摄像机外参矩阵E.
由于制造工艺等因素影响，摄像机成像过程会存在光学畸变.畸变会导致P点在摄像机图像坐标系下的位置与理想状态下的位置存在偏差，见图3的Pr(xr,yr).由相机的畸变校正模型，两个位置之间的计算关系为
(3)
由式(2)～(3)可知:为确定物理空间物体表面某点的三维几何位置与其在摄像机成像平面上对应点之间的关系，必须建立起摄像机成像的几何数学模型，这些数学模型的参数就是摄像机的参数.本文采用张正友标定算法，具体步骤如下.
步骤1打印标定图案贴到一个平面上或者购买标准的标定板，利用标定板的特征点坐标作为已知条件.
步骤2移动平面或者摄像机，并从至少两个角度拍摄标定板图像，建立和已知条件的映射.
步骤3在计算机上用算法实现对图像中特征点的亚像素级检测.
步骤4在不考虑摄像机畸变的情况下，估计出内参矩阵和全部外参.
步骤5考虑摄像机畸变，通过最小二乘法计算畸变参数.
步骤6通过收敛函数改进全部参数，使得误差最小.
1.2 特征匹配与三维重建
目标物体的高精度摄影测量需要高质量的点云数据，这些数据来源于立体匹配.基于特征点的立体匹配具有良好的尺度不变性和较快的匹配速度，应用较为成熟.
在特征提取阶段，AKAZE算法采用与尺度不变换特征(scale invariant feature transform, SIFT)算法类似的模型，即在同一金字塔层内的不同尺度空间的一组图像中寻找最大特征点.不同的是，AKAZE算法采用的是非线性扩散滤波来构建尺度空间，从而保留图像更多的边缘特征.在特征描述子生成阶段，采用M-LDB特征描述算法[8]来生成二进制描述子，使得最终得到的特征具有旋转不变性，提高了鲁棒性.AKAZE算法的核心步骤[9]如下.
1) 非线性扩散滤波 非线性扩散滤波方法是将图像亮度L在不同尺度空间上的变化作为热扩散函数的散度因子，通常采用非线性偏微分方程进行求解，其表达式为
(4)
式中：c(x,y,t)为传导函数；
t为图像尺度参数；
div和的分别表示散度和梯度.
2) 构造非线性尺度空间 AKAZE算法特征的尺度空间构造与SIFT算法类似，都有O组图像，每组图像有S个子层级，尺度级别按对数递增，各个层级采用与原始图像相同的分辨率，并通过下面的式子与尺度参数σi相对应.
(5)
式中：o∈[0,…,O-1]；
s∈[0,…,S-1]；
i∈[0,N-1]；
N为图像总数；
σ0为图像初始尺度参数.
综合式(4)～式(5)，对于一副输入图像，通过加性算子分裂算法(additive operator splitting，AOS)可以得到非线性尺度空间的所有图像，图像高度L的尺度空间表达式为
Li+1=(I+τA(Li))Li
(6)
式中：τ为步长；
I为单位矩阵；
A(Li)为图像在维度上的传导矩阵.
3) 特征点检测 图像的非线性尺度空间构造后，AKAZE算法通过寻找不同尺度归一化后各像素点对应的Hessian矩阵局部极大值，来实现特征点的检测.比较Hessian矩阵值与同层八个领域点，以及上下邻层18个像素点的响应值大小，得到的极小值点就是图像在该尺度下的特征点，见图4.二维图像中Hessian矩阵表达式为
(7)
图4 特征点检测
4) M-LDB描述子 AKAZE算法根据特征点的局部图像结构来确定其主方向，保证了图像的旋转不变性.图5为特征点主方向，特征点的尺度参数为σi，则以6σi为半径搜索，以60°扇形为一个滑动窗口，按照近大远小(即离特征点越近，贡献值越大，反之则越小)的原则，对像素点进行向量加权叠加，向量之和最大的向量方向即为特征点主方向.
图5 特征点主方向
确定了特征点的主方向，AKAZE算法采用M-LDB描述子对特征点进行描述.对网格像素进行尺度自适应，比例采样后以平均值来近似划分网格像素，进行尺度采样的描述子对尺度变化具有更好的鲁棒性[10].
当两幅图像的特征点和特征向量生成后，采用KNN算法对AKAZE特征点的描述子进行匹配，再使用RANSAC算法进行误匹配点消除，得到理想的特征匹配.
1) KNN特征匹配算法 KNN算法采用测量不同特征向量之间的欧式距离进行相似判断度量，因此具有精度高、对异常值不敏感等优点.文中取K=2，即对于特征点p，在另一幅图像中搜索匹配点q、q的最邻近点q1和q的次邻近点q2.设p与q1的欧式距离为d1，与q2的欧式距离为d2.当远d1小于d2时，则认为q为有效匹配点；
当d1与d2接近时，设置一个阈值范围，若d1/d2在阈值范围内，则认为q为有效匹配点，否则不是有效匹配点.实验表明，阈值范围为0.6～0.8时,算法效果最好.
2) RANSAC算法 采用KNN算法进行粗匹配会后出现特征点的误匹配，采用RANSAC算法进行匹配优化.RANSAC算法的原理是依据迭代算法从参考模型中选取样本集合，然后计算所有样本与该模型的偏差，再使用一个预先设定好的阈值与偏差比较，当偏差小于阈值时，该样本点属于模型内点，否则属于模型外点，迭代这一过程，最终找出性能更加优异的参数模型，实现精确匹配.
图像的特征匹配只能获取稀疏的点云，无法完成物体的三维重建.准稠密匹配算法是将特征匹配得到的稀疏匹配点作为种子点集合，再次利用RANSAC算法设置阈值从集合中得到种子匹配点，这些种子匹配点在邻域内进行扩散来获取更多新匹配点，最终得到稠密匹配结果.利用视差法对所有匹配点进行计算，还原出匹配点在世界坐标系下的三维坐标，得到目标物体的点云数据，实现三维重建.此部分原理在1节已经阐述，这里不再赘述.
2.1 相机标定实验
按照上述摄像机标定步骤，在位置1和位置2分别对标定板进行拍摄.标定板规格为6×9，单格尺寸为40 mm.为保证拍摄图像最清晰，先使用自动对焦模式，使得图像最清晰.然后将摄像机定焦，通过移动标定板，在不同角度一共拍摄14组照片，图像分辨率为2 736×1 842，标定板棋盘格见图6.
图6 标定板
基于OpenCV库在Visual Studio 2017上实现张正友标定法，对摄像机进行标定.标定过程分为两步，第一步为对两个位置的摄像机分别进行单目标定，得到摄像机的内参数，见表1.第二步进行摄像机双目标定，得到摄像机的外参数，即旋转矩阵R和平移向量T.
表1 摄像机标定参数
旋转矩阵R和平移向量T的标定数值为
T=[-848.693 1 1 892.277 6 216.282 1]
实验表明，本标定算法可行、有效，重投影误差RMS达到0.05个像素左右，算法具有较高的标定精度，可以进行下一步摄影测量实验.
2.2 船体曲板测量实验
为了从二维图像中获得被测物体特征的三维坐标，基于结构光的双目摄影测量系统至少需要从不同位置获取包含物体特征点的两幅图像.事实上，获取两幅图像并不一定需要两个摄像机，由一个摄像机通过运动，在不同位置观测同一个静止的物体，都可以满足要求.采用单摄像机，降低了系统的成本，根据摄像机的移动的位置的不同，很容易构成不同基线距的双目视觉测量系统，具有很大的灵活性.为增加船体曲板表面的特征点，用投影仪投射结构光，本次投射的结构光为棋盘格，见图7.使用标定好的摄像机在位置1和位置2分别对船体曲板进行拍摄.
图7 摄影测量系统
摄像机位置1和摄像机位置2为相对固定，保证船体曲板在视场中即可，本次实验中B=1 000 mm；
摄像机与船体曲板的距离值为r=1 500 mm；
投影与船体曲板的距离不设值，根据船体曲板大小，保证投射的结构光完全覆盖船体曲板，并且投影到船体曲板上显示最清晰.此次拍摄的船体曲板尺寸大约为1 500 mm×2 000 mm，用摄像机分别在位置1和位置2拍摄船体曲板，图像的分辨率为2 736×1 824，见图8a).图像采集完成之后，在Visual Studio 2017上进行图像处理，见图8b)～d).
图8 船体曲板摄影测量结果
2.3 船体曲板测量精度分析
标准块是边长为100 mm的标准立方体，实际尺寸由游标卡尺测得，精度为±0.02 mm.放在与船体曲板同一个视场中进行测量.用标定后的相机对标准块进行拍摄，标准块模型图、实物图见图9a).按照文中的摄影测量方法对标准块进行测量，得到标准块的点云数据，见图9b).
图9 标准块测量
从图9中标准块的点云数据中采样a、b、c、d点的坐标，见表2.并计算ab、bc、bd的距离，与标准块的实际边长比较，计算其绝对误差，见表3，以此来分析本文船体曲板摄影测量的精度.
表2 标准块角点坐标 单位：mm
表3 标准块边长误差 单位：mm
由表2～3可知，测量结果的误差在1 mm以内，测量精度满足船舶行业工程制造的精度要求.
文中针对大尺寸船体曲板测量提出了一种非接触式基于计算机视觉的船体曲板摄影测量方法，搭建了摄影测量系统进行船体曲板的测量实验.利用张正友标定法对摄像机进行了标定.在图像处理方面，利用AKAZE算法进行船体曲板的特征点提取；
利用KNN机器学习算法与RANSAC算法相结合进行特征匹配算法；
用获得匹配点对，基于特征匹配点进行匹配扩散，获得稠密匹配，实现对船体曲板的三维重建. 试验结果显示，在测量场景不小于3 000 mm×3 000 mm，船体曲板尺寸为1 500 mm×2 000 mm范围内，摄影测量误差在1 mm以内，满足船舶制造的精度要求.在实验过程中发现了在结构光投影到目标物体显示有重影处，特征提取和匹配的误差相对较大，为了进一步提高测量精度，后续工作将围绕如何将这类点进行精确提取与匹配展开.