ImFusion SDK 4.3
DeformableIcp Class Reference

#include <ImFusion/Reg/DeformableIcp.h>

Implements a regularized-least-squares-based iterative closest point (ICP) procedure for deformable registration of two sets of points. More...

+ Inheritance diagram for DeformableIcp:

Detailed Description

Implements a regularized-least-squares-based iterative closest point (ICP) procedure for deformable registration of two sets of points.

The Frobenius norm of a user-defined regularization matrix (e.g. a Laplacian) applied to point displacements is used as regularization term. The regularization is gradually reduced by multiplying the regularization parameter \(\lambda\) by a factor smaller than 1 each iteration. Point correspondences are re-established each iteration using bi- or unidirectional nearest neighbor search. In each iteration, a sparse linear system is solved to compute the closed-form solution of:

\[\mathbf D_i = \arg\min_{\mathbf D} \left\{ \frac{1}{n} \|\mathbf S_i^{y\leftarrow x} \mathbf D - (\mathbf Y-\mathbf S_i^{y\leftarrow x}\mathbf X_0)\|_F^2 + \frac{1}{m} \|\mathbf D - (\mathbf S_i^{x\leftarrow y}\mathbf Y-\mathbf X_0)\|_F^2 + \lambda_i \|\mathbf L \mathbf D\|_F^2 \right\} \]

where

  • \(i\) is the iteration number ( \(i=0\) refers to the initial configuration),
  • \(\mathbf X_i \in \mathbb R^{m\times 3}\) are the \(m\) moving points,
  • \(\mathbf Y \in \mathbb R^{n\times 3}\) are the \(n\) fixed points,
  • \(\mathbf D_i \in \mathbb R^{m\times 3}\) are displacements of the moving points,
  • \(\mathbf S_i^{y\leftarrow x} \in \mathbb R^{n\times m}\) extracts for each point in \(\mathbf Y\) the row of the nearest point in \(\mathbf X_{i-1}\),
  • \(\mathbf S_i^{x\leftarrow y} \in \mathbb R^{m\times n}\) extracts for each point in \(\mathbf X_{i-1}\) the row of the nearest point in \(\mathbf Y\),
  • \(\mathbf L \in \mathbb R^{p\times m}\) is the regularization matrix (typically \(p=m\)), and
  • \(\lambda_i = r^{i-1} \lambda_0\) is the regularization parameter ( \(r\) is the \(\lambda\) reduction factor).

The moving points are then updated with \(\mathbf X_i = \mathbf X_0 + \mathbf D_i\).

Iterations are stopped if the relative RMSE or \(\lambda\) falls below a specified threshold or the maximum number of iterations has been reached. The same stopping criteria apply for an optional iterative rigid or rigid-then-rigid-and-scaling pre-alignment. Additionally, the pre-alignment is also stopped if the difference of relative RMSEs between successive iterations falls below a threshold.

The method is similar to [1] in that it also performs "optimal" steps in each iteration and uses a decreasing regularization parameter, but it differs in the type of regularization (Laplacian or similar applied to displacements vs. difference of locally affine transformations between neighboring points). To speed up the registration, this method also uses a more greedy iteration scheme by omitting the inner loop described in [1]. Furthermore, this method uses bi-directional point correspondences by default, while [1] establishes correspondences with respect to the moving points only.

[1] https://doi.org/10.1109/CVPR.2007.383165

Classes

struct  IterationInfo
 Used to report information about an iteration. More...
 
struct  RansacProperties
 

Public Types

enum class  PreAlignment { None , Rigid , RigidThenRigidAndScaling }
 Pre-alignment modes. More...
 
enum class  NearestNeighborDirection { WrtFixed , WrtMoving , WrtBoth }
 Direction in which point correspondences are established using nearest neighbor search. More...
 
enum class  StoppingReason {
  None , Lambda , RelRmse , RelRmseDiff ,
  MaxIter , Callback
}
 Used to report why iterations have been stopped. Stopping criteria are evaluated in this order. More...
 
enum class  Phase {
  None , Centroids , Rigid , RigidAndScaling ,
  Deformable
}
 Used to report in which phase the algorithm is in. More...
 

Public Member Functions

 DeformableIcp ()
 Constructor.
 
 ~DeformableIcp ()=default
 Destructor.
 
Algorithm::Status compute (const matX3 &pointsFixed, matX3 &pointsMoving)
 Registers the moving points in-place to the fixed points.
 
void setNormals (const matX3 *normalsFixed)
 Sets normals for the fixed points. If set, point-to-plane instead of point-to-point ICP will be used.
 
void setRegularizationMatrix (const Eigen::SparseMatrix< double > *matrix)
 Sets the regularization matrix (typically a differential operator like the Laplacian).
 
void setRegularizationWeights (const vecX *weights)
 Sets regularization weights.
 
void setExtractionMatrix (const Eigen::SparseMatrix< double > *matrix)
 Sets a matrix that is applied to the moving points in all computations except the regularization term.
 
void setRegions (const vecXi *regionsFixed, const vecXi *regionsMoving)
 Assigns a region number to each fixed and each moving point.
 
void setNearestNeighborDirectionPerRegion (const std::map< int, NearestNeighborDirection > *directionPerRegion)
 Sets region-specific directions in which point correspondences are established using nearest neighbor search.
 
void setWeightPerRegion (const std::map< int, double > *weightPerRegion)
 Sets region-specific weights.
 
void setIterationCallback (std::function< bool(const IterationInfo &)> callback)
 Sets a function that is called at the end of each iteration.
 
void setProgressTask (Progress::Task *task)
 Sets the progress task to be used.
 
void setCustomTerm (std::function< const Eigen::SparseMatrix< double > &(const IterationInfo &, const matX3 &pointsMoving)> lhs, std::function< const matX3 &(const IterationInfo &, const matX3 &pointsMoving)> rhs)
 For power users only.
 
void setPreIterationFixedPointsProcessing (std::function< matX3(const IterationInfo &, const matX3 &pointsFixed, const matX3 &pointsMoving, const vecXi *idsMoving)> processing)
 For power users only.
 
IterationInfo lastIterationInfo () const
 Returns information about the most recent completed iteration.
 
- Public Member Functions inherited from Configurable
virtual void configure (const Properties *p)
 Configure this object instance by de-serializing the given Properties.
 
virtual void configuration (Properties *p) const
 Serialize the current object configuration into the given Properties object.
 
virtual void configureDefaults ()
 Retrieve the properties of this object, replaces values with their defaults and sets it again.
 
void registerParameter (ParameterBase *param)
 Register the given Parameter or SubProperty, so that it will be configured during configure()/configuration().
 
void unregisterParameter (const ParameterBase *param)
 Remove the given Parameter or SubProperty from the list of registered parameters.
 
 Configurable (const Configurable &rhs)
 
 Configurable (Configurable &&rhs) noexcept
 
Configurableoperator= (const Configurable &)
 
Configurableoperator= (Configurable &&) noexcept
 

Public Attributes

Parameter< bool > p_firstAlignCentroids = {"firstAlignCentroids", true, this}
 Whether to first translate the moving points so that their centroid aligns with the fixed points' centroid.
 
Parameter< PreAlignmentp_preAlignment = {"preAlignment", PreAlignment::RigidThenRigidAndScaling, this}
 Which pre-alignment to perform.
 
Parameter< bool > p_onlyPreAlignment = {"onlyPreAlignment", false, this}
 Whether to stop after pre-alignment.
 
Parameter< double > p_initLogLambda = {"initLogLambda", 2., this}
 Initial value of the regularization parameter's logarithm.
 
Parameter< double > p_lambdaReductionFactor = {"lambdaReductionFactor", 0.8, this}
 Factor smaller than 1 by which the regularization parameter (not its logarithm) is reduced each iteration.
 
Parameter< double > p_stoppingLogLambda = {"stoppingLogLambda", -10., this}
 If the regularization parameter's logarithm is less than or equal to this value, iterations are stopped.
 
Parameter< double > p_stoppingRelRmse = {"stoppingRelRmse", 1., this}
 If the relative RMSE is less than or equal to this value, iterations are stopped.
 
Parameter< double > p_stoppingRelRmseDiff = {"stoppingRelRmseDiff", 0.001, this}
 If the difference of relative RMSE between successive iterations is less than or equal to this value, pre-alignment iterations are stopped.
 
Parameter< int > p_maxIter = {"maxIter", 200, this}
 Maximum number of ICP iterations.
 
Parameter< NearestNeighborDirectionp_nearestNeighborDirection = {"nearestNeighborDirection", NearestNeighborDirection::WrtBoth, this}
 Direction in which point correspondences are established using nearest neighbor search.
 
Parameter< double > p_conjGradTolerance = {"conjGradTolerance", 1e-3, this}
 The tolerance threshold for stopping the conjugate gradient iterations used to solve the linear system (upper bound of relative residual error).
 
Parameter< int > p_conjGradMaxIterations = {"conjGradMaxIterations", 4, this}
 Conjugate gradient (CG) with a complete Cholesky preconditioner is used to solve the linear system.
 
Parameter< bool > p_updateWhileIterating = {"updateWhileIterating", false, this}
 Whether to update the display via the progress task while iterating.
 
Parameter< bool > p_verbose = {"verbose", false, this}
 Whether to print additional information (used regions, current phase, timings) to the log.
 
SubProperty< std::optional< RansacProperties > > p_preAlignmentRansac = {"preAlignmentRansac", std::nullopt, this}
 If set, the algorithm will use RANSAC to make the pre-alignment more robust.
 
- Public Attributes inherited from Configurable
Signal signalParametersChanged
 Emitted whenever one of the registered Parameters' or SubPropertys' signalValueChanged signal was emitted.
 

Additional Inherited Members

- Protected Attributes inherited from Configurable
std::vector< Paramm_params
 List of all registered Parameter and SubProperty instances.
 

Member Enumeration Documentation

◆ PreAlignment

enum class PreAlignment
strong

Pre-alignment modes.

Enumerator
None 

No pre-alignment.

Rigid 

ICP pre-alignment using a rigid transformation.

RigidThenRigidAndScaling 

Two subsequent runs of ICP pre-alignment:

  1. Using a rigid transformation;
  2. Using a rigid-and-scaling transformation.

◆ NearestNeighborDirection

enum class NearestNeighborDirection
strong

Direction in which point correspondences are established using nearest neighbor search.

Enumerator
WrtFixed 

For each point in fixed mesh, find nearest point in moving mesh.

WrtMoving 

For each point in moving mesh, find nearest point in fixed mesh.

WrtBoth 

For each point in both meshes, find nearest point in other mesh.

◆ StoppingReason

enum class StoppingReason
strong

Used to report why iterations have been stopped. Stopping criteria are evaluated in this order.

Enumerator
None 

None of the stopping criteria is fulfilled.

Lambda 

Logarithm of regularization parameter lambda is less than or equal to p_stoppingLogLambda.

RelRmse 

Relative RMSE is less than or equal to p_stoppingRelRmse.

RelRmseDiff 

Difference of relative RMSE between iterations is less than or equal to p_stoppingRelRmseDiff.

MaxIter 

Number of iterations is greater than or equal to p_maxIter.

Callback 

The user-defined callback has returned false.

◆ Phase

enum class Phase
strong

Used to report in which phase the algorithm is in.

Enumerator
None 

The algorithm has not started with any of the phases below.

Centroids 

The algorithm is aligning the centroids.

Rigid 

The algorithm is performing the rigid pre-alignment.

RigidAndScaling 

The algorithm is performing the rigid-and-scaling pre-alignment.

Deformable 

The algorithm is performing the deformable registration.

Member Function Documentation

◆ compute()

Algorithm::Status compute ( const matX3 & pointsFixed,
matX3 & pointsMoving )

Registers the moving points in-place to the fixed points.

If normals for the fixed points are provided, point-to-plane ICP is performed instead of point-to-point ICP.

◆ setRegularizationMatrix()

void setRegularizationMatrix ( const Eigen::SparseMatrix< double > * matrix)

Sets the regularization matrix (typically a differential operator like the Laplacian).

The number of columns of this matrix must be the same as the number of moving points.

◆ setRegularizationWeights()

void setRegularizationWeights ( const vecX * weights)

Sets regularization weights.

The number of weights must be the same as the number of moving points.

◆ setExtractionMatrix()

void setExtractionMatrix ( const Eigen::SparseMatrix< double > * matrix)

Sets a matrix that is applied to the moving points in all computations except the regularization term.

This matrix can be used to extract a subset of points and/or to linearly interpolate new points to be registered while still regularizing the original set of moving points. The number of columns of this matrix must be the same as the number of moving points.

◆ setRegions()

void setRegions ( const vecXi * regionsFixed,
const vecXi * regionsMoving )

Assigns a region number to each fixed and each moving point.

Each region is treated separately when establishing point correspondences. Negative region numbers can be used to define regions to be ignored.

◆ setNearestNeighborDirectionPerRegion()

void setNearestNeighborDirectionPerRegion ( const std::map< int, NearestNeighborDirection > * directionPerRegion)

Sets region-specific directions in which point correspondences are established using nearest neighbor search.

Keys of the map have to be region numbers. The map does not have to include all regions. p_nearestNeighborDirection is used for the remaining regions.

◆ setWeightPerRegion()

void setWeightPerRegion ( const std::map< int, double > * weightPerRegion)

Sets region-specific weights.

Keys of the map have to be region numbers. The map does not have to include all regions. A weight of 1 is used for the remaining regions.

◆ setIterationCallback()

void setIterationCallback ( std::function< bool(const IterationInfo &)> callback)

Sets a function that is called at the end of each iteration.

If this function returns false, the iterations are stopped.

◆ setCustomTerm()

void setCustomTerm ( std::function< const Eigen::SparseMatrix< double > &(const IterationInfo &, const matX3 &pointsMoving)> lhs,
std::function< const matX3 &(const IterationInfo &, const matX3 &pointsMoving)> rhs )

For power users only.

Sets a custom term to add to the left- and right-hand side of the normal equation. pointsMoving will contain all moving points before applying the extraction matrix.

◆ setPreIterationFixedPointsProcessing()

void setPreIterationFixedPointsProcessing ( std::function< matX3(const IterationInfo &, const matX3 &pointsFixed, const matX3 &pointsMoving, const vecXi *idsMoving)> processing)

For power users only.

Sets a function that is applied to the fixed points before each iteration. pointsMoving will contain all moving points after applying the extraction matrix. If regions are used, the function is called for each region, pointsFixed will contain only the fixed points for this region, idsMoving will contain the indices of pointsMoving corresponding to pointsFixed and the function must return only the processed fixed points for this region.

Member Data Documentation

◆ p_lambdaReductionFactor

Parameter<double> p_lambdaReductionFactor = {"lambdaReductionFactor", 0.8, this}

Factor smaller than 1 by which the regularization parameter (not its logarithm) is reduced each iteration.

A larger factor should yield better point correspondences, but also increases the computation time. 0.8 is usually a good compromise.

◆ p_stoppingRelRmse

Parameter<double> p_stoppingRelRmse = {"stoppingRelRmse", 1., this}

If the relative RMSE is less than or equal to this value, iterations are stopped.

The relative RMSE is the absolute RMSE divided by the RMS distance of fixed points from their centroid.

◆ p_stoppingRelRmseDiff

Parameter<double> p_stoppingRelRmseDiff = {"stoppingRelRmseDiff", 0.001, this}

If the difference of relative RMSE between successive iterations is less than or equal to this value, pre-alignment iterations are stopped.

This criterion is not used for the deformable registration.

◆ p_nearestNeighborDirection

Parameter<NearestNeighborDirection> p_nearestNeighborDirection = {"nearestNeighborDirection", NearestNeighborDirection::WrtBoth, this}

Direction in which point correspondences are established using nearest neighbor search.

If regions are used and individual nearest neighbor directions have been set for some but not all regions, this value is still used for the other regions.

◆ p_conjGradTolerance

Parameter<double> p_conjGradTolerance = {"conjGradTolerance", 1e-3, this}

The tolerance threshold for stopping the conjugate gradient iterations used to solve the linear system (upper bound of relative residual error).

If set to 0, the system is solved directly using a Cholesky factorization.

◆ p_conjGradMaxIterations

Parameter<int> p_conjGradMaxIterations = {"conjGradMaxIterations", 4, this}

Conjugate gradient (CG) with a complete Cholesky preconditioner is used to solve the linear system.

The Cholesky factorization is reused across ICP iterations and is only recomputed if the number of CG iterations needed during the previous ICP iteration is greater than or equal to this value. If set to 0, the system is solved directly using the Cholesky factorization.


The documentation for this class was generated from the following file:
Search Tab / S to search, Esc to close