Geometry & 3D Vision

Point Cloud Registration

Estimating the transformation that aligns two or more point clouds of the same object or scene, by local refinement such as ICP, by global methods built on 3D feature correspondences, by globally optimal search, or by probabilistic and non-rigid models, with the benchmarks and metrics used to compare them.

intermediate

Point cloud registration is the task of finding the transformation that brings one point cloud into the coordinate frame of another, so that points sampled from the same surface coincide. The clouds may be two range scans of an object taken from different sides, consecutive sweeps of a LiDAR on a moving vehicle, depth frames of a handheld camera, or a model and a scan of a scene that contains it. In the common rigid case the answer is a rotation and a translation; in the non-rigid case it is a deformation field. The task is also called scan matching, scan registration or point-set registration, and it is the point-cloud counterpart of image registration.

Two families of methods dominate. Local methods, above all iterative closest point (ICP), refine an initial alignment by alternating nearest-neighbour matching with a closed-form fit; they are accurate but only from a good start. Global methods need no initial pose: most match local 3D descriptors such as FPFH between the clouds and fit a transform to the matches robustly, with RANSAC, robust optimization (fast global registration) or certifiable estimators (TEASER); others search the whole space of rotations and translations (Go-ICP). Probabilistic methods such as coherent point drift treat one cloud as a Gaussian mixture and extend to non-rigid registration, and learned descriptors and networks have taken over the correspondence step on benchmarks such as 3DMatch and KITTI. The closed-form alignment of matched points that ICP, RANSAC pipelines and many learned methods call internally is the Kabsch algorithm.

Problem Definition

Given a source cloud P={pi}i=1NP = \{\mathbf{p}_i\}_{i=1}^{N} and a target cloud Q={qj}j=1MQ = \{\mathbf{q}_j\}_{j=1}^{M} that sample overlapping parts of the same surface, rigid registration seeks R∈SO(3)R \in SO(3) and t∈R3\mathbf{t} \in \mathbb{R}^3 such that Rpi+tR\mathbf{p}_i + \mathbf{t} lies on the surface sampled by QQ for every pi\mathbf{p}_i in the overlap. Neither the correspondences nor the overlap are given. A generic objective makes this explicit: with an unknown assignment j(i)j(i) and a robust penalty ρ\rho,

min⁡R, t, j(⋅)  ∑i=1Nρ ⁣(∥Rpi+t−qj(i)∥).\min_{R,\,\mathbf{t},\,j(\cdot)} \; \sum_{i=1}^{N} \rho\!\left( \lVert R\mathbf{p}_i + \mathbf{t} - \mathbf{q}_{j(i)} \rVert \right).

ICP fixes the assignment to nearest neighbours under the current pose and takes ρ(r)=r2\rho(r) = r^2 with pairs beyond a distance rejected. Correspondence-based global methods fix the assignment once, by descriptor matching, and make ρ\rho robust to the wrong matches. Branch-and-bound methods optimize over pose and nearest-neighbour assignment together. When the assignment is known and correct, the problem reduces to the least-squares fit that the Kabsch algorithm solves in closed form; the difficulty of registration lies entirely in the unknown, partly wrong correspondences.

Three assumptions are usually made and are worth stating. The scene is static between the two captures, or the moving parts are treated as outliers. The overlap is large enough to determine the transform: the 3DMatch benchmark scores only pairs with more than 30% overlap, and PREDATOR’s 3DLoMatch benchmark targets pairs with 10–30% (Huang et al., 2021). And the overlapping geometry constrains all six degrees of freedom; a plane, a cylinder or a long corridor leaves some motions undetermined, as the ICP article explains.

Input and Output

  • Input: two point clouds, sometimes with normals, colours or intensities; optionally an initial transform (from odometry, an IMU, the previous frame or a user), and parameters such as a voxel size, a descriptor radius and an inlier distance.
  • Output: a rigid transform (R,t)(R, \mathbf{t}), or a similarity with scale, or a non-rigid displacement for every source point; often the inlier correspondences and a quality score (an inlier fraction or an RMSE) that the caller uses to accept or reject the result.

Variants

  • Rigid, similarity and non-rigid. Most 3D scanning and robotics registration is rigid. A scale factor is needed when one cloud comes from monocular reconstruction; TEASER’s formulation estimates scale, rotation and translation (Yang and Carlone, 2019). Non-rigid registration aligns deforming shapes such as faces; coherent point drift was introduced for it (Myronenko, Song and Carreira-Perpiñán, 2006).
  • Local and global. A local method needs an initial estimate inside its basin of convergence; a global method does not. Zhou, Park and Koltun (2016) describe the usual practical pipeline as global alignment by RANSAC over candidate correspondences followed by local refinement with ICP.
  • Pairwise and multi-way. Reconstructing a scene from many fragments requires poses for all of them. Existing pipelines first compute pairwise alignments and then synchronize them into globally consistent poses; fast global registration also has a multi-way form that optimizes a single objective over all poses directly (Zhou et al., 2016).
  • Scan-to-scan, scan-to-map and model-to-scene. Odometry registers consecutive scans; mapping registers each scan to an accumulated map; object pose estimation registers a model to a cluttered scene, as in the UWA benchmark of objects in scenes used by Zhou et al. (2016) and the cap-localization experiment of Go-ICP (Yang, Li and Jia, 2013).

Challenges

  • Unknown correspondences. Independent scans never sample the same points, so even at the true pose a point’s nearest neighbour is offset along the surface.
  • Partial overlap. Points seen by only one sensor have no counterpart. A squared-error fit is pulled by them, and the fraction of matched points in the overlap falls as the overlap shrinks; PREDATOR was designed for pairs below 30% overlap, where earlier learned pipelines lost much of their registration recall (Huang et al., 2021).
  • Outlier correspondences. Yang and Carlone (2019) note that 3D keypoint matching is less accurate than its 2D counterparts such as SIFT and ORB, that 95% wrong correspondences are considered common, and that the run time of RANSAC grows exponentially with the outlier ratio.
  • Local minima. ICP’s objective has many local minima; ICP converges to one of them and, as Yang et al. (2013) point out, cannot tell whether it is the global one.
  • Noise and sensor artefacts. Range data carry noise and systematic distortion; the scene benchmark used by Zhou et al. (2016) has fragments with high-frequency noise and low-frequency distortion that simulate scans from consumer depth cameras.
  • Symmetry, repetition and degeneracy. Repeated structures produce look-alike descriptors and alternative alignments that explain the data almost as well as the true one, as in feature matching between images; flat or extruded geometry leaves motions undetermined.
  • Scale. Large clouds are downsampled before registration, often with a voxel grid, as PREDATOR and FCGF do; Go-ICP was run on 1,000-point subsets because full-resolution clouds took hours (Zhou et al., 2016).

Approaches

Local registration: ICP and its variants

ICP, as the ICP article describes Besl and McKay’s (1992) statement of it, alternates closest-point matching with a closed-form rigid fit and converges monotonically to a local minimum of the mean squared distance; Chen and Medioni (1992) minimized distances to tangent planes, the basis of point-to-plane ICP. Neither paper could be read for this article. The ICP article covers the iteration, the point-to-point and point-to-plane metrics, Generalized-ICP, the normal distributions transform and projective association in depth. For the task, three properties matter. First, ICP is a refinement method: according to the ICP article, Besl and McKay obtained global registration by running it from a set of initial rotations and translations and keeping the best. Second, its basin of convergence is limited; in Zhou et al.’s (2016) controlled comparison on synthetic range images, the accuracy of the PCL and Sparse ICP variants they tested degraded once the initialization deviated by 5° in rotation or by 5–10% of the cloud’s diameter in translation. Third, partial overlap and outliers are handled by rejecting distant pairs, trimming or robust losses (ICP article), and in the example below the rejection distance trades accuracy against basin width.

Tracking systems run local registration at every frame from the previous estimate: KinectFusion aligns each depth frame to its model starting from the previous pose, and LiDAR SLAM registers successive scans by scan matching.

Global registration from feature correspondences

The standard global pipeline computes a local descriptor at each (downsampled) point, matches descriptors between the clouds, and estimates the transform robustly from the putative matches.

FPFH descriptors. Point Feature Histograms describe the geometry around a point by histograms of relations between the surface normals of point pairs in its neighbourhood (PCL documentation). For points ps,pt\mathbf{p}_s, \mathbf{p}_t with normals ns,nt\mathbf{n}_s, \mathbf{n}_t and d=∥pt−ps∥d = \lVert\mathbf{p}_t - \mathbf{p}_s\rVert, a Darboux frame is set at the source point,

u=ns,v=u×pt−psd,w=u×v,\mathbf{u} = \mathbf{n}_s, \qquad \mathbf{v} = \mathbf{u} \times \frac{\mathbf{p}_t - \mathbf{p}_s}{d}, \qquad \mathbf{w} = \mathbf{u} \times \mathbf{v},

and three angular features are measured (PCL’s tutorial notation; its code also normalizes v\mathbf{v} to unit length and computes it with the opposite sign, which mirrors the α\alpha and θ\theta histograms):

α=v⋅nt,ϕ=u⋅pt−psd,θ=arctan⁡ ⁣(w⋅nt,  u⋅nt).\alpha = \mathbf{v} \cdot \mathbf{n}_t, \qquad \phi = \mathbf{u} \cdot \frac{\mathbf{p}_t - \mathbf{p}_s}{d}, \qquad \theta = \arctan\!\left(\mathbf{w} \cdot \mathbf{n}_t,\; \mathbf{u} \cdot \mathbf{n}_t\right).

They depend only on the relative position and orientation of the two oriented points, so they do not change when the cloud is rotated or translated. PFH bins these features over all pairs within a radius, which costs O(k2)O(k^2) per point for kk neighbours. The Fast Point Feature Histogram of Rusu, Blodow and Beetz (2009) reduces this to O(k)O(k): a simplified histogram (SPFH) of each point is built from its pairs with its own neighbours only, and the FPFH of a query point adds its neighbours’ SPFHs weighted by the inverse of their distance,

FPFH⁡(pq)=SPFH⁡(pq)+1k∑i=1k1ωi SPFH⁡(pi),\operatorname{FPFH}(\mathbf{p}_q) = \operatorname{SPFH}(\mathbf{p}_q) + \frac{1}{k}\sum_{i=1}^{k} \frac{1}{\omega_i}\,\operatorname{SPFH}(\mathbf{p}_i),

which reaches points up to twice the radius away. The three features are binned separately into 11 bins each and concatenated, giving 33 values (PCL 1.15.1 tutorials). The PCL implementation differs from the tutorial’s formula in details: its weighting loop skips the query point itself, so the query’s own SPFH is not added, its weights are the inverse squared distances that PCL’s radius search returns, each 11-bin block is scaled to sum to 100, and the pair-feature code chooses as source the point whose normal makes the smaller angle with the line joining the two points. The ICRA 2009 paper could not be read for this article; the description follows PCL’s documentation and source code. Fast global registration chose FPFH because it can be computed in a fraction of a millisecond and matches well across a broad range of datasets (Zhou et al., 2016); it depends on the quality of the normals, and in FCGF’s comparison on 3DMatch its feature-match recall was 0.359, against 0.952 for the learned FCGF descriptor (Choy, Park and Koltun, 2019, Table 1).

Matching and pre-filtering. Each source descriptor is paired with its nearest neighbour among the target descriptors. Zhou et al. (2016) keep only reciprocal matches, in which each descriptor is the other’s nearest neighbour, and then apply a tuple test: three random matches are kept only if every pairwise distance among the three source points is within a factor τ=0.9\tau = 0.9 to 1/τ1/\tau of the corresponding target distance, a constraint that every rigid motion preserves. 3DMatch’s registration pipeline likewise keeps mutually closest keypoint descriptors before RANSAC (Zeng et al., 2017).

RANSAC. Three correspondences determine a rigid motion through the Kabsch fit; RANSAC draws many such triples, fits each, and keeps the best hypothesis by some score. Two scores are in common use, and the example shows they behave differently. One counts the putative correspondences that the hypothesis maps to within a distance. The other validates the hypothesis on the whole cloud: Open3D’s RegistrationRANSACBasedOnCorrespondence (version 0.20.0) transforms the entire source cloud for each hypothesis that passes cheap checks, such as the edge-length test above, and keeps the one with the largest fraction of source points within the distance of the target, breaking ties by the smaller inlier RMSE; PCL’s SampleConsensusPrerejective, which credits Buch et al. (2013), likewise evaluates each surviving hypothesis on the nearest-neighbour distances of all source points and keeps the one with the lowest mean squared error of its inliers among those whose inlier fraction reaches a set minimum (PCL 1.15.1 source). Both libraries stop after a maximum number of iterations. Open3D can also stop early once the iterations reach log⁡(1−confidence)/log⁡(1−wn)\log(1-\text{confidence})/\log(1-w^{n}) for nn-point samples; its header documents ww as the fitness, while its code uses the fraction of the putative correspondences that the best hypothesis so far maps within the distance. A RANSAC result is coarse, since it comes from three noisy points, and is refined by ICP. Zhou et al. (2016) criticize this two-stage design: most of RANSAC’s effort goes into validating hypotheses that are discarded, and both stages need nearest-neighbour queries inside their loops.

Fast global registration. Zhou, Park and Koltun (2016) replace sampling with a single robust optimization over a fixed set of FPFH correspondences K\mathcal{K}:

E(T)=∑(p,q)∈Kρ ⁣(∥p−Tq∥),ρ(x)=μx2μ+x2,E(T) = \sum_{(\mathbf{p},\mathbf{q}) \in \mathcal{K}} \rho\!\left(\lVert \mathbf{p} - T\mathbf{q} \rVert\right), \qquad \rho(x) = \frac{\mu x^2}{\mu + x^2},

a scaled Geman–McClure penalty that is quadratic for small residuals and flat for large ones, so wrong matches lose their influence. Using the Black–Rangarajan duality they introduce a weight lp,ql_{\mathbf{p},\mathbf{q}} per correspondence and alternate two steps that decrease the same objective: with the weights fixed, a weighted least-squares step on TT, linearized as a six-vector and solved by Gauss–Newton; with TT fixed, the closed-form weights l=(μ/(μ+∥p−Tq∥2))2l = \left(\mu / (\mu + \lVert\mathbf{p} - T\mathbf{q}\rVert^2)\right)^2. To avoid local minima they use graduated non-convexity: μ\mu starts at D2D^2, with DD the diameter of the largest surface, and is halved every four iterations down to δ2\delta^2, with δ\delta the distance threshold for genuine correspondences. Correspondences are never recomputed, so no nearest-neighbour search runs in the loop. On their synthetic range images with Gaussian noise of 0.005 surface diameters, the average RMSE of their method was 0.008 diameters against at least 0.017 for the other global methods they tested, at an average of 0.22 s per pair against 11.1 s for the fastest of them (Tables 1 and 2); on the UWA object-in-scene benchmark its 0.05-recall was 84%, comparable to PCL’s RANSAC pipeline (82%); on the scene-fragment benchmark of Choi et al. its recall was 51.1% against 59.2% for Choi et al.’s own pipeline, with the highest precision (23.2%) and at least an order of magnitude less time (Table 5). The TEASER authors report that FGR tends to fail when more than 80% of the correspondences are outliers (Yang and Carlone, 2019).

TEASER. Yang and Carlone (2019) formulate registration with a truncated least squares cost, ∑imin⁡ ⁣(∥bi−sRai−t∥2/βi2, cˉ2)\sum_i \min\!\left(\lVert \mathbf{b}_i - sR\mathbf{a}_i - \mathbf{t}\rVert^2/\beta_i^2,\ \bar{c}^2\right), in which a correspondence whose residual exceeds a noise bound contributes a constant and has no influence. They decouple the unknowns: differences of correspondences, bj−bi\mathbf{b}_j - \mathbf{b}_i, do not depend on the translation (translation-invariant measurements), and their norms do not depend on rotation either, so the ratio of norms measures the scale. Scale and each component of translation are then estimated exactly by an adaptive voting scheme, outliers are pruned by finding a maximal clique in a graph whose edges join mutually consistent correspondences, and rotation is solved through a semidefinite relaxation that they found tight in practice, with per-instance bounds on its performance. The noise is assumed unknown but bounded, ∥ϵi∥≤βi\lVert\boldsymbol{\epsilon}_i\rVert \le \beta_i. The paper reports tolerance of up to 99% outlier correspondences. The journal version (Yang, Shi and Carlone, 2021) could not be read for this article; in the TEASER++ library (version 2.0) the default rotation solver is graduated non-convexity with the truncated least squares cost, and the default inlier selection is an exact maximum clique.

Correspondence-free global search: Go-ICP

Go-ICP (Yang, Li and Jia, 2013) finds the global minimum of ICP’s own objective, the sum of squared closest-point distances, by branch and bound. Rotations are parameterized by angle-axis vectors in the cube [−π,π]3[-\pi, \pi]^3 that encloses the ball of radius π\pi, translations by a bounded cube, and both are subdivided as octrees. For a sub-cube, the error at its centre is an upper bound; a lower bound follows from an uncertainty radius, the furthest that any rotation or translation in the cube can move a transformed point, by subtracting that radius from each closest-point distance. An outer search over rotations calls an inner search over translations, and whenever the search finds a better pose, local ICP is run from it to lower the upper bound. Closest-point distances come from a precomputed distance transform on a 3003300^3 grid. With 1,000 data points registered to models of 30,000–40,000 points, a run took about 30 s for the Stanford bunny and 15 s for a hand on a 3.4 GHz CPU; a trimmed variant handles outliers. Zhou et al. (2016) ran Go-ICP on 1,000-point subsets, as its authors suggest, because full-resolution clouds took hours. The journal version (Yang, Li, Campbell and Jia, 2016) could not be read for this article.

Probabilistic and non-rigid registration: coherent point drift

Coherent point drift (Myronenko, Song and Carreira-Perpiñán, 2006) treats the template points YY as the centroids of an equally weighted Gaussian mixture with isotropic variance σ2\sigma^2 and fits it to the reference points XX by maximum likelihood, with a prior that the displacement field is smooth. The regularizer penalizes the high-frequency content of the velocity field through a Gaussian low-pass filter, which the authors show is equivalent to the motion coherence theory’s penalty on derivatives of all orders; the optimal displacement then has the kernel form v(z)=∑mwmG(z−y0m)v(\mathbf{z}) = \sum_m \mathbf{w}_m G(\mathbf{z} - \mathbf{y}_{0m}). An EM algorithm alternates soft correspondences (posterior probabilities of each mixture component for each data point, the E-step) with a linear system for the kernel weights (the M-step), and deterministic annealing shrinks the capture range σ\sigma from large to small, a coarse-to-fine strategy. A uniform component in the mixture absorbs outliers. Each iteration costs O(M3)O(M^3) for MM template points, reducible to O(M2)O(M^2) with truncated kernels or conjugate gradients. On 2D fish shapes with missing parts and noise, CPD recovered the deformation where ICP stayed in a local minimum. A journal version (Myronenko and Song, 2010) followed; it could not be read for this article. Soft assignment of this kind avoids ICP’s hard nearest-neighbour decisions, but the authors note that their EM algorithm converges to a local optimum, so, in our reading, it depends on its starting position as ICP does.

Learned registration

Learning entered registration first through descriptors. 3DMatch (Zeng et al., 2017) trains a Siamese 3D ConvNet on 30×30×3030 \times 30 \times 30 truncated-distance-function patches around keypoints, using correspondences mined automatically from existing RGB-D reconstructions of 62 scenes (54 for training, 8 for testing), and registers fragments by RANSAC over mutually closest descriptor matches. FCGF (Choy et al., 2019) computes 32-dimensional descriptors for all points of a scan in one pass of a sparse fully convolutional network; on 3DMatch its registration recall with RANSAC averaged 0.82 over the test scenes, against 0.40 for FPFH and 0.67 for the 3DMatch descriptor (Table 5). PREDATOR (Huang et al., 2021) adds an overlap-attention module that exchanges information between the two clouds’ encodings, predicts for each point how likely it is to lie in the overlap and to be matchable, and samples RANSAC keypoints in proportion to those scores; in Table 1 of the paper this sampling gives a higher registration recall than taking the top-scoring points in nine of ten settings, which the authors attribute to the top-scoring points’ inliers clustering too closely to fix the transform.

A parallel line embeds pose estimation in the network. As PREDATOR’s related work summarizes it, PointNetLK combines a PointNet global descriptor with a Lucas–Kanade-like iteration; DCP computes soft correspondences from learned features and solves for the pose with the Kabsch algorithm; PRNet adds keypoint detection for partial correspondence; RPM-Net updates the similarity matrix with a differentiable Sinkhorn layer and a weighted Kabsch fit; and other methods learn to weight putative correspondences before a weighted Kabsch step. The PREDATOR authors could not train PRNet or RPM-Net to reasonable results on 3DMatch, and RPM-Net with FCGF features reached 54.9% registration recall there against 85.1% for FCGF with RANSAC; they concluded that such direct methods suited geometrically simple objects in controlled settings such as ModelNet40.

Example: local against global registration by initial rotation

The ICP article’s example compares two ICP metrics at one initial pose on fully overlapping samples of a height field. This example asks a different question: how the success of local and global registration depends on the initial rotation when two scans overlap only partly. A small object made of two boxes, a sphere and a capped cylinder, about 1 m across, is ray cast by two range sensors (160×120160 \times 120 pixels, f=200f = 200 px) about 2.2 m away, whose viewing directions differ by 40–60° in azimuth; each scan is expressed in its own sensor frame (camera looking along +z+z, xx right, yy down) and is downsampled by a 2 cm voxel grid that keeps the centroid of each voxel. Depth noise is either absent or Gaussian with a standard deviation of 5 mm. The transform to recover maps the source scan B to the target scan A, xA=RxB+t\mathbf{x}_A = R\mathbf{x}_B + \mathbf{t}.

The global method follows the pipeline above: normals by principal component analysis of the neighbours within 4 cm, oriented toward the sensor; FPFH with a 10 cm radius, written in NumPy after PCL’s tutorial formula (including the query’s own SPFH, with weights the inverse of the distance) and its pair-feature code, with each 11-bin block scaled to sum to 100 as in PCL’s code; reciprocal matches; and 5,000 three-point RANSAC hypotheses filtered by the edge-length test with τ=0.9\tau = 0.9. Each hypothesis is scored both ways described above, since the choice changes the result: by the number of putative matches within 3 cm (with a final Kabsch refit on them), and by the dense fraction of all source points within 3 cm of the target as in Open3D 0.20.0 (no refit). The dense result is then refined by point-to-plane ICP with a 0.03 m rejection distance. ICP alone, point-to-point or point-to-plane as in the ICP article, starts from the true pose turned by θ\theta about a random axis through the source centroid, so the initial error is a pure rotation about the object. Pairs farther apart than a rejection distance are discarded; the ICP article’s example used 0.1 m on fully overlapping clouds, and because the result depends on this choice, both 0.1 m and 0.03 m are run. Errors are the relative rotation error RRE=arccos⁡ ⁣((tr⁡(R^⊤R)−1)/2)\text{RRE} = \arccos\!\left((\operatorname{tr}(\hat{R}^\top R) - 1)/2\right) and the relative translation error RTE=∥t^−t∥\text{RTE} = \lVert \hat{\mathbf{t}} - \mathbf{t} \rVert, the translation of the sensor origin, as FCGF and PREDATOR define them. A registration counts as a success when RRE < 5° and RTE < 10 cm; the 2 m of KITTI’s criterion (below) is for driving scenes, not a 1 m object. Each of the 20 trials per noise level (seeds 0–19) draws new sensor positions and noise. The script uses NumPy 2.3.5 and SciPy 1.16.3 and takes one to one and a half minutes.

This set-up was revised twice before its first run, each time after the editor reviewed a pilot on noise-free scans with a 5° initial rotation that every compared method had to pass. In the first pilot, ICP with the planned 0.1 m rejection distance and RANSAC with the planned correspondence-count score did not pass, and diagnostic runs looked for a broken component without finding one. The 0.03 m rejection distance and the dense score were then added, with both planned variants kept and printed; the printout of each chosen score against the true pose’s score was added; and the planned ICP refinement of the count-scored result was replaced by point-to-plane ICP at 0.03 m on the dense result. In the second pilot, dense scoring alone did not pass, and the editor set dense scoring followed by that ICP refinement as the global method’s base check, with dense scoring alone still printed. The script has not changed since its first run.

import numpy as np
from scipy.spatial import cKDTree
from scipy.spatial.transform import Rotation

# Scene: two boxes, a sphere and a capped cylinder (metres, z up), about 1 m across.
BOXES = [((-0.45, -0.30, 0.00), (0.45, 0.30, 0.30)), ((-0.37, -0.02, 0.30), (-0.13, 0.22, 0.60))]
SPHERE = ((0.20, -0.08, 0.45), 0.15)
CYL = ((0.30, 0.20), 0.08, 0.30, 0.65)                  # axis (x, y), radius, z0, z1
F, W, H = 200.0, 160, 120                               # range sensor: focal length (px), image size
LOOK_AT = np.array([0.0, 0.0, 0.3])
VOX, NRAD, FRAD, THR = 0.02, 0.04, 0.10, 0.03           # voxel, normal radius, FPFH radius, inlier distance

def camera(az, el, dist):
    """Camera-to-world rotation (columns: x right, y down, z forward) and camera centre."""
    c = LOOK_AT + dist * np.array([np.cos(el) * np.cos(az), np.cos(el) * np.sin(az), np.sin(el)])
    z = (LOOK_AT - c) / np.linalg.norm(LOOK_AT - c)
    x = np.cross(z, [0, 0, 1.0]); x /= np.linalg.norm(x)
    return np.column_stack([x, np.cross(z, x), z]), c

def first_hit(o, d):
    """Ray parameter of the first intersection of rays o + s d with the scene (inf if none)."""
    s = np.full(len(d), np.inf)
    for lo, hi in BOXES:                                # slab test
        with np.errstate(divide="ignore", invalid="ignore"):
            t1, t2 = (np.array(lo) - o) / d, (np.array(hi) - o) / d
        tmin, tmax = np.nanmax(np.minimum(t1, t2), 1), np.nanmin(np.maximum(t1, t2), 1)
        s = np.where((tmax >= tmin) & (tmin > 0), np.minimum(s, tmin), s)
    (ctr, r) = SPHERE
    oc = o - np.array(ctr); a = np.sum(d * d, 1); b = 2 * d @ oc
    disc = b * b - 4 * a * (oc @ oc - r * r)
    t = (-b - np.sqrt(np.maximum(disc, 0))) / (2 * a)
    s = np.where((disc >= 0) & (t > 0), np.minimum(s, t), s)
    (cx, cy), r, z0, z1 = CYL
    ox, oy = o[0] - cx, o[1] - cy
    a = d[:, 0] ** 2 + d[:, 1] ** 2; b = 2 * (ox * d[:, 0] + oy * d[:, 1])
    disc = b * b - 4 * a * (ox ** 2 + oy ** 2 - r * r)
    t = (-b - np.sqrt(np.maximum(disc, 0))) / (2 * a); z = o[2] + t * d[:, 2]
    s = np.where((disc >= 0) & (t > 0) & (z >= z0) & (z <= z1), np.minimum(s, t), s)
    with np.errstate(divide="ignore", invalid="ignore"):
        t = (z1 - o[2]) / d[:, 2]                       # top cap
    s = np.where((t > 0) & ((ox + t * d[:, 0]) ** 2 + (oy + t * d[:, 1]) ** 2 <= r * r), np.minimum(s, t), s)
    return s

def render(R_wc, c, sigma, rng):
    """Range scan in the sensor's own frame, Gaussian noise on depth."""
    u, v = np.meshgrid(np.arange(W), np.arange(H))
    d_c = np.column_stack([(u.ravel() - (W - 1) / 2) / F, (v.ravel() - (H - 1) / 2) / F, np.ones(W * H)])
    s = first_hit(c, d_c @ R_wc.T)                      # unit z component, so s is the depth
    ok = np.isfinite(s)
    return d_c[ok] * (s[ok] + sigma * rng.standard_normal(ok.sum()))[:, None]

def voxel_down(p, size):
    _, inv = np.unique(np.floor(p / size).astype(np.int64), axis=0, return_inverse=True)
    inv = inv.ravel(); out = np.zeros((inv.max() + 1, 3)); np.add.at(out, inv, p)
    return out / np.bincount(inv)[:, None]

def normals(p, radius):
    n = np.zeros_like(p)
    for i, nb in enumerate(cKDTree(p).query_ball_point(p, radius)):
        q = p[nb] - p[nb].mean(0); n[i] = np.linalg.eigh(q.T @ q)[1][:, 0]
    n[np.sum(n * p, 1) > 0] *= -1                       # orient toward the sensor (origin)
    return n

def fpfh(p, n, radius, bins=11):
    """FPFH after PCL 1.15.1: SPFH of (alpha, phi, theta) pairs, then distance-weighted sum."""
    nbrs = cKDTree(p).query_ball_point(p, radius)
    i = np.concatenate([np.full(len(nb), k) for k, nb in enumerate(nbrs)]); j = np.concatenate(nbrs)
    i, j = i[i != j], j[i != j]
    d = p[j] - p[i]; dist = np.linalg.norm(d, axis=1)
    n1, n2 = n[i].copy(), n[j].copy()
    a1, a2 = np.sum(n1 * d, 1) / dist, np.sum(n2 * d, 1) / dist
    sw = np.arccos(np.clip(np.abs(a1), 0, 1)) > np.arccos(np.clip(np.abs(a2), 0, 1))
    n1[sw], n2[sw] = n[j][sw], n[i][sw]; d[sw] *= -1    # source = point with smaller angle
    f3 = np.where(sw, -a2, a1)
    v = np.cross(d, n1); vn = np.linalg.norm(v, axis=1); good = vn > 0
    v /= np.where(good, vn, 1)[:, None]; w = np.cross(n1, v)
    f2 = np.sum(v * n2, 1)
    f1 = np.arctan2(np.sum(w * n2, 1), np.sum(n1 * n2, 1))
    k = np.bincount(i, minlength=len(p)).astype(float)
    spfh = np.zeros((len(p), 3 * bins))
    for c, (f, lo, hi) in enumerate([(f1, -np.pi, np.pi), (f2, -1, 1), (f3, -1, 1)]):
        b = np.clip(np.floor(bins * (f - lo) / (hi - lo)).astype(int), 0, bins - 1)
        np.add.at(spfh, (i[good], c * bins + b[good]), 100.0 / k[i[good]])
    out = spfh.copy()
    np.add.at(out, i, spfh[j] / dist[:, None] / k[i][:, None])   # SPFH(p) + 1/k sum SPFH(p_i)/w_i
    for c in range(3):
        blk = out[:, c * bins:(c + 1) * bins]; blk /= np.maximum(blk.sum(1, keepdims=True), 1e-12) / 100
    return out

def kabsch(p, q):
    p0, q0 = p.mean(0), q.mean(0)
    U, _, Vt = np.linalg.svd((p - p0).T @ (q - q0))
    R = Vt.T @ np.diag([1, 1, np.sign(np.linalg.det(Vt.T @ U.T))]) @ U.T
    return R, q0 - R @ p0

def kabsch_batch(P, Q):
    p0, q0 = P.mean(1, keepdims=True), Q.mean(1, keepdims=True)
    U, _, Vt = np.linalg.svd(np.einsum("bki,bkj->bij", P - p0, Q - q0))
    V, Ut = np.swapaxes(Vt, 1, 2), np.swapaxes(U, 1, 2)
    D = np.zeros((len(P), 3, 3)); D[:, 0, 0] = D[:, 1, 1] = 1; D[:, 2, 2] = np.sign(np.linalg.det(V @ Ut))
    R = V @ D @ Ut
    return R, q0[:, 0] - np.einsum("bij,bj->bi", R, p0[:, 0])

def mutual_matches(fs, ft):
    _, s2t = cKDTree(ft).query(fs); _, t2s = cKDTree(fs).query(ft)
    src = np.arange(len(fs)); keep = t2s[s2t] == src
    return src[keep], s2t[keep]

def ransac(ps, pt, src, tree, n_hyp, rng, tau=0.9, batch=500):
    """3-point hypotheses with an edge-length test, scored two ways.
    count: putative matches within THR (best hypothesis refitted on them).
    dense (Open3D 0.20 rule): fraction of all source points within THR of the target,
    ties broken by lower inlier RMSE; the sample's own transform is returned."""
    idx_all = np.stack([rng.choice(len(ps), 3, replace=False) for _ in range(n_hyp)])
    best_c, best_d = (-1, None), ((-1.0, -np.inf), None)
    for start in range(0, n_hyp, batch):
        idx = idx_all[start:start + batch]; P, Q = ps[idx], pt[idx]
        ok = np.ones(len(idx), bool)
        for a, b in [(0, 1), (0, 2), (1, 2)]:
            r = np.linalg.norm(P[:, a] - P[:, b], axis=1) / np.maximum(np.linalg.norm(Q[:, a] - Q[:, b], axis=1), 1e-12)
            ok &= (r > tau) & (r < 1 / tau)
        if not ok.any():
            continue
        R, t = kabsch_batch(P[ok], Q[ok])
        cnt = (np.linalg.norm(np.einsum("bij,mj->bmi", R, ps) + t[:, None] - pt[None], axis=2) < THR).sum(1)
        if cnt.max() > best_c[0]:
            best_c = (cnt.max(), (R[cnt.argmax()], t[cnt.argmax()]))
        dd = tree.query((np.einsum("bij,mj->bmi", R, src) + t[:, None]).reshape(-1, 3))[0].reshape(len(R), -1)
        inl = dd < THR; fit = inl.mean(1)
        rmse = np.sqrt(np.where(inl, dd ** 2, 0).sum(1) / np.maximum(inl.sum(1), 1))
        for k in np.flatnonzero(fit == fit.max()):
            if (fit[k], -rmse[k]) > best_d[0]:
                best_d = ((fit[k], -rmse[k]), (R[k], t[k]))
    R, t = best_c[1]
    inl = np.linalg.norm(ps @ R.T + t - pt, axis=1) < THR
    return kabsch(ps[inl], pt[inl]), best_d[1]

def icp(src, tgt, tree, tn, R, t, metric, reject, max_iter=100, tol=1e-7):
    prev = np.inf
    for _ in range(max_iter):
        moved = src @ R.T + t
        d, j = tree.query(moved); keep = d < reject
        if keep.sum() < 6:
            break
        if metric == "p2p":
            dR, dt = kabsch(moved[keep], tgt[j[keep]]); err = np.mean(d[keep] ** 2)
        else:                                           # linearized point-to-plane step
            p, q, nn = moved[keep], tgt[j[keep]], tn[j[keep]]
            b = np.sum((q - p) * nn, 1)
            x = np.linalg.lstsq(np.hstack([np.cross(p, nn), nn]), b, rcond=None)[0]
            dR, dt = Rotation.from_rotvec(x[:3]).as_matrix(), x[3:]; err = np.mean(b ** 2)
        R, t = dR @ R, dR @ t + dt
        if prev - err < tol:
            break
        prev = err
    return R, t

def errors(R, t, Rg, tg):
    """RRE (degrees) and RTE (metres) against the true transform."""
    return np.degrees(Rotation.from_matrix(R @ Rg.T).magnitude()), np.linalg.norm(t - tg)

def ok(e):
    return e[0] < 5 and e[1] < 0.10                     # success: RRE < 5 deg and RTE < 10 cm

THETAS = [5, 15, 30, 45, 60, 90, 120, 180]
ICPS = [(m, r) for m in ["p2p", "p2plane"] for r in [0.10, 0.03]]

def trial(seed, sigma):
    rng = np.random.default_rng(seed)
    az = rng.uniform(0, 2 * np.pi)
    RA, cA = camera(az, np.radians(30), 2.2)            # sensor A: target scan
    RB, cB = camera(az + np.radians(50 + rng.uniform(-10, 10)), np.radians(30 + rng.uniform(-5, 5)),
                    2.2 + rng.uniform(-0.1, 0.1))       # sensor B: source scan
    tgt = voxel_down(render(RA, cA, sigma, rng), VOX)
    src = voxel_down(render(RB, cB, sigma, rng), VOX)
    Rg, tg = RA.T @ RB, RA.T @ (cB - cA)                # truth: x_A = Rg x_B + tg
    tn, sn = normals(tgt, NRAD), normals(src, NRAD)
    tree = cKDTree(tgt)
    fit = lambda R, t: np.mean(tree.query(src @ R.T + t)[0] < THR)
    si, ti = mutual_matches(fpfh(src, sn, FRAD), fpfh(tgt, tn, FRAD))
    ps, pt = src[si], tgt[ti]
    count = lambda R, t: np.sum(np.linalg.norm(ps @ R.T + t - pt, axis=1) < THR)
    (Rc, tc), (Rd, td) = ransac(ps, pt, src, tree, 5000, rng)
    out = {"overlap": fit(Rg, tg), "matches": len(si), "inlier ratio": count(Rg, tg) / len(si),
           "count": (errors(Rc, tc, Rg, tg), count(Rc, tc) >= count(Rg, tg)),
           "dense": (errors(Rd, td, Rg, tg), fit(Rd, td) >= fit(Rg, tg)),
           "dense+ICP": errors(*icp(src, tgt, tree, tn, Rd, td, "p2plane", 0.03), Rg, tg)}
    c = src.mean(0)
    for th in THETAS:                                   # initial estimate: truth rotated by th about
        axis = rng.standard_normal(3); axis /= np.linalg.norm(axis)   # a random axis through the source centroid
        R0 = Rotation.from_rotvec(np.radians(th) * axis).as_matrix() @ Rg
        t0 = Rg @ c + tg - R0 @ c
        for m, r in ICPS:
            out[(m, r, th)] = errors(*icp(src, tgt, tree, tn, R0, t0, m, r), Rg, tg)
    return out

N = 20
for sigma in [0.0, 0.005]:
    res = [trial(seed, sigma) for seed in range(N)]
    print(f"depth noise {1000 * sigma:.0f} mm, {N} scan pairs (seeds 0-{N - 1})")
    ov = [r["overlap"] for r in res]
    print(f"  overlap {min(ov):.2f}-{max(ov):.2f}; mutual FPFH matches median {np.median([r['matches'] for r in res]):.0f},"
          f" inlier ratio median {np.median([r['inlier ratio'] for r in res]):.2f}")
    for key, label in [("count", "RANSAC, count score "), ("dense", "RANSAC, dense score ")]:
        e = np.array([r[key][0] for r in res])
        print(f"  {label}: success {sum(ok(x) for x in e):2d}/{N}  median RRE {np.median(e[:, 0]):5.2f} deg"
              f"  RTE {100 * np.median(e[:, 1]):5.1f} cm  chosen score >= true pose's in {sum(r[key][1] for r in res):2d}/{N}")
    e = np.array([r["dense+ICP"] for r in res])
    print(f"  dense RANSAC + ICP  : success {sum(ok(x) for x in e):2d}/{N}  median RRE {np.median(e[:, 0]):5.2f} deg"
          f"  RTE {100 * np.median(e[:, 1]):5.1f} cm")
    print("  ICP alone, successes by initial rotation (reject = pair rejection distance)")
    print("   theta  " + "  ".join(f"{m:>7s} {r:.2f} m" for m, r in ICPS))
    for th in THETAS:
        print(f"   {th:4d}   " + "  ".join(f"{sum(ok(r[(m, rr, th)]) for r in res):11d}  " for m, rr in ICPS))

Output:

depth noise 0 mm, 20 scan pairs (seeds 0-19)
  overlap 0.58-0.92; mutual FPFH matches median 480, inlier ratio median 0.11
  RANSAC, count score : success  9/20  median RRE  5.07 deg  RTE  14.0 cm  chosen score >= true pose's in 18/20
  RANSAC, dense score : success 14/20  median RRE  3.23 deg  RTE   7.8 cm  chosen score >= true pose's in  6/20
  dense RANSAC + ICP  : success 20/20  median RRE  0.12 deg  RTE   0.4 cm
  ICP alone, successes by initial rotation (reject = pair rejection distance)
   theta      p2p 0.10 m      p2p 0.03 m  p2plane 0.10 m  p2plane 0.03 m
      5            15             20             15             20  
     15            15              8             14             18  
     30            13              0             12              5  
     45             4              0              4              0  
     60             1              0              1              0  
     90             0              0              0              0  
    120             0              0              0              0  
    180             0              0              0              0  
depth noise 5 mm, 20 scan pairs (seeds 0-19)
  overlap 0.59-0.91; mutual FPFH matches median 387, inlier ratio median 0.08
  RANSAC, count score : success  8/20  median RRE  6.87 deg  RTE  25.9 cm  chosen score >= true pose's in 17/20
  RANSAC, dense score : success 12/20  median RRE  3.18 deg  RTE   9.4 cm  chosen score >= true pose's in  5/20
  dense RANSAC + ICP  : success 19/20  median RRE  0.11 deg  RTE   0.4 cm
  ICP alone, successes by initial rotation (reject = pair rejection distance)
   theta      p2p 0.10 m      p2p 0.03 m  p2plane 0.10 m  p2plane 0.03 m
      5            14             20             15             20  
     15            14              7             15             18  
     30            13              0             12              4  
     45             5              0              4              0  
     60             1              0              2              0  
     90             0              0              0              0  
    120             0              0              0              0  
    180             0              0              0              0  

Local registration has a narrow basin, and its width depends on pair rejection. With a 0.03 m rejection distance both ICP metrics recover every pose from 5°; from 15° point-to-plane still succeeds in 18 of 20 trials but point-to-point in only 7 or 8, and from 30° point-to-plane succeeds a few times and point-to-point never. With 0.1 m the basin reaches further (12–13 successes from 30°, a few from 45°), but even from 5° five or six trials fail. The printed overlap, the fraction of source points within 3 cm of the target at the true pose, is between 0.58 and 0.92, so every pair has points seen by one sensor only. In our reading of the construction, many of those points still find some target point within 10 cm on an object 1 m across, enter the fit with full weight and pull the pose away from the truth; at 3 cm they are rejected, but far from the solution too few correct pairs survive to pull the pose in. A fixed rejection distance thus trades accuracy for reach, and noise changes little of this. No ICP variant recovers any pose from 90° or more.

The global method does not depend on the initial rotation, since FPFH is invariant to rigid motion and RANSAC never uses the initial estimate; it is computed once per trial. Its success is limited instead by the matches: the median inlier ratio of the reciprocal FPFH matches is 0.11 without noise and 0.08 with it. The two scores then differ. Counting matched inliers succeeds in fewer than half the trials, and in 17–18 of 20 the chosen pose (after its refit) has at least as many inliers as the true pose, so the failure lies in the score, not in the sampling: in our reading, the flat faces and right-angled edges of the boxes give many look-alike FPFHs, so poses that map one face onto another explain as many putative matches as the correct pose. Scoring every source point against the target succeeds more often (14 and 12 of 20), but with a median rotation error near 3° it is a coarse stage: in only 5–6 trials does the chosen hypothesis score at least as well as the true pose, so in the rest no sampled three-point fit came as close to the truth as the score could have rewarded. Refined by point-to-plane ICP, the dense result succeeds in 20 of 20 trials without noise and 19 of 20 with it, with a median rotation error near 0.1°. This is the division of labour that Zhou et al. (2016) describe: the global stage brings the clouds into ICP’s basin, and ICP supplies the accuracy.

Datasets and Benchmarks

  • 3DMatch (Zeng et al., 2017). Fragments built from RGB-D reconstructions of 62 indoor scenes (from Analysis-by-Synthesis, 7-Scenes, SUN3D, RGB-D Scenes v.2 and Halber et al.), 54 used for training and 8 for testing in the original paper. Its geometric registration benchmark adopts the evaluation of Choi et al.: pairs of non-consecutive fragments are registered, and an estimate is a true positive if the pair overlaps by more than 30% and the RMSE of ground-truth correspondences after applying the estimate is below τ=0.2\tau = 0.2 (the paper prints no unit; FCGF states it as 0.2 m). The paper also evaluated on Choi et al.’s synthetic benchmark of 207 fragments from four ICL-NUIM scenes, where its descriptor with RANSAC reached 65.1% recall and 25.2% precision against 46.1% and 19.1% for FPFH (“Rusu et al.”) with the same RANSAC pipeline (Table 2); on the real-scene benchmark the figures were 66.8% and 40.1% against 44.2% and 30.7% (Table 3). Later work, such as PREDATOR, uses 46 scenes for training, 8 for validation and 8 for testing.
  • 3DLoMatch (Huang et al., 2021). The 3DMatch test scans restricted to pairs with 10–30% overlap. The authors note that a bug in the official 3DMatch overlap computation puts a few (under 7%) pairs in both sets. With 5,000 sampled points, PREDATOR’s registration recall was 88.3% on 3DMatch and 54.2% on 3DLoMatch, against 85.1% and 40.1% for FCGF (Table 2).
  • KITTI odometry (KITTI). LiDAR scans from the odometry training sequences 00–10, whose GPS/IMU poses are refined with ICP. FCGF’s protocol pairs each scan with the first scan at least 10 m away, removes pairs where ICP fails or fewer than 1,000 voxels overlap, and yields 1,358 training, 180 validation and 555 test pairs; its text assigns sequences 0–5 to training, “7 to 8” to validation and “8 to 10” to testing, which overlap at sequence 8, while PREDATOR, following FCGF, gives 0–5, 6–7 and 8–10. PREDATOR’s text describes its evaluation pairs as at most 10 m apart, where FCGF’s are at least 10 m apart. Success is RTE < 2 m and RRE < 5°, following 3DFeat-Net. With 30 cm voxels, FCGF with RANSAC reached a 98.92% success rate, with an RTE of 6.5 cm and an RRE of 0.23° (Table 6); PREDATOR reported 99.8% registration recall, 6.8 cm RTE and 0.27° RRE (Table 5), under D3Feat’s evaluation protocol, in a table that lists FCGF at 96.6%, 9.5 cm and 0.30°, so the two papers’ FCGF figures are not comparable.
  • Synthetic and object benchmarks. Choi et al.’s Augmented ICL-NUIM fragments, used for both pairwise and multi-way registration (Zhou et al., 2016); the UWA object-in-scene benchmark, with 188 registration tests and overlaps as low as 21% (Zhou et al., 2016); synthetic range images of standard models with Gaussian noise (Zhou et al., 2016); and partial scans of ModelNet40 CAD models for learned direct registration, with PREDATOR’s lower-overlap variant ModelLoNet (Huang et al., 2021).

Evaluation Metrics

  • Relative rotation and translation error. RRE=arccos⁡ ⁣((tr⁡(R^⊤R∗)−1)/2)\text{RRE} = \arccos\!\left((\operatorname{tr}(\hat{R}^\top R^*) - 1)/2\right), the geodesic angle between estimated and true rotations, and RTE=∥t^−t∗∥\text{RTE} = \lVert\hat{\mathbf{t}} - \mathbf{t}^*\rVert (Choy et al., 2019; Huang et al., 2021). RTE depends on where the frame origin is: for sensor-centred scans it is the error in the sensor’s position, which a small rotation error about a distant object also inflates.
  • Success rate. The fraction of pairs whose RRE and RTE are below thresholds, such as 5° and 2 m on KITTI.
  • Registration recall (RR). The fraction of pairs, among those with enough overlap, for which the RMSE of the ground-truth correspondences after the estimated transform,
RMSE=1∣H∗∣∑(p,q)∈H∗∥T^(p)−q∥2,\text{RMSE} = \sqrt{\frac{1}{\lvert\mathcal{H}^*\rvert}\sum_{(\mathbf{p},\mathbf{q})\in\mathcal{H}^*} \lVert \hat{T}(\mathbf{p}) - \mathbf{q}\rVert^2},

is below 0.2 m. PREDATOR’s supplementary material calls it the most reliable metric because it measures the end task, and notes that immediately adjacent fragments are excluded, as in the original 3DMatch evaluation script; FCGF credits it to Choi, Zhou and Koltun’s 2015 work on indoor reconstruction and argues that recall matters more than precision, since precision can be improved by pruning.

  • Feature-match recall (FMR) and inlier ratio (IR). IR is the fraction of putative matches within τ1=10\tau_1 = 10 cm of each other under the true transform; FMR is the fraction of pairs whose IR exceeds τ2=5%\tau_2 = 5\% (Choy et al., 2019). They measure descriptors, not registration: PREDATOR’s authors note that FMR does not check whether the transformation can be recovered from the matches, and that a high inlier ratio need not give a high registration recall when the inliers are clustered too closely.
  • RMSE recall and precision. Zhou et al. (2016) plot “α-recall”, the fraction of tests with RMSE below α in units of the surface diameter; Choi et al.’s scene benchmark reports precision as well as recall.
  • Chamfer distance between registered scans, used on synthetic object data (Huang et al., 2021).

Applications

  • 3D scanning and reconstruction. Registering range scans or RGB-D fragments, then fusing them with methods such as TSDF fusion; 3DMatch’s descriptors were trained on, and used to improve, such reconstructions (Zeng et al., 2017).
  • LiDAR odometry and mapping. Scan-to-scan and scan-to-map registration in LiDAR SLAM; the SLAM article describes LOAM’s split into fast odometry and slower map registration.
  • RGB-D tracking. KinectFusion tracks a depth camera by point-to-plane ICP against a ray-cast model.
  • Object pose estimation and robotics. Registering a model to a scene for grasping and robotic manipulation; Go-ICP localized a cap in a cluttered RGB-D scene using geometry alone (Yang et al., 2013).
  • Evaluation and ground truth. KITTI’s disparity and flow ground truth accumulates laser scans registered with ICP, and reconstruction benchmarks align reference and reconstruction with ICP, as the KITTI and 3D reconstruction articles describe.
  • Medical and deformable shapes. Non-rigid point-set registration of surfaces such as faces, on which CPD was demonstrated; CPD’s authors name image-guided surgery among the uses of point-set registration, and deformable registration of images is part of image registration.

Open Problems

  • Low overlap. PREDATOR’s registration recall on 3DLoMatch with 5,000 points, 54.2% (58.0% for its wider variant), is far below its 88% on 3DMatch, and its authors show a case whose overlap is insufficient for registration even with their attention mechanism (Huang et al., 2021).
  • Generalization of learned methods. Learned descriptors are trained on particular sensors and scenes; FCGF used rotation augmentation to make its features invariant to the relative pose, and the direct registration networks that PREDATOR’s authors tried did not train to reasonable results on real indoor scans.
  • Certifiability at scale. Branch and bound guarantees the global optimum but, as Yang and Carlone (2019) note, runs in exponential time in the worst case; TEASER is presented as a polynomial-time algorithm with per-instance bounds under its bounded-noise assumption, although its maximal-clique pruning is exponential in general and, the authors note, fast in practice because the graphs are sparse under high outlier rates; they set this against RANSAC, whose performance guarantees they consider weaker.
  • Symmetry and repetition. As the example shows, scores built on local evidence can prefer a wrong alignment when structures repeat, and no method can resolve a truly symmetric object without extra information such as colour.
  • Non-rigid and dynamic scenes. Coherent point drift scales cubically with the number of points without approximations, and real scenes mix rigid background, moving objects and deforming bodies.
  • Evaluation protocols. Benchmark details such as overlap thresholds, pair selection and train/validation/test splits are described inconsistently across papers, as the KITTI protocol above shows, which makes reported numbers hard to compare.

Related

  • Point Cloud

    A set of 3D points sampled from the surfaces of a scene, possibly with color, intensity, or normals per point, and the most direct output of 3D sensors and reconstruction.

  • RANSAC

    Random sample consensus, a randomized method that fits a model to data containing many outliers by repeatedly fitting minimal random samples and keeping the model that the most data points agree with.

  • Feature Descriptors

    Vectors that summarize the image appearance around a keypoint so that the same scene point seen in different images gets similar vectors, from raw patches and gradient histograms to binary strings and learned embeddings.

  • LiDAR

    An active sensor that measures distance by timing laser light reflected from surfaces, and sweeps its beams across the scene to produce a point cloud with range, intensity and a timestamp for every return.

  • 3D Reconstruction

    Recovering the 3D shape of objects and scenes from images, depth maps, or range scans, as point clouds, meshes, volumes, or implicit and radiance fields, and measuring how accurate and complete the result is.

  • Image Registration

    Aligning two images of the same scene or object, taken at different times, from different viewpoints or with different sensors, by finding the geometric transformation that maps points in one onto the corresponding points in the other.

References

  1. Besl, P. J. & McKay, N. D. (1992). A Method for Registration of 3-D Shapes. IEEE Transactions on Pattern Analysis and Machine Intelligence, 14(2), 239–256 (not read).
  2. Chen, Y. & Medioni, G. (1992). Object Modelling by Registration of Multiple Range Images. Image and Vision Computing, 10(3), 145–155 (not read).
  3. Rusu, R. B., Blodow, N. & Beetz, M. (2009). Fast Point Feature Histograms (FPFH) for 3D Registration. IEEE International Conference on Robotics and Automation (ICRA) (not read).
  4. Point Cloud Library (PCL) 1.15.1: FPFH and PFH estimation tutorials (doc/tutorials/content/fpfh_estimation.rst, pfh_estimation.rst) and source code (features/src/pfh.cpp, features/include/pcl/features/impl/fpfh.hpp, registration/include/pcl/registration/sample_consensus_prerejective.h and impl/sample_consensus_prerejective.hpp) (2025).
  5. Zhou, Q.-Y., Park, J. & Koltun, V. (2016). Fast Global Registration. European Conference on Computer Vision (ECCV), Lecture Notes in Computer Science (read in the authors' copy).
  6. Open3D 0.20.0 source code: cpp/open3d/pipelines/registration (Registration.cpp, Registration.h, CorrespondenceChecker.cpp, Feature.cpp, FastGlobalRegistration.cpp).
  7. Yang, H. & Carlone, L. (2019). A Polynomial-time Solution for Robust Registration with Extreme Outlier Rates. Robotics: Science and Systems (RSS) XV.
  8. Yang, H., Shi, J. & Carlone, L. (2021). TEASER: Fast and Certifiable Point Cloud Registration. IEEE Transactions on Robotics, 37(2), 314–333 (not read).
  9. MIT-SPARK. TEASER++ v2.0 source code: README.md and teaser/include/teaser/registration.h.
  10. Yang, J., Li, H. & Jia, Y. (2013). Go-ICP: Solving 3D Registration Efficiently and Globally Optimally. IEEE International Conference on Computer Vision (ICCV).
  11. Yang, J., Li, H., Campbell, D. & Jia, Y. (2016). Go-ICP: A Globally Optimal Solution to 3D ICP Point-Set Registration. IEEE Transactions on Pattern Analysis and Machine Intelligence, 38(11), 2241–2254 (not read).
  12. Myronenko, A., Song, X. & Carreira-Perpiñán, M. Á. (2006). Non-rigid Point Set Registration: Coherent Point Drift. Advances in Neural Information Processing Systems (NIPS), 19.
  13. Myronenko, A. & Song, X. (2010). Point Set Registration: Coherent Point Drift. IEEE Transactions on Pattern Analysis and Machine Intelligence, 32(12), 2262–2275 (not read).
  14. Zeng, A., Song, S., Nießner, M., Fisher, M., Xiao, J. & Funkhouser, T. (2017). 3DMatch: Learning Local Geometric Descriptors from RGB-D Reconstructions. IEEE Conference on Computer Vision and Pattern Recognition (CVPR).
  15. Choy, C., Park, J. & Koltun, V. (2019). Fully Convolutional Geometric Features. IEEE/CVF International Conference on Computer Vision (ICCV).
  16. Huang, S., Gojcic, Z., Usvyatsov, M., Wieser, A. & Schindler, K. (2021). PREDATOR: Registration of 3D Point Clouds with Low Overlap. IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR) (with its supplementary material).