Colmap中 patchmatch源码

---

惊世骇俗的Colmap PatchMatch源码

---

500多行，简明易懂

先了解PatchMatch的options：

1. max_image_size: 任一维度的最大图像尺寸
2. gpu_index： gpu的下标
3. depth_min，depth_max： 随机采样的深度范围
4. window_radius：用于计算NCC（两图片的归一化互相关），NCC的计算需要用到 win_size. 图像匹配 | NCC 归一化互相关损失 | 代码 + 讲解 - 知乎 (zhihu.com)
5. window_step: 用于计算NCC时的步长
6. sigma_spatial ，sigma_color： 用于NCC的双边权重，即 灰度差和 像素空间距离的权重
7. num_samples： 蒙特卡洛采样数目
8. ncc_sigma：NCC似然函数传播参数
9. min_triangulation_angle： 最小三角测量角
10. incident_angle_sigma： 入射角似然函数的传播参数
11. num_iterations： 座标扫描迭代次数
12. geom_consistency：是否为代价函数添加正则化的几何一致性项。 如果为 true，则 depth_maps 和 normal_maps 不能为 null
13. geom_consistency_max_cost： NCC的 几何一致性到 图像（像素）一致性的权重
14. filter： 是否启用过滤参数
15. filter_min_ncc：像素具有光一致性的最小 NCC 系数。
16. filter_min_triangulation_angle：达到重建稳定的最小角度
17. filter_min_num_consistent：源图像的最小数量必须一致才能使像素不被过滤。
18. filter_geom_consistency_max_cost：能够保证几何一致性的最大重投影误差
19. cache_size： 单位G，代表缓存大小，越大计算越快。如果一致性图太大，那么会很耗费内存大小。
20. write_consistency_graph： 是否写（保存）一致性图

---

在上述的optios之后，再接着看 PatchMatch 的封装类,并进行代码解释

class PatchMatch {
 public:
  struct Problem {
    // 参考图像的下标
    int ref_image_idx = -1;
    // 源图像们 的所有下标
    std::vector<int> src_image_idxs;
	// 输入的图像数据，用于计算灰度一致性
    std::vector<Image>* images = nullptr;
	// 输入的深度图数据，用于计算几何一致性
    std::vector<DepthMap>* depth_maps = nullptr;
    // 输入的法向图数据，用于计算几何一致性
    std::vector<NormalMap>* normal_maps = nullptr;
    // Print the configuration to stdout.
    void Print() const;
  };
// 构造方法，传入 patch match 的options和 problem结构体
  PatchMatch(const PatchMatchOptions& options, const Problem& problem);
  ~PatchMatch();
// 检查 options 和 problem 参数的正确性
// Check the options and the problem for validity.
  void Check() const;
// 执行算法
// Run the patch match algorithm.
  void Run();

  // 一些get方法，返会 Run得到的结果
  // Get the computed values after running the algorithm.
  DepthMap GetDepthMap() const;
  NormalMap GetNormalMap() const;
  ConsistencyGraph GetConsistencyGraph() const;
  Mat<float> GetSelProbMap() const;
  // 加 _ 的私有变量
 private:
  const PatchMatchOptions options_;
  const Problem problem_;
  std::unique_ptr<PatchMatchCuda> patch_match_cuda_;
};

回顾之前的稀疏重建，其实重建的结果中，我们有

1. images
2. camera.txt 也就是相机的参数
3. points3D 也就是稀疏的3D点

那么到了现在的PatchMatch部分的话， 文件夹应为如下的结构：

images/*
sparse/{cameras.txt, images.txt, points3D.txt}
stereo/
depth_maps/*
normal_maps/*
consistency_graphs/*
patch-match.cfg

前几个似乎就是文件夹用来存放结果的，最后一个文件是什么？
patch-match.cfg 是配置文件， 两行连续指定用于计算一个补丁匹配问题的图像。 第一行指定参考图像，第二行指定源图像。 图像名称相对于images 目录。 在本例中，第一幅参考图像使用所有其他图像作为源图像，第二幅参考图像使用连接最多的 20 张图像作为源图像，第三幅参考图像使用第一幅和第二幅作为源图像。 请注意，必须在 sparse 文件夹中提供的 COLMAP 重建中重建所有指定的图像。 具有如下的形式:

    image_name1.jpg
    __all__
    image_name2.jpg
    __auto__, 20
    image_name3.jpg
    image_name1.jpg, image_name2.jpg

接着是处理PatchMatch任务的这个Thread的定义：

class PatchMatchController : public Thread {
 public:
    // 处理 PatchMatch 任务的 Thread
  PatchMatchController(const PatchMatchOptions& options,
                       const std::string& workspace_path,
                       const std::string& workspace_format,
                       const std::string& pmvs_option_name);

 private:
    // 线程的方法
  void Run();
  void ReadWorkspace();
  void ReadProblems();
  void ReadGpuIndices();
  void ProcessProblem(const PatchMatchOptions& options,
                      const size_t problem_idx);
	//
  const PatchMatchOptions options_;
  const std::string workspace_path_;
  const std::string workspace_format_;
  const std::string pmvs_option_name_;
	// 锁及一些变量
  std::unique_ptr<ThreadPool> thread_pool_;
  std::mutex workspace_mutex_;
  std::unique_ptr<Workspace> workspace_;
  std::vector<PatchMatch::Problem> problems_;
  std::vector<int> gpu_indices_;
  std::vector<std::pair<float, float>> depth_ranges_;
};

说完了方法接口的定义，让我们来看真正的实现把!

以下代码均来自 patch_match.cc 文件

// 一个Check函数
void PatchMatch::Check() const {
  CHECK(options_.Check());
	//check是否支持GPU
  CHECK(!options_.gpu_index.empty());
  const std::vector<int> gpu_indices = CSVToVector<int>(options_.gpu_index);
  //check 检查GPU下标
  CHECK_EQ(gpu_indices.size(), 1);
  CHECK_GE(gpu_indices[0], -1);
  // check 保存的图片非空
  CHECK_NOTNULL(problem_.images);
  // check几何一致性
  if (options_.geom_consistency) {
    // check 深度图，法向图，不为空
      CHECK_NOTNULL(problem_.depth_maps);
    CHECK_NOTNULL(problem_.normal_maps);
    CHECK_EQ(problem_.depth_maps->size(), problem_.images->size());
    CHECK_EQ(problem_.normal_maps->size(), problem_.images->size());
  }

  CHECK_GT(problem_.src_image_idxs.size(), 0);
	// check图像不能有重复
  // Check that there are no duplicate images and that the reference image
  // is not defined as a source image.
  std::set<int> unique_image_idxs(problem_.src_image_idxs.begin(),
                                  problem_.src_image_idxs.end());
  unique_image_idxs.insert(problem_.ref_image_idx);
  CHECK_EQ(problem_.src_image_idxs.size() + 1, unique_image_idxs.size());
	// check 输入的图片是合法的
  // Check that input data is well-formed.
  for (const int image_idx : unique_image_idxs) {
    CHECK_GE(image_idx, 0) << image_idx;
    CHECK_LT(image_idx, problem_.images->size()) << image_idx;
    // 检查图片！
    const Image& image = problem_.images->at(image_idx);
    CHECK_GT(image.GetBitmap().Width(), 0) << image_idx;
    CHECK_GT(image.GetBitmap().Height(), 0) << image_idx;
    CHECK(image.GetBitmap().IsGrey()) << image_idx;
    CHECK_EQ(image.GetWidth(), image.GetBitmap().Width()) << image_idx;
    CHECK_EQ(image.GetHeight(), image.GetBitmap().Height()) << image_idx;

    // Make sure, the calibration matrix only contains fx, fy, cx, cy.
    CHECK_LT(std::abs(image.GetK()[1] - 0.0f), 1e-6f) << image_idx;
    CHECK_LT(std::abs(image.GetK()[3] - 0.0f), 1e-6f) << image_idx;
    CHECK_LT(std::abs(image.GetK()[6] - 0.0f), 1e-6f) << image_idx;
    CHECK_LT(std::abs(image.GetK()[7] - 0.0f), 1e-6f) << image_idx;
    CHECK_LT(std::abs(image.GetK()[8] - 1.0f), 1e-6f) << image_idx;

    if (options_.geom_consistency) {
      CHECK_LT(image_idx, problem_.depth_maps->size()) << image_idx;
      const DepthMap& depth_map = problem_.depth_maps->at(image_idx);
      CHECK_EQ(image.GetWidth(), depth_map.GetWidth()) << image_idx;
      CHECK_EQ(image.GetHeight(), depth_map.GetHeight()) << image_idx;
    }
  }
// check 参考图像数目要和 法向图像数目一致
  if (options_.geom_consistency) {
    const Image& ref_image = problem_.images->at(problem_.ref_image_idx);
    const NormalMap& ref_normal_map =
        problem_.normal_maps->at(problem_.ref_image_idx);
    CHECK_EQ(ref_image.GetWidth(), ref_normal_map.GetWidth());
    CHECK_EQ(ref_image.GetHeight(), ref_normal_map.GetHeight());
  }
}

Run 函数

void PatchMatch::Run() {
  PrintHeading2("PatchMatch::Run");

  Check();

  patch_match_cuda_.reset(new PatchMatchCuda(options_, problem_));
  patch_match_cuda_->Run();
}

• 可以看到patch_match.cc 依赖于 patch_match_cuda 进行匹配，所以我们还是先看一下 patch_match cuda吧！

接下来逐行解释Cuda的代码(遇到一些事情，这个后续马上更新)

#ifndef COLMAP_SRC_MVS_PATCH_MATCH_CUDA_H_
#define COLMAP_SRC_MVS_PATCH_MATCH_CUDA_H_

#include <iostream>
#include <memory>
#include <vector>

#include <cuda_runtime.h>

#include "mvs/cuda_array_wrapper.h"
#include "mvs/depth_map.h"
#include "mvs/gpu_mat.h"
#include "mvs/gpu_mat_prng.h"
#include "mvs/gpu_mat_ref_image.h"
#include "mvs/image.h"
#include "mvs/normal_map.h"
#include "mvs/patch_match.h"

namespace colmap {
namespace mvs {

class PatchMatchCuda {
 public:
  PatchMatchCuda(const PatchMatchOptions& options,
                 const PatchMatch::Problem& problem);
  ~PatchMatchCuda();

  void Run();

  DepthMap GetDepthMap() const;
  NormalMap GetNormalMap() const;
  Mat<float> GetSelProbMap() const;
  std::vector<int> GetConsistentImageIdxs() const;

 private:
  template <int kWindowSize, int kWindowStep>
  void RunWithWindowSizeAndStep();

  void ComputeCudaConfig();

  void InitRefImage();
  void InitSourceImages();
  void InitTransforms();
  void InitWorkspaceMemory();

  // Rotate reference image by 90 degrees in counter-clockwise direction.
  void Rotate();

  const PatchMatchOptions options_;
  const PatchMatch::Problem problem_;

  // Dimensions for sweeping from top to bottom, i.e. one thread per column.
  dim3 sweep_block_size_;
  dim3 sweep_grid_size_;
  // Dimensions for element-wise operations, i.e. one thread per pixel.
  dim3 elem_wise_block_size_;
  dim3 elem_wise_grid_size_;

  // Original (not rotated) dimension of reference image.
  size_t ref_width_;
  size_t ref_height_;

  // Rotation of reference image in pi/2. This is equivalent to the number of
  // calls to `rotate` mod 4.
  int rotation_in_half_pi_;

  std::unique_ptr<CudaArrayWrapper<uint8_t>> ref_image_device_;
  std::unique_ptr<CudaArrayWrapper<uint8_t>> src_images_device_;
  std::unique_ptr<CudaArrayWrapper<float>> src_depth_maps_device_;

  // Relative poses from rotated versions of reference image to source images
  // corresponding to _rotationInHalfPi:
  //
  //    [S(1), S(2), S(3), ..., S(n)]
  //
  // where n is the number of source images and:
  //
  //    S(i) = [K_i(0, 0), K_i(0, 2), K_i(1, 1), K_i(1, 2), R_i(:), T_i(:)
  //            C_i(:), P(:), P^-1(:)]
  //
  // where i denotes the index of the source image and K is its calibration.
  // R, T, C, P, P^-1 denote the relative rotation, translation, camera
  // center, projection, and inverse projection from there reference to the
  // i-th source image.
  std::unique_ptr<CudaArrayWrapper<float>> poses_device_[4];

  // Calibration matrix for rotated versions of reference image
  // as {K[0, 0], K[0, 2], K[1, 1], K[1, 2]} corresponding to _rotationInHalfPi.
  float ref_K_host_[4][4];
  float ref_inv_K_host_[4][4];

  // Data for reference image.
  std::unique_ptr<GpuMatRefImage> ref_image_;
  std::unique_ptr<GpuMat<float>> depth_map_;
  std::unique_ptr<GpuMat<float>> normal_map_;
  std::unique_ptr<GpuMat<float>> sel_prob_map_;
  std::unique_ptr<GpuMat<float>> prev_sel_prob_map_;
  std::unique_ptr<GpuMat<float>> cost_map_;
  std::unique_ptr<GpuMatPRNG> rand_state_map_;
  std::unique_ptr<GpuMat<uint8_t>> consistency_mask_;

  // Shared memory is too small to hold local state for each thread,
  // so this is workspace memory in global memory.
  std::unique_ptr<GpuMat<float>> global_workspace_;
};

}  // namespace mvs
}  // namespace colmap

#endif  // COLMAP_SRC_MVS_PATCH_MATCH_CUDA_H_