# Original code from Stefan Toeberg (2023)
# Copyright (c) 2023, Stefan Toeberg.
#
# Modified by Maxime Paschoud (2025)
# Copyright (c) 2025, Maxime Paschoud.
# All rights reserved.
#
# This source code is licensed under the BSD-style license found in the
# LICENSE file in the root directory of this source tree.
#
# (http://opensource.org/licenses/BSD-3-Clause)
#
# __author__ = "Stefan Toeberg, LUH: IMR; Maxime Paschoud, ETHZ"
#
from pipeline import DAGNode
import numpy as np
import open3d as o3d
import os
[docs]
class Reconstruction(DAGNode):
"""
Reconstructs 3D point clouds from matched 2D image points using camera intrinsics and extrinsics.
Attributes
----------
pcds : list of list of np.ndarray or None
A nested list of reconstructed point clouds.
- Outer list corresponds to each face.
- Inner list contains one point cloud per image pair.
by_pair : bool
Determines the format of input matches. If True, matches are expected as pairs.
"""
def __init__(self, by_pair=False):
"""
Initializes the reconstruction pipeline.
Parameters
----------
by_pair : bool, optional
If True, expects dense_matches to be formatted as pairs of point sets.
If False, expects a flat list of match tuples. Default is False.
"""
super().__init__()
self.pcds = None
self.by_pair = by_pair
[docs]
def eval(self, camera_position, matches):
"""
Computes the 3D point cloud for each image pair using triangulation.
Parameters
----------
camera_position : list of dict
A list of dictionaries containing camera rotation matrices and intrinsics.
Each dictionary represents one face and must contain:
- 'rotations': np.ndarray of shape (N, 3, 3)
- 'intrinsics': list with two elements:
- intrinsics[0]: list of scale factors (k1, k2, ...)
- intrinsics[1]: list containing [c, s] values (affine scale & shear)
matches : list
2D point matches between images. Format depends on `by_pair`.
Returns
-------
list of list of np.ndarray
Reconstructed 3D point clouds per face and image pair.
"""
pcds = self.compute_point_cloud(camera_position, matches)
self.pcds = pcds
return [pcds]
[docs]
def compute_point_cloud(self, camera_position, dense_matches):
"""
Computes 3D point clouds via affine triangulation from dense 2D matches.
Parameters
----------
camera_position : list of dict
List of camera parameters for each face.
dense_matches : list
List of matched 2D point arrays. Format depends on `by_pair`.
Returns
-------
list of list of np.ndarray
Nested list of point clouds. One list per face, one cloud per image pair.
"""
# iterate over all entries of self.config['selectedPairs']
dense_pcds = []
for j, matches_face in enumerate(dense_matches):
dense_pcds.append([])
l = len(matches_face)
if self.by_pair:
l -= 1
for i in range(l):
idx = (i, i + 1)
# intrinsics
k1 = camera_position[j]['intrinsics'][0][idx[0]]
k2 = camera_position[j]['intrinsics'][0][idx[1]]
c = camera_position[j]['intrinsics'][1][0]
s = camera_position[j]['intrinsics'][1][1]
# get truncated rotation matrices
rot1 = camera_position[j]['rotations'][idx[0]][:2]
rot2 = camera_position[j]['rotations'][idx[1]][:2]
# build affine camera matrices
P1 = k1*np.array([[c, 0], [s, 1]]).dot(rot1)
P2 = k2*np.array([[c, 0], [s, 1]]).dot(rot2)
# center matches / no need to further consider the translations
if self.by_pair:
centroid1 = np.mean(matches_face[idx[0]], axis=0)[:, None]
centroid2 = np.mean(matches_face[idx[1]], axis=0)[:, None]
matches1 = matches_face[idx[0]].T - centroid1
matches2 = matches_face[idx[1]].T - centroid2
else:
centroid1 = np.mean(matches_face[i][0], axis=0)[:, None]
centroid2 = np.mean(matches_face[i][1], axis=0)[:, None]
matches1 = matches_face[i][0].T - centroid1
matches2 = matches_face[i][1].T - centroid2
# perform linear triangulation (inhomogeneous equation system)
matches12 = np.vstack((matches1, matches2))
P12 = np.vstack((P1, P2))
pcd = self.inhom_triangulation(matches12, P12)
dense_pcds[-1].append(pcd)
return dense_pcds
[docs]
def inhom_triangulation(self, b, A):
"""
Solves an inhomogeneous linear system `Ax = b` using the Moore-Penrose pseudo-inverse.
Parameters
----------
b : np.ndarray
Right-hand side matrix (stacked 2D correspondences).
A : np.ndarray
Projection matrix (camera model matrix).
Returns
-------
np.ndarray
Triangulated 3D point cloud of shape (N, 3).
"""
A_pinv = np.linalg.pinv(A)
x = A_pinv.dot(b)
return x.T
[docs]
def save(self, path):
"""
Saves the reconstructed 3D point clouds to a `.npy` file.
Parameters
----------
path : str
Path to save the numpy file. Must end with `.npy`.
Raises
------
ValueError
If path does not end in `.npy` or if no point clouds exist.
"""
if not path.endswith('.npy'):
raise ValueError(f"path should end with an .npy extension. path is {path}.")
if self.pcds is None:
raise ValueError(f"self.pcds is None. self.pcds is set running self.eval().")
np.save(path, self.pcds)
[docs]
def show(self, face=0):
"""
Displays the reconstructed point cloud for a given face using Open3D.
Parameters
----------
face : int, optional
Index of the face to display (default: 0).
"""
pcd_face = self.pcds[face]
pcd_list = []
for p in pcd_face:
pcd_list.append(o3d.geometry.PointCloud(o3d.utility.Vector3dVector(np.asarray(p))))
o3d.visualization.draw_geometries(pcd_list)