张宝祥，玉振明，杨秋慧

（1. 广西大学 计算机与电子信息学院，广西 南宁 530004；2. 梧州学院 机器视觉与智能控制广西重点实验室，广西 梧州 543002；3. 梧州学院 大数据与软件工程学院，广西 梧州 543002）

1 引言

近年来，随着人工智能和机器视觉技术的快速发展，基于图像的三维重建方法越来越受到关注，但在面向显微镜应用领域，其主要应用于平面标本的观测，鲜有针对立体标本成像以及建模的显微设备。研发新一代智能显微设备的同时，开发出面向立体标本的三维模型重建和立体标本几何参数测量的技术具有重要意义。

作为机器视觉的热门方向，已有诸多方法和理论用于基于图像的三维重建，其中，运动恢复结构法（Structure from Motion，SfM）［1］是应用最广泛的经典方法之一。基于SfM 设计出了多种三维建模系统，如Wu 开发的开源软件Visual SfM［2］，Moulon 等 人 开 发 的 开 源 软 件Open-MVG［3］，Schönberger 等 人 开 发 的 开 源 三 维 建 模系 统COLMAP［4，5］。蒋 华 强 等［6］针 对 现 有SfM算法重建模型存在的点云稀疏等问题，提出一种利用不同匹配数据进行模型重建的算法。但由于SfM 仅计算图像之间成功匹配的特征点的三维信息并将其恢复为三维坐标形成三维点云，而图像中含有的特征点信息相对较少［7］，因此SfM 计算得到的点云模型比较稀疏，重建后的模型精度较低。多视图立体方法（Multi-View Stereo，MVS）［8］可以从物体多个视角图像计算得到场景的密集三维点云。PMVS（Patchbased MVS）［5］以SfM 重建的稀疏点云为输入，利用图像面邻域信息迭代，使用点云扩张策略进行点云扩张和过滤，最后重建出稠密点云。Schönberger 等人［9］以SfM 得到的稀疏重建模型和相机姿态为输入，采用深度图融合的方法来恢复出稠密点云。基于学习的MVS 方法如Point-MVSNet［10］，R-MVSNet［11］，Cascade-MVSNet［12］等采用深度图融合的方法恢复较高精度的场景稠密点云也取得不错的结果。

以上三维重建算法对拍摄设备要求较高，适用于对大场景成像的纹理等细节丰富的场景，而在实体显微镜应用场景中，由于显微镜工作环境、珠宝等立体标本在显微镜局部放大观测时存在弱纹理、细节不丰富和反光等问题，运用上述算法得到的三维重建效果并不理想。本文提出了一种适用于显微镜应用场景下基于特征提取的多视图立体三维重建算法。以基恩士VHX-6000 高级体视显微镜作为实验平台，以集成电路板元器件或珠宝首饰作为实验对象，通过显微镜采集物体不同角度的图像序列（大于10 张），对立体标本显微镜图像进行三维重建，取得较好效果。

2 相关工作

2.1 图像采集与预处理

显微镜具有高放大率和窄场的特点，对立体标本的图像采集往往只能获得物体较小的局部放大图像，所采集的显微图像包含的纹理信息较少属于弱纹理图像。因此，为了增强图像的细节信息和纹理信息以便后续过程可以提取到图像更多的特征信息，使用直方图均衡化方法增强图像的细节信息以使图像灰度分布更均匀［13］，使用一种融合光照和彩色信息的灰度图像补偿算法增强图像的纹理信息对图像进行灰度化［14］。本文的实验设备为基恩士VHX-6000 高端数码显微系统，如图1 所示。

图1 基恩士VHX-6000Fig.1 KEYENCE VHX-6000

2.2 改进的Harris-SIFT 匹配算法

特征点提取与匹配是基于图像三维重建算法的关键一步，特征提取与匹配的精度直接影响重建效果。Harris［15］角点检测算法提取的特征角点比较全面，具有平移和旋转不变性，对参数变化、噪声、图像对比度和光照变化等有很好的鲁棒性，但对图像尺度变化敏感，不具有尺度不变特性。尺度不变特征变换（Scale-invariant Feature Transform，SIFT）［16］算法对旋转、尺度、亮度等变化具有不变性，对视角转换、噪声、仿射变换有很高的容忍度，但由于算法实现过程中需要通过高斯模糊、降采样等方式构建高斯金字塔，得到的特征点精度不是很高，对边缘光滑的目标无法准确提取特征点，且计算量大，效率低。本文采用一种将Harris 与SIFT 算法相结合的方法［17］，首先，用Harris 进行图像中的特征角点提取，使用特征角点包含的信息熵值对特征角点进一步筛选；其次，使用马氏距离度量特征点之间的相似性；再次，用随机抽样一致性（RANSAC）方法对匹配模型进行优化，筛选匹配点对，提高匹配的正确率；最后使用显微镜采集的显微图像进行对比实验，验证本算法的匹配准确率。

2.2.1 Harris 角点检测

Harrs 角点检测算法是利用移动的局部窗口在图像中计算灰度变化值，若在任意方向上移动都存在较大的灰度变化就认为该窗口中存在角点。通过设置合理的阈值来决定检测到的角点数。Harris 角点检测变化如图2 所示。

图2 Harris 角点检测变化Fig.2 Change of Harris corner detection

窗口由像素点(x，y)移动(u，v)产生的灰度变化E可以表示为：

其中：I是图像灰度函数，W是窗口函数，将式（1）用泰勒级数化简可得：

为了在显微图像中可能存在的弱纹理区域检测到尽可能多的特征点，对图像逐像素点计算灰度变化。 使自相关矩阵M的权重， 则 有H(x，y)，表示图像每个像素点在x，y方向的二阶偏导数和x、y方向的导数的Hessian 矩阵。从图像中提取到的角点可表示为：

其中：R(x，y)为角点响应函数；det(H)为矩阵H的行列式；tr(H)为矩阵H的迹；λ为常数（取值范围0.04～0.06），在本文实验中函数R的局部极大值点即为特征角点，此时取λ=0.04。

通过式（3）初步检测出图像特征点。但由于采用逐一计算每一个像素点的局部窗口内灰度变化，因此检测结果仅反映图像局部变化情况，并会将背景中噪声点或纹理不清晰的像素点判断为特征点。熵值［18］的大小与局部区域包含的信息量成正比，背景中噪声点或纹理不清晰的像素点包含信息量较少导致计算出的熵值低，熵值大的特征点比熵值小的特征点匹配的准确率高。用以下两式进行熵计算：

其中：pij为图像某像素处的灰度值与周围像素的灰度分布的综合特征，N为图像尺度，Hf为图像二维熵，n，i，j为像素灰度值。分别对图像中初步检测出的特征角点利用式（5）在以特征角点为中心的n×n邻域内计算图像块的局部熵，统计每幅图像中特征点熵的均值，取熵均值的80%作为阈值，删除低于阈值的特征点，对特征点进行提纯。对于图像弱纹理区域统计的熵均值较低，可以删除较少的点进而在弱纹理区域保留较多的特征点，从而提高算法在弱纹理区域匹配的正确率，同时提高算法的效率。

2.2.2 生成SIFT 特征描述子

Harris 算法提取到特征角点，还需要采用SIFT 算法来描述并建立特征向量。SIFT 通过利用特征角点邻域像素的梯度方向分布特性为每个关键点指定方向参数，使算子具备旋转不变性。

其中：m(x，y)为(x，y)处梯度的幅度，θ(x，y)为(x，y)处的梯度方向，L(x，y)为所有关键点的尺度。在特征点周围n×n邻域内以尺度大小建立一个大采样区域，用式（6）、（7）计算窗口内像素点的梯度方向和大小，然后在该邻域内构建一个小采样区域，并用梯度直方图统计该小采样区域中每个“格子”内像素的梯度方向和大小，形成种子点，从而构建出特征描述子。

在本文实验中，SIFT 首先将采样区域旋转至特征点方向上，然后在其周围16×16 的邻域内计算梯度方向和大小。然后在每个4×4 的采样窗口中计算8 个方向的梯度直方图，绘制每个方向的累加值，即可形成一个种子点。如图3 所示，每个特征角点有4×4 个种子点，每个种子点有8 个方向，即每个特征角点由4×4×8=128维的描述子组成。

图3 SIFT 描述特征点Fig.3 SIFT description of feature points

2.2.3 特征匹配

为了提高对显微弱纹理图像的匹配精度，用马氏距离［19］替代原始SIFT 算法中的欧氏距离作为两幅图像中特征角点之间的的相似性判定度量。欧式距离将特征点的不同属性之间的差别等同看待，有时不能满足实际要求。马氏距离与欧式距离的不同之处在于：它考虑到特征点各种特性之间的联系；并且能对两个样本集的相似性进行有效度量，不受尺度变化及变量之间相关性的影响；它还可以排除变量之间的相关性的干扰。马氏距离可以表示为：

其中：X，Y为两个特征向量，S为两组特征向量的协方差矩阵，E为均值。从式（9）中可以发现，马氏距离计算是建立在总体样本的基础上，即计算的是两个向量间的整体距离。相比欧式距离，马氏距离多了一个协方差矩阵的逆矩阵，能更准确的度量特征点之间的相似性。虽然增加了计算量，但马氏距离把方差归一化，使得特征之间的相似关系更加合理，匹配精度更高。为进一步提高算法在弱纹理显微图像上的匹配准确率，需要剔除掉误匹配点，因此采用随机抽样一致性（RANSAC）［20］算 法 替 代 原SIFT 算 法 中 的 比 例阈值法，选择一组局域熵值高的特征点作为初始模型，用初始模型去测试所有的其它数据，如果某个点适用于估计的模型则认为它也是局内点，否则是局外点；然后，用所有假设的局内点去重新估计模型。算法通过K次迭代，可精确估计匹配模型参数。

2.3 实体显微镜场景下物体的三维重建

2.3.1 稀疏点云重建

利用基恩士VHX-6000 可变倍镜头可将智能显微镜镜头调整到一个固定的角度，然后上下、前后、左右等多方向移动载物台，从不同视角拍摄立体样本获得图像序列。本文采用如文献［4］中的增量SfM 方法重建物体的三维结构。增量SfM 是具有迭代重建组件的顺序处理方法。首先，用本文算法对图像进行特征点提取与匹配，采用包含正确匹配特征点对最多的两视图作为重建种子，通过它们之间的匹配点对估计出本征矩阵，由本征矩阵采用SVD 分解法计算出相机位姿。然后，依次将下一最佳匹配图像Inext（匹配点对数量最多）添加进模型迭代，重复估计相机位姿并对相机进行标定；同时将场景中匹配的点通过三角化来生成三维点，并计算出相机的运动轨迹和场景中物体的三维信息。最后，对所有已生成的三维点和已估计出的相机位姿使用光束平差法（Bundle Adjustment，BA）进行优化。BA 结束，得到相机估计参数和场景几何信息，即稀疏的3D 点云。

2.3.1.1 有序图像输入

由于增量SfM 需要对所有图像依次进行匹配，用无序图像进行模型重建时由于图像的无序导致图像位姿的先验条件未知，在匹配环节存在多张图片之间的交叉匹配和多次交叉验证匹配点，所以计算量大，实时性差。为提高图像的匹配效率，采用有序图像序列作为输入。通常相邻拍摄的两张图像所包含的共同特征信息最多，也更加容易正确匹配。在采集图像时根据显微镜载物台移动的轨迹以及图像拍摄时间对所采集的物体不同视角图像添加标签，将图像根据载物台运行轨迹和时间的先后排列成一组顺序固定的图像序列（类似视频中的帧）。显微镜载物台工作距离较小，拍摄的相邻图片中两两之间物体信息的重合很多，包含的特征点信息最多；在特征匹配时使相邻的图片先进行特征匹配。在SfM 算法选好初始的两张匹配图象时默认把与这两张图像相邻的图像作为下一最佳匹配图像Inext进行重建，依次将与上一张图像邻近的图像作为下一最佳输入进行迭代，这样所有图像只需添加匹配一次即可完成所有图像的重建。

2.3.1.2 三角化

SfM 在图像正确匹配和估计得到相机参数之后，就可以将图像场景之间正确匹配的特征点通过三角化方法恢复为空间三维点，进而得到三维点云。三角化的方法如图4 所示。

图4 三角化Fig.4 Triangulation

在图4 中，I1，I2为两幅图像，P1，P2分别为两幅图像上匹配的特征点；O1，O2分别为两幅图像拍摄视角时的相机中心，将每幅图像的相机中心与其特征点的连线延长在空间中交于点P，则P点为图像上的点所恢复的世界坐标系下的三维点。得到三维点之后通过BA 调整，过滤离群值，估计得到物体的三维点云。

2.3.2 稠密点云重建

SfM 重建得到的点云模型为稀疏点云，只有大致的模型结构，缺少了很多细节部分，所以还需要将稀疏点云进行稠密化。多视图立体（MVS）是一种将SfM 输出的稀疏点云稠密化的方法，MVS 重建利用SfM 算法输出的稀疏点云、相机参数及相机位姿，然后融合图像的深度图获得物体的稠密三维点云。但对于显微镜工作场景下拍摄的物体局部存在弱纹理的图像，由于表面灰度信息有限，传统MVS 方法图像深度信息估计误差较大。因此，针对上述在显微镜应用场景中的物体稠密重建问题，利用CNN 网络与MVS 结合的方法来估计和预测采集的显微图像的深度图，进而使用深度图进行稠密重建。

为保证模型收敛，本文采用文献［21］提出的基于深度残差网络（ResNet-50）［22-23］的全卷积神经网络对采集的显微图像序列进行深度估计和预测。在ResNet-50 的基础上设置了反卷积层和上采样层逐步增加每层的输出大小，并在反卷积层输出后增加了激活函数（ReLu）［24］，对输入图像进行同样比例大小的深度预测。通过比较全卷积神经网络预测的深度与MVS 所得深度的绝对差值并设置阈值，将大于阈值的预测深度融合到MVS 估计的深度图中得到改善的深度图，然后利用SfM 得到的稀疏点云和相机的位姿信息与改善后的深度图恢复出物体的稠密点云。

用于预测图像深度的全卷积网络结构如图5所示。为了增加系统的实时性并减少计算时内存的占用量，去掉了原ResNet-50 中的平均池化层和全连接层部分，仅保留ResNet-50 的卷积层结构来提取图像中的特征信息，从而可以输出尺寸为7×7，2 048 通道的特征图，减少了计算参数，加速模型的训练，提高预测效率。

图5 深度预测网络结构Fig.5 Depth prediction network structure

本文的研究旨在恢复显微镜应用场景中所拍摄物体的近似真实的三维结构。重建物体近似真实的三维结构是以正确估计图像深度为前提的。所以，为了尽可能正确地预测图像的深度信息，约束预测的图像深度值和真实值之间的误差。结合文献［21］采用平均平方误差作为损失函数，损失函数L定义为：

其中：在图像(i，j)处预测的深度值为yij，在图像(i，j)处真实深度值为tij，图像包含的像素点数为m×n，ε为常数。L将预测深度图像素点的预测值和真实值之间的平均误差定义为损失，该损失函数对预测的深度值y惩罚较严格，可以很好地减少显微镜载物台背景区域带来的深度误差影响，进一步约束了前景物体的深度误差。

3 实验结果与分析

本文的实验平台为基恩士VHX-6000 高性能实体显微镜系统，以电子电路板上元件为实体标本采集的立体标本图像序列来验证本算法的可行性。实验PC 平台配置为：系统为Ubuntu20.04.3LTS，CPU 为 英 特 尔Xeon E5-1630 V3，主频为3.70 GHz，内存为16 G，显卡为英伟达Quadro K2200，显存为4 G，硬盘为2 T 固态硬盘。如图1 所示的基恩士VHX-6000 智能显微镜可以通过人为转动手柄灵活调整显微镜物镜的观测角度，以物镜垂直载物台为原点可以任意左右倾斜0°～90°，并通过上下、旋转以及水平任意角度移动载物台实现对立体标本的多视角观测。在本文实验中将显微镜物镜向左调整倾斜60°，放大倍数设置为20 倍，然后移动载物台采集了36 张电路板上的电容的多视角图像序列。对图像根据载物台移动轨迹和采集时间进行排序，使其成为有序图像序列。采集的图像如图6 所示。

图6 电容的多视角图像Fig.6 Multi-view image of capacitance

随机从拍摄的36 张图像中选取2 张图像进行图像匹配实验，选取的图片如图7 所示。将选取的2 张图像经过本文的图像灰度化补偿算法进行预处理，然后与Harris 和SIFT 匹配算法经行实验比较，实验数据如表1 所示，三种匹配算法的实验结果图如图8 所示。

图7 选取的待匹配图像Fig.7 Selected image to be matched

表1 不同匹配算法实验数据Tab.1 Experimental data of different matching algorithms

可以发现，尽管其他算法能提取出较多的特征点，但本文算法在针对显微图像中弱纹理区域的匹配上具有明显优势。如图8 所示，本文算法 在对电容壁和电容顶部弱纹理区域的特征点匹配中的误匹配较少，而其他算法存在明显的误匹配现象。如表1，本算法也具有更好的性能指标，原因是本文算法在用Harris 提取出特征角点之后又利用熵值法对特征点进一步筛选，减少了后续特征点匹配过程的计算量，减少了运算时间；筛选出的包含丰富图像特征信息的点，以及在匹配环节采用马氏距离代替欧氏距离，可更好地计算弱纹理区域特征点之间的相似性，提高了匹配的准确率。匹配准确率相比SIFT 算法提高了3.29%，时间比Harris 算法减少了23.42%。

图8 三种匹配算法实验结果Fig.8 Experimental results of three matching algorithms

本文在稀疏重建环节改进了特征提取与匹配，将使用无序图像进行增量SfM 重建改为使用有序图像重建。对图6 所示的36 张图像进行稀疏重建并与SfM 算法和OpenMVG 算法进行实验比较，实验结果如图9 所示，实验数据如表2所示。

表2 SfM 稀疏重建实验数据Tab.2 Experimental data of SfM sparse reconstruction

图9 稀疏点云重建Fig.9 Sparse point cloud reconstruction

可以看到改进后的SfM 重建结果具有更好的结构完整性，本文算法在电容壁处的弱纹理区域重建后的点云比较密集，本文算法相比SfM 算法整体点云数量增加16.69%，相比OpenMVG算法重建时间节省了10.29%。

在使用本文深度预测网络模型进行深度预测前先对模型进行训练。本文使用基恩士-VHX6000 超景深显微镜拍摄的宝石、电子元器件等不同颜色、形状、材质的立体标本的多视角显微图像构建的本地图像数据集，随机从中取出部分图像构成训练集和验证集，对深度预测网络进行监督训练。基恩士-VHX6000 超景深显微镜具有3D 测量功能可以同步扫描立体标本的深度信息生成深度图，即在拍摄立体标本多视角图像时可以同步扫描同一视角的深度图作为模型训练的输入图像，图像对应的深度图作为训练标签。 模型的前$ 部分即卷积层部分使用ResNet50 在 公 开 数 据 集ImageNet［25］上 预 训 练 好的权重参数作为初始化参数，模型训练时用AdaDelta 优化训练效率［26］，通过批处理的随机梯度下降法优化网络参数，本文训练中batch size 设置为4，迭代过程中采用了“早停”机制，即在迭代过程中观察每次迭代输出预测值与真实值的误差，当在验证集上预测值与真实值的误差不再明显变化时停止训练，此时损失函数L在0.217 附近小幅波动、没有明显变化，说明模型收敛，完成本文模型的训练。

图10 展示了部分预测结果，从预测的深度图可以看出立体标本“电容器”的具体形状，深度变化均匀，很好地预测了立体标本的深度。但不足之处在于，在场景背景处与“电容器”壁身弱纹理处存在噪声以及部分缺失，这在之后用MVS 深度图进行修正并利用修正后的深度图恢复三维结构时进行了弥补。

图10 深度预测结果Fig.10 Depth prediction result

通过以SfM 稀疏重建阶段获得的物体稀疏点云和估计的场景几何参数为基础，用本文全卷积神经网络方法预测的图像深度和MVS 估计的深度图相融合，以增强显微镜场景下对弱纹理区域的表征效果，重建出物体的稠密点云。将本文的稠密重建方法的重建效果与开源的三维建模系统COLMAP 中的MVS 方法相比较，估计的深度图如图11 所示，重建后的密集点云如图12所示。

图11 深度图Fig.11 Depth map

图12 稠密重建后的点云Fig.12 Densely reconstructed point cloud

在深度图估计方面，COLMAP 中MVS 方法估计出的深度图会受显微镜聚焦位置的影响，与场景中物体真实结构深度信息不一致。本文方法用全卷积神经网络预测的图像深度对MVS 方法估计的深度图进行矫正，得到的深度图更加精确。在稠密点云方面，可以明显看到本文方法比COLMAP 得到的稠密点云具有更好的结构完整性，在图12 中电容周身弱纹理区域的重建结构更加完整，重建得到的点云数量明显多于COLMAP。稠密重建实验数据如表3 所示。

表3 稠密重建实验数据Tab.3 Dense reconstruction experimental data

用均方根误差（VRMSE）、平均相对误差（VREL）和平均log10 误差（VLOG10）这3 种常见的预测评价指标验证本文使用的全卷积神经网络模型预测深度的精确性，表示为：

其 中：yi为 第i个 像 素 点 的 真 实 深 度，di为i个 像素点的深度预测值，N为图像像素点总数。可以发现与COLMAP 的MVS 算法相比，本算法在对包含了弱纹理区域的显微图像深度预测误差评价指标均更小。点云数目增加31.25%，整体重建时间节省了21.16%。

对重建得到的稠密点云进行模型化，使用泊松表面重建方法将点云网格化，用MeshLab 软件对重建后的模型进行显示，实验结果如图13 所示。由于本文算法重建得到点云数量更多、点云更加稠密，最终得到的物体模型结构更加完整，细节效果更优。使用结构相似度（VSSIM）［27］评价指标验证重建模型的结构完整性，表示为：

图13 泊松表面重建结果Fig.13 Poisson surface reconstruction results

由表3 所示，本文算法重建后模型的SSIM提升了39.06%。

为了对本文算法的鲁棒性［28］经行验证，将显微镜物镜以垂直载物台为原点分别倾斜0°、30°、60°，并选取不同材质和形状的物体进行稠密点云重建实验，实验结果如图14 所示。

从图14 中可以看出，在调整显微镜物镜不同的观测角度，以及采用不同材质的物体进行三维重建实验，本文算法都能够重建出结构完整性较好的稠密点云，具有很好的鲁棒性。

图14 鲁棒性测试Fig.14 Robustness test

4 结论

本文针对实体显微镜应用场景提出了一种对立体标本显微图像进行三维重建的算法。首先，将显微镜镜头固定角度，移动载物台拍摄立体标本的多视角图像，根据拍摄时间以及载物台移动坐标使得采集的图像为有序图像序列。然后，通过使用改进的Harris 与SIFT 相结合的特征提取与匹配算法，提升了算法对于立体标本显微图像在弱纹理区域特征信息的提取，提升了增量SfM 算法性能，得到结构完整性更好的稀疏点云。最后，通过使用与深度残差网络相结合的全卷积神经网络对输入的图像进行深度估计和预测，将预测得到的深度信息通过阈值法与MVS深度图相融合，对MVS 深度图进行修正，重建出物体的稠密点云，提升了重建结构完整性和点云数量。实验结果表明，本算法在显微镜应用场景下无需事先对显微镜镜头进行标定，只需固定显微镜镜头移动载物台采集立体标本的多视角图像进行重建。有效地解决了在显微镜观测场景下对立体标本局部放大观测时，存在的弱纹理区域三维重建效果不佳的问题。相比一些文献的方法，本文方法可以在实体显微镜应用场景下重建出物体结构完整性更好、更加稠密的三维点云，且具有很好的鲁棒性。