use the id mapping from imagecols to do colmap point triangulation.

cfgs/fitnmerge/default.yaml

@@ -24,7 +24,6 @@ sfm:
     ranges:
         range_robust: [0.05, 0.95]
         k_stretch: 1.25
-    fbase: "sift" # ["sift", "hloc"]
     hloc:
         descriptor: "superpoint_aachen"
         matcher: "NN-superpoint"

cfgs/fitnmerge/hypersim.yaml

@@ -10,6 +10,4 @@ input_stride: 1
 line2d:
     detector: 
         method: "sold2"
-sfm:
-    fbase: "hloc"
 

cfgs/fitnmerge/scannet.yaml

@@ -8,6 +8,4 @@ stride: 5
 line2d:
     detector: 
         method: "sold2"
-sfm:
-    fbase: "hloc"
 

cfgs/triangulation/default.yaml

@@ -28,7 +28,6 @@ sfm:
     ranges:
         range_robust: [0.05, 0.95]
         k_stretch: 1.25
-    fbase: "sift" # ["sift", "hloc"]
     hloc:
         descriptor: "superpoint_aachen"
         matcher: "NN-superpoint"

cfgs/triangulation/hypersim.yaml

@@ -10,6 +10,4 @@ input_stride: 1
 line2d:
     detector: 
         method: "sold2"
-sfm:
-    fbase: "hloc"
 

cfgs/triangulation/scannet.yaml

@@ -8,6 +8,4 @@ stride: 5
 line2d:
     detector: 
         method: "sold2"
-sfm:
-    fbase: "hloc"
 

limap/base/bindings.cc

@@ -346,6 +346,7 @@ void bind_camera(py::module& m) {
         .def(py::init<const std::string&, M3D, int, std::pair<int, int>>(), py::arg("model_name"), py::arg("K"), py::arg("cam_id")=-1, py::arg("hw")=std::make_pair<int, int>(-1, -1))
         .def(py::init<py::dict>())
         .def(py::init<const Camera&>())
+        .def(py::init<int, int, std::pair<int, int>>(), py::arg("model_id"), py::arg("cam_id"), py::arg("hw")=std::make_pair<int, int>(-1, -1)) // empty camera
         .def(py::pickle(
             [](const Camera& input) { // dump
                 return input.as_dict();
@@ -396,6 +397,7 @@ void bind_camera(py::module& m) {
         .def(py::init<const Camera&, const CameraPose&, const std::string&>(), py::arg("camera"), py::arg("pose"), py::arg("image_name") = "none")
         .def(py::init<py::dict>())
         .def(py::init<const CameraImage&>())
+        .def(py::init<const int&, const std::string&>(), py::arg("cam_id"), py::arg("image_name") = "none") // empty image
         .def(py::pickle(
             [](const CameraImage& input) { // dump
                 return input.as_dict();

limap/base/camera.cc

@@ -10,6 +10,14 @@ Camera::Camera(const colmap::Camera& cam) {
     SetWidth(cam.Width());
 }
 
+// empty camera
+Camera::Camera(int model_id, int cam_id, std::pair<int, int> hw) {
+    SetModelId(model_id);
+    SetCameraId(cam_id);
+    SetHeight(hw.first);
+    SetWidth(hw.second);
+}
+
 Camera::Camera(int model_id, const std::vector<double>& params, int cam_id, std::pair<int, int> hw) {
     SetModelId(model_id);
     THROW_CHECK_EQ(params.size(), NumParams());

limap/base/camera.h

@@ -44,6 +44,7 @@ class Camera: public colmap::Camera {
     Camera(const std::string& model_name, M3D K, int cam_id=-1, std::pair<int, int> hw=std::make_pair<int, int>(-1, -1));
     Camera(py::dict dict);
     Camera(const Camera& cam);
+    Camera(int model_id, int cam_id, std::pair<int, int> hw=std::make_pair<int, int>(-1, -1)); // empty camera
     bool operator ==(const Camera&);
 
     py::dict as_dict() const;
@@ -66,8 +67,8 @@ class CameraPose {
     CameraPose(py::dict dict);
     CameraPose(const CameraPose& campose): qvec(campose.qvec), tvec(campose.tvec) {}
 
-    V4D qvec;
-    V3D tvec;
+    V4D qvec = V4D::Zero();
+    V3D tvec = V3D::Zero();
 
     py::dict as_dict() const;
     M3D R() const { return colmap::QuaternionToRotationMatrix(qvec); }

limap/base/camera_view.h

@@ -24,6 +24,7 @@ class CameraImage {
     CameraImage(const Camera& input_cam, const CameraPose& input_pose, const std::string& image_name = "none"): cam_id(input_cam.CameraId()), pose(input_pose), image_name_(image_name) {}
     CameraImage(py::dict dict);
     CameraImage(const CameraImage& camimage): cam_id(camimage.cam_id), pose(camimage.pose) {SetImageName(camimage.image_name());}
+    CameraImage(const int& input_cam_id, const std::string& image_name = "none"): cam_id(input_cam_id), image_name_(image_name) {} // empty image
 
     int cam_id;
     CameraPose pose;

limap/pointsfm/bindings.cc

@@ -29,13 +29,13 @@ void bind_pointsfm(py::module &m) {
     // bind the new sfm model.
     py::class_<SfmModel>(m, "SfmModel")
         .def(py::init<>())
-        .def("addImage", &SfmModel::addImage)
+        .def("addImage", &SfmModel::addImage, py::arg("image"), py::arg("img_id") = -1)
         .def("addPoint", &SfmModel::addPoint)
         .def("ReadFromCOLMAP", &SfmModel::ReadFromCOLMAP)
         .def("GetImageNames", &SfmModel::GetImageNames)
         .def("ComputeNumPoints", &SfmModel::ComputeNumPoints)
         .def("ComputeSharedPoints", &SfmModel::ComputeSharedPoints)
-        .def("GetMaxOverlappingImages", &SfmModel::GetMaxOverlappingImages)
+        .def("GetMaxOverlapImages", &SfmModel::GetMaxOverlapImages)
         .def("GetMaxIoUImages", &SfmModel::GetMaxIoUImages)
         .def("GetMaxDiceCoeffImages", &SfmModel::GetMaxDiceCoeffImages)
         .def("ComputeNumPoints", &SfmModel::ComputeNumPoints)

limap/pointsfm/colmap_sfm.py

@@ -25,7 +25,31 @@ def run_colmap_sift_matches(image_path, db_path, use_cuda=False):
            '--SiftMatching.use_gpu', str(int(use_cuda))]
     subprocess.run(cmd, check=True)
 
-def run_hloc_matches(cfg, image_path, db_path, keypoints=None):
+def import_images_with_known_cameras(image_dir, database_path, imagecols):
+    if imagecols is None:
+        return reconstruction.import_images(image_dir=image_path, database_path=db_path)
+    image_ids = imagecols.get_img_ids()
+    image_name_list = sorted(os.listdir(image_dir))
+    assert len(image_name_list) == len(image_ids)
+
+    # connect to the database
+    db = database.COLMAPDatabase(database_path)
+    db.create_tables()
+    # add camera
+    for cam_id in imagecols.get_cam_ids():
+        cam = imagecols.cam(cam_id)
+        import pdb
+        db.add_camera(cam.model_id(), cam.w(), cam.h(), cam.params(), camera_id=cam_id)
+    # add image
+    for img_name, img_id in zip(image_name_list, image_ids):
+        cam_id = imagecols.camimage(img_id).cam_id
+        db.add_image(img_name, cam_id, image_id=img_id)
+    db.commit()
+
+def run_hloc_matches(cfg, image_path, db_path, keypoints=None, imagecols=None):
+    '''
+    use the id mapping from _base.ImageCollection to do the match
+    '''
     image_path = Path(image_path)
     from hloc import extract_features, match_features, reconstruction, triangulation
     outputs = Path(os.path.join(os.path.dirname(db_path), 'hloc_outputs'))
@@ -41,17 +65,20 @@ def run_hloc_matches(cfg, image_path, db_path, keypoints=None):
 
     sfm_dir.mkdir(parents=True, exist_ok=True)
     reconstruction.create_empty_db(db_path)
-    reconstruction.import_images(image_dir=image_path, database_path=db_path)
+    if imagecols is None:
+        reconstruction.import_images(image_dir=image_path, database_path=db_path)
+    else:
+        import_images_with_known_cameras(image_path, db_path, imagecols) # use cameras and id mapping
     image_ids = reconstruction.get_image_ids(db_path)
     reconstruction.import_features(image_ids, db_path, feature_path)
     reconstruction.import_matches(image_ids, db_path, sfm_pairs, match_path, None, None)
     triangulation.estimation_and_geometric_verification(db_path, sfm_pairs)
 
-def run_colmap_sfm_with_known_poses(cfg, imagecols, output_path='tmp/tmp_colmap', keypoints=None, use_cuda=False, skip_exists=False, map_to_original_image_names=False):
+def run_colmap_sfm_with_known_poses(cfg, imagecols, output_path='tmp/tmp_colmap', keypoints=None, skip_exists=False, map_to_original_image_names=True):
     ### set up path
     db_path = os.path.join(output_path, 'db.db')
     image_path = os.path.join(output_path, 'images')
-    model_path = os.path.join(output_path, 'sparse', 'model')
+    model_path = os.path.join(output_path, 'sparse', 'reference_model')
     point_triangulation_path = os.path.join(output_path, 'sparse')
 
     ### initialize sparse folder
@@ -78,11 +105,7 @@ def run_colmap_sfm_with_known_poses(cfg, imagecols, output_path='tmp/tmp_colmap'
             keypoints_in_order.append(keypoints[img_id])
 
     # sift feature extraction and matching
-    if cfg["fbase"] == "sift":
-        run_colmap_sift_matches(image_path, db_path, use_cuda=use_cuda)
-    elif cfg["fbase"] == "hloc":
-        assert (use_cuda == True)
-        run_hloc_matches(cfg["hloc"], image_path, Path(db_path), keypoints=keypoints_in_order)
+    run_hloc_matches(cfg["hloc"], image_path, Path(db_path), keypoints=keypoints_in_order, imagecols=imagecols)
 
     ### write cameras.txt
     colmap_cameras = {}
@@ -105,18 +128,17 @@ def run_colmap_sfm_with_known_poses(cfg, imagecols, output_path='tmp/tmp_colmap'
     db = database.COLMAPDatabase(db_path)
     rows = db.execute("SELECT * from images")
     colmap_images = {}
-    for idx, img_id in enumerate(imagecols.get_img_ids()):
+    for img_id in imagecols.get_img_ids():
         kk = next(rows)
-        assert (kk[0] == idx + 1)
+        assert (kk[0] == img_id)
         imname = kk[1]
-        img_id_orig = int(imname[5:-4])
-        camimage = imagecols.camimage(img_id_orig)
+        camimage = imagecols.camimage(img_id)
         cam_id = camimage.cam_id
         qvec = camimage.pose.qvec
         tvec = camimage.pose.tvec
-        colmap_images[idx+1] = colmap_utils.Image(id=idx+1, qvec=qvec, tvec=tvec,
-                                                  camera_id=cam_id, name=imname,
-                                                  xys=[], point3D_ids=[])
+        colmap_images[img_id] = colmap_utils.Image(id=img_id, qvec=qvec, tvec=tvec,
+                                                   camera_id=cam_id, name=imname,
+                                                   xys=[], point3D_ids=[])
     fname = os.path.join(model_path, 'images.txt')
     colmap_utils.write_images_text(colmap_images, fname)
 
@@ -137,7 +159,7 @@ def run_colmap_sfm_with_known_poses(cfg, imagecols, output_path='tmp/tmp_colmap'
     if map_to_original_image_names:
         fname_images_bin = os.path.join(point_triangulation_path, "images.bin")
         colmap_images = colmap_utils.read_images_binary(fname_images_bin)
-        for idx, img_id in enumerate(imagecols.get_img_ids()):
-            colmap_images[idx + 1] = colmap_images[idx + 1]._replace(name = imagecols.image_name(img_id))
+        for img_id in imagecols.get_img_ids():
+            colmap_images[img_id] = colmap_images[img_id]._replace(name = imagecols.image_name(img_id))
         colmap_utils.write_images_binary(colmap_images, fname_images_bin)
 

limap/pointsfm/functions.py

@@ -14,44 +14,18 @@ def filter_by_cam_id(cam_id, prev_imagecols, prev_neighbors):
     neighbors = imagecols.update_neighbors(prev_neighbors)
     return imagecols, neighbors
 
-def ComputeNeighborsVector(model, n_neighbors, min_triangulation_angle=1.0, neighbor_type="iou"):
+def ComputeNeighbors(model, n_neighbors, min_triangulation_angle=1.0, neighbor_type="iou"):
     '''
-    Returns: vector<vector<int>>
+    Returns: map<int, vector<int>>
     '''
     if neighbor_type == "iou":
-        neighbors_vec = model.GetMaxIoUImages(n_neighbors, min_triangulation_angle)
+        neighbors = model.GetMaxIoUImages(n_neighbors, min_triangulation_angle)
     elif neighbor_type == "overlap":
-        neighbors_vec = model.GetMaxOverlappingImages(n_neighbors, min_triangulation_angle)
+        neighbors = model.GetMaxOverlapImages(n_neighbors, min_triangulation_angle)
     elif neighbor_type == "dice":
-        neighbors_vec = model.GetMaxDiceCoeffImages(n_neighbors, min_triangulation_angle)
+        neighbors = model.GetMaxDiceCoeffImages(n_neighbors, min_triangulation_angle)
     else:
         raise NotImplementedError
-    return neighbors_vec
-
-def ComputeNeighbors(model, n_neighbors, min_triangulation_angle=1.0, neighbor_type="iou"):
-    '''
-    Returns: map<int, vector<int>>
-    '''
-    neighbors_vec = ComputeNeighborsVector(model, n_neighbors, min_triangulation_angle=min_triangulation_angle, neighbor_type=neighbor_type)
-    neighbors = {}
-    for img_id in range(len(neighbors_vec)):
-        neighbors[img_id] = neighbors_vec[img_id]
-    return neighbors
-
-# compute neighborhood for a image list sorted as 'image{0:08d}.png'
-def ComputeNeighborsSorted(model, n_neighbors, min_triangulation_angle=1.0, neighbor_type="iou"):
-    '''
-    Returns: map<int, vector<int>>
-    '''
-    # get neighbors
-    neighbors_vec = ComputeNeighborsVector(model, n_neighbors, min_triangulation_angle=min_triangulation_angle, neighbor_type=neighbor_type)
-
-    # map indexes
-    image_names = model.GetImageNames()
-    image_id_list = [int(name[5:-4]) for name in image_names]
-    neighbors = {}
-    for idx, img_id in enumerate(image_id_list):
-        neighbors[img_id] = [image_id_list[k] for k in neighbors_vec[idx]]
     return neighbors
 
 def compute_metainfos(cfg, model, n_neighbors=20):

limap/pointsfm/sfm_model.cc

@@ -1,6 +1,9 @@
 #include "pointsfm/sfm_model.h"
+#include "util/log_exceptions.h"
 
 #include <colmap/util/math.h>
+#include <colmap/base/reconstruction.h>
+#include <colmap/util/misc.h>
 
 namespace limap {
 
@@ -25,8 +28,28 @@ void SfmModel::addPoint(double x, double y, double z, const std::vector<int>& im
     points.push_back(p);
 }
 
-void SfmModel::addImage(const colmap::mvs::Image& image) {
+void SfmModel::ReadFromCOLMAP(const std::string& path,
+                              const std::string& sparse_path,
+                              const std::string& images_path) {
+    colmap::mvs::Model::ReadFromCOLMAP(path, sparse_path, images_path);
+
+    // store image ids
+    colmap::Reconstruction reconstruction;
+    reconstruction.Read(colmap::JoinPaths(path, sparse_path));
+    for (size_t i = 0; i < reconstruction.NumRegImages(); ++i) {
+        reg_image_ids.push_back(reconstruction.RegImageIds()[i]);
+    }
+}
+
+void SfmModel::addImage(const colmap::mvs::Image& image, const int img_id) {
     images.push_back(image);
+    if (img_id == -1) {
+        if (!reg_image_ids.empty())
+            THROW_CHECK_EQ(reg_image_ids.back(), reg_image_ids.size() - 1);
+        reg_image_ids.push_back(reg_image_ids.size());
+    }
+    else
+        reg_image_ids.push_back(img_id);
 }
 
 std::vector<int> SfmModel::ComputeNumPoints() const {
@@ -48,8 +71,26 @@ std::vector<std::string> SfmModel::GetImageNames() const {
     return image_names;
 }
 
+std::map<int, std::vector<int>> SfmModel::neighbors_vec_to_map(const std::vector<std::vector<int>>& neighbors) const {
+    std::map<int, std::vector<int>> m_neighbors;
+    for (size_t i = 0; i < reg_image_ids.size(); ++i) {
+        int img_id = reg_image_ids[i];
+        m_neighbors.insert(std::make_pair(img_id, std::vector<int>()));
+        for (auto it = neighbors[i].begin(); it != neighbors[i].end(); ++it) {
+            m_neighbors.at(img_id).push_back(reg_image_ids[*it]);
+        }
+    }
+    return m_neighbors;
+}
+
+std::map<int, std::vector<int>> SfmModel::GetMaxOverlapImages(
+        const size_t num_images, const double min_triangulation_angle) const {
+    std::vector<std::vector<int>> overlapping_images = colmap::mvs::Model::GetMaxOverlappingImages(num_images, min_triangulation_angle);
+    return neighbors_vec_to_map(overlapping_images);
+}
+
 // The original SfmModel::GetMaxOverlappingImages in COLMAP uses the number of intersection rather than IoU
-std::vector<std::vector<int>> SfmModel::GetMaxIoUImages(
+std::map<int, std::vector<int>> SfmModel::GetMaxIoUImages(
     const size_t num_images, const double min_triangulation_angle) const {
   std::vector<std::vector<int>> overlapping_images(images.size());
 
@@ -105,10 +146,10 @@ std::vector<std::vector<int>> SfmModel::GetMaxIoUImages(
     }
   }
 
-  return overlapping_images;
+  return neighbors_vec_to_map(overlapping_images);
 }
 
-std::vector<std::vector<int>> SfmModel::GetMaxDiceCoeffImages(
+std::map<int, std::vector<int>> SfmModel::GetMaxDiceCoeffImages(
     const size_t num_images, const double min_triangulation_angle) const {
   std::vector<std::vector<int>> overlapping_images(images.size());
 
@@ -164,7 +205,7 @@ std::vector<std::vector<int>> SfmModel::GetMaxDiceCoeffImages(
     }
   }
 
-  return overlapping_images;
+  return neighbors_vec_to_map(overlapping_images);
 }
 
 std::pair<float, float> SfmModel::get_robust_range(std::vector<float>& data, const std::pair<double, double>& range_robust, const double& kstretch) const {

limap/pointsfm/sfm_model.h

@@ -27,22 +27,32 @@ class SfmModel: public colmap::mvs::Model {
 
     void addPoint(double x, double y, double z, const std::vector<int>& image_ids);
 
-    void addImage(const colmap::mvs::Image& image);
+    void addImage(const colmap::mvs::Image& image, const int img_id = -1);
+
+    void ReadFromCOLMAP(const std::string& path, 
+                        const std::string& sparse_path = "sparse",
+                        const std::string& images_path = "images");
 
     std::vector<std::string> GetImageNames() const;
 
     std::vector<int> ComputeNumPoints() const;
 
-    std::vector<std::vector<int>> GetMaxIoUImages(
+    std::map<int, std::vector<int>> GetMaxOverlapImages(
+            const size_t num_images, const double min_triangulationo_angle) const;
+
+    std::map<int, std::vector<int>> GetMaxIoUImages(
             const size_t num_images, const double min_triangulationo_angle) const;
 
-    std::vector<std::vector<int>> GetMaxDiceCoeffImages(
+    std::map<int, std::vector<int>> GetMaxDiceCoeffImages(
             const size_t num_images, const double min_triangulationo_angle) const;
 
     std::pair<V3D, V3D> ComputeRanges(const std::pair<double, double>& range_robust, const double& kstretch) const;
 
 private:
     std::pair<float, float> get_robust_range(std::vector<float>& data, const std::pair<double, double>& range_robust, const double& kstretch) const;
+
+    std::vector<int> reg_image_ids;
+    std::map<int, std::vector<int>> neighbors_vec_to_map(const std::vector<std::vector<int>>& neighbors) const;
 };
 
 } // namespace pointsfm