Colmap使用笔记

Colmap命令行稀疏重建与位姿导出

Colmap概述

Colmap是一个开源的基于图像的稀疏重建工具，目前的sfm的标杆工具箱。不过虽然Colmap是基于C++写的，但是主要使用方式还是通过命令行来使用，C++没有API文档，对应的pycolmap也没有详细文档。

因为一般Colmap是在服务器上运行的，要先按照文档装好CUDA环境，然后需要在对应的CmakeLists.txt中关闭GUI选项，然后编译安装即可。

最好启用GPU，因为在特征提取时用的GPU-SIFT可以显著加速特征提取。

Colmap命令行方式来进行稀疏重建

Colmap稀疏重建分为下面几个步骤：

0. 准备需要重建的图像和相机参数
1. 特征提取：使用SIFT提取特征，可以配置相机参数
2. 特征匹配：可以配置多种匹配方式
3. 建图：可以配置直接稀疏建图或分段稀疏建图
4. 特征点三角化（可选）：如果是分段重建那么往往在建图后三角化
5. 全局BA（可选）：如果是分段重建那么往往在三角化后进行全局BA

特征提取

特征提取非常简单：

colmap feature_extractor     \
--database_path ./database.db     \
--image_path ./frames     \
--ImageReader.camera_model OPENCV     \
--ImageReader.single_camera 1 \
--ImageReader.camera_params "804.5569,803.4388,967.5456,544.1265,-0.004259,-0.000188,0.000657,0.001718"     

因为一般相机都是通过OPENCV标定的，所以直接在Colmap里设置成OPENCV即可；然后需要设置相机参数，这个参数是相机内参，分别是fx,fy,cx,cy,k1,k2,p1,p2，这些内参直接用字符串的方式传进去即可。也就是设置--ImageReader.camera_params "fx,fy,cx,cy,k1,k2,p1,p2"和--ImageReader.camera_model OPENCV。

而因为是一个相机录制一段视频，所以设置--ImageReader.single_camera 1即可。

特征匹配

特征匹配有多种算法，这里列出几个常用的，分别是exhaustive_matcher、sequential_matcher和vocab_tree_matcher。

实际测试下来2000张左右的图像，直接用exhaustive_matcher就行，大概10分钟就能完成匹配。

exhaustive_matcher

暴力匹配，也就是把每个图像与其他图像进行匹配，速度比较慢，但是对于几百张图片效果很好。

colmap exhaustive_matcher  \
--database_path ./database.db 

sequential_matcher

序列匹配，针对顺序采集的视频图像，由于相邻帧存在视觉上的重叠且没有必要进行完全匹配，它只匹配视频流中的相邻帧。同时，这种匹配方式能够基于vocabulary tree进行回环检测。

colmap sequential_matcher  \
--database_path ./database.db 

vocab_tree_matcher

针对大量图像（几千帧量级），可以通过提供vocabulary tree从而快速检索视觉上最相近的图像进行匹配。

colmap vocab_tree_matcher  \
--database_path ./database.db 

建图

建图有两种方式，一种是直接建图，一种是分段建图。

分段建图是有一个bug的，如果你的帧太少，colmap总共只分了一个段出来，那么在merge的时候会报错，这时候需要修改对应的colmap的源码来解决。

所以少量图像直接建图即可。

通常都是使用直接建图，除非你的图像特别多才会使用分段建图。这里给出几个参考（可能由于配置不同而存在差异）：

1. 300张图像左右，直接建图大概耗时10分钟
2. 600张图像左右，直接建图大概耗时30分钟
3. 1400张图像左右，直接建图大概耗时540分钟

直接建图

colmap mapper \
--database_path ./database.db \
--image_path ./frames \
--output_path ./sparse \
# --Mapper.ba_refine_principal_point true

最后这个--Mapper.ba_refine_principal_point是用来优化主点的，可以按需添加。

分段建图

实际测试下来，用sequential_matcher和hierarchical_mapper来对视频进行重建效果不咋地，可能是拍摄的时候运动模糊太大了导致的。

并且官方文档上推荐在执行完分段建图后最好还额外执行几轮三角化和全局BA。

colmap hierarchical_mapper \
--database_path ./database.db \
--image_path ./frames \
--output_path ./sparse 

特征点三角化

colmap point_triangulator \
--database_path ./database.db \
--image_path ./frames/ \
--input_path ./sparse \
--output_path ./triangulator

全局BA

colmap bundle_adjuster \
--input_path ./triangulator \
--output_path  ./ba

稀疏重建后位姿导出与使用

上面已经完成了稀疏重建，可以在sparse/x找到对应的bin文件，这个x是一个数字，代表的是重建的序号，一般是0，但是如果一次重建效果很差，colmap会自动重试最多三次。

这个bin文件是一个二进制文件，里面存储了重建的信息，可以通过下面的命令来导出：

colmap model_converter 、
--input_path ./sparse/1 \
--output_path ./output \
--output_type TXT

这样就能得到三个txt文件，其中保存相机位姿和点云的文件为images.txt。这个文件从第五行开始，奇数行代表对应的帧和位姿，偶数行表示上一个奇数行的帧观测到的点云。

不过注意，这里的位姿并不能直接是相机位姿，而是世界坐标系到相机的位姿，所以需要进行变换才能使用。

也就是

下面给出一个简单的python代码来读取这个文件并使用open3d来绘制相机轨迹图

import time
import numpy as np
import open3d as o3d


class DrawOpen3d:
    def __init__(self) -> None:

        self.frame = o3d.geometry.TriangleMesh.create_coordinate_frame(size=0.6)
        self.vis = o3d.visualization.Visualizer()
        self.vis.create_window()

        self.CAM_POINTS = np.array(
            [
                [0, 0, 0],
                [-1, -1, 1.5],
                [1, -1, 1.5],
                [1, 1, 1.5],
                [-1, 1, 1.5],
                [-0.5, 1, 1.5],
                [0.5, 1, 1.5],
                [0, 1.2, 1.5],
            ]
        )
        self.CAM_LINES = np.array(
            [
                [1, 2],
                [2, 3],
                [3, 4],
                [4, 1],
                [1, 0],
                [0, 2],
                [3, 0],
                [0, 4],
                [5, 7],
                [7, 6],
            ]
        )

    def quat2rotm(self, q: np.array) -> np.array:
        q2, q3, q4, q1 = q
        R = np.array(
            [
                [
                    q1 * q1 + q2 * q2 - q3 * q3 - q4 * q4,
                    2 * q2 * q3 - 2 * q1 * q4,
                    2 * q2 * q4 + 2 * q1 * q3,
                ],
                [
                    2 * q2 * q3 + 2 * q1 * q4,
                    q1 * q1 - q2 * q2 + q3 * q3 - q4 * q4,
                    2 * q3 * q4 - 2 * q1 * q2,
                ],
                [
                    2 * q2 * q4 - 2 * q1 * q3,
                    2 * q3 * q4 + 2 * q1 * q2,
                    q1 * q1 - q2 * q2 - q3 * q3 + q4 * q4,
                ],
            ]
        )
        return R

    def create_camera_actor(self, g, scale=0.05):
        camera_actor = o3d.geometry.LineSet(
            points=o3d.utility.Vector3dVector(scale * self.CAM_POINTS),
            lines=o3d.utility.Vector2iVector(self.CAM_LINES),
        )

        color = (g * 1.0, 0.5 * (1 - g), 0.9 * (1 - g))
        camera_actor.paint_uniform_color(color)
        return camera_actor

    def parse_colmap_image_txt(self, file):
        data_l = []

        with open(file, "r") as f:
            f.readline()
            f.readline()
            f.readline()
            f.readline()

            while line := f.readline():
                f.readline()
                frame, pose = line.strip("\n").split(" ", 1)
                frame = pose.split("_")[-1].split(".")[0]
                data_l.append((int(frame), pose))

        data_l = sorted(data_l, key=lambda x: x[0])
        rot = []
        trans = []

        for frame, pose in data_l:
            pose = pose.split(" ")
            qw = float(pose[0])
            qx = float(pose[1])
            qy = float(pose[2])
            qz = float(pose[3])

            tx = float(pose[4])
            ty = float(pose[5])
            tz = float(pose[6])

            trans.append([tx, ty, tz])
            rot.append([qx, qy, qz, qw])
        return np.array(rot), np.array(trans)

    def draw_file(self, file):
        count = 0
        rot, trans = self.parse_colmap_image_txt(file)
        pose = []
        for roti, transi in zip(rot, trans):
            T_base = np.eye(4)
            rot_m = self.quat2rotm(roti)
            T_base[:3, :3] = rot_m.T
            T_base[:3, 3] = -rot_m.T @ transi
            pose.append(T_base)

        for i in range(len(pose)):
            T = pose[i]
            cam_actor = self.create_camera_actor(0)
            cam_actor.transform(T)
            self.vis.add_geometry(cam_actor)
        self.vis.poll_events()
        self.vis.update_renderer()
        time.sleep(0.05)
        print(count)

    def run(self):
        self.vis.run()

file = f"images.txt"
draw = DrawOpen3d()
draw.draw_file(file)
draw.run()

注意，Colmap生成的IMAGE_ID不一定是真实的顺序，所以如果是视频的话最好按照拆帧时的序号来对位姿进行排序。

比如，一个视频拆分成形如frame_052.jpg的格式，就如上述代码中的

frame, pose = line.strip("\n").split(" ", 1)
frame = pose.split("_")[-1].split(".")[0]
data_l.append((int(frame), pose))

来对所有位姿进行重排

这样就能得到colmap重建后的相机轨迹图了。