Rigid Shape Registration Based on Extended Hamiltonian Learning.

Shape registration, finding the correct alignment of two sets of data, plays a significant role in computer vision such as objection recognition and image analysis. The iterative closest point (ICP) algorithm is one of well known and widely used algorithms in this area. The main purpose of this paper is to incorporate ICP with the fast convergent extended Hamiltonian learning (EHL), so called EHL-ICP algorithm, to perform planar and spatial rigid shape registration. By treating the registration error as the potential for the extended Hamiltonian system, the rigid shape registration is modelled as an optimization problem on the special Euclidean group S E ( n ) ( n = 2 , 3 ) . Our method is robust to initial values and parameters. Compared with some state-of-art methods, our approach shows better efficiency and accuracy by simulation experiments.

1. Introduction

Point set registration is to find an optimal transformation that aligns one set to the other. This problem occurs in many applications, such as motion tracking [1,2], target localization [3], super-resolution [4], mosaicing [5] and medical image analysis [6,7,8,9,10]. Depending on different goals, the methods can be categorized as coarse registration and fine registration [11,12]. The former is to find an initial estimation between two point sets while the latter is to obtain a more accurate solution. In general, the popular method used in coarse registration is point-to-point and iterative, such as point signature and spin image. Some other methods may include principle component analysis (PCA), algebraic surface model and principle curvature. While methods like iterative closed point (ICP), Chen’s method, signed distance fields and genetic algorithms, are more commonly used in fine registration.

Literately, the ICP algorithm is efficient and accurate among all point set registration techniques. However, it requires that the initial configuration of the two point sets is sufficiently close and it is sensitive to noise and outliers. Moreover, the ICP algorithm is unable to handle affine case. To improve the robustness and cope with more situations, many researchers extended the ICP method. In particular, Ying et al. applied the theory of differential geometry and Lie group to ICP algorithm and formed a unified framework to design fast and robust shape registration algorithms [13,14,15,16,17,18,19,20,21].

From the geometric point of view, given the correspondence between two sets of points, the affine registration is to find an element in the general linear group that transfers one set to the other optimally. In this paper, we deal with the rigid registration, which uses the geometry of the special Euclidean group. From this perspective, differential geometry is a powerful tool not only in registration problems but also in computer vision extensively. There is a vast community dealing with shape registration/analysis/comparison from geometrical viewpoints. For example, E. Celledoni etc. applied the theory of Lie groups and homogeneous manifolds to the problem of shape analysis [22,23]. A.L. Brigant introduced an algorithm that finds an optimal matching between two curves by computing the geodesic of the infinite-dimensional manifold of curves [24]. X. Pennec etc. concentrated on feature-based approaches for rigid registrations using differential geometry of surfaces [25]. Other methods, using harmonic analysis and statistical optimizations, can be found in [26,27,28,29]. Experiments showed that these geometry-based approaches performed better than other state-of-art methods.

There is also an immense literature in applying Hamiltonian dynamical systems to learning theory. F. Barbaresco studied the symplectic extension of Souriau Lie groups thermodynamics and used this model for data analysis and machine learning on Lie groups [30,31]. S. Fiori proposed extended Hamiltonian learning (EHL) on Riemannian manifolds motivated by Hamilton variational principle [32,33]. Compared with the classical gradient-based learning, EHL has the advantages of averaging out the oscillations and mitigating the plateau effect. In this paper, we investigate the 2D/3D rigid registration problem from the view of Hamiltonian learning. The innovation of our proposed method is treating the registration error as the potential of the extended Hamiltonian system on the underlying space of transformation, the special Euclidean group. Under the system of extended Hamiltonian, an EHL-based algorithm is designed to achieve the optimal transformation corresponding to the minimum of the registration error. Some experiments are carried out to validate the efficiency and robustness of our algorithm.

The structure of the paper is organized as follows. In Section 2, some basic background about geometric structures of the special Euclidean group and is reviewed, including the Riemannian metric, geodesic, and Riemannian gradient. In Section 3, extended Hamiltonian learning on a general Riemannian manifold is discussed. And in the following section, we formulate the method for 2D/3D rigid registration problem by extended Hamiltonian learning on and . In Section 5, some numerical experiments prove the efficiency of this proposed method, comparing with Du’s SVD method (SVD) [16], Ying’s Iwasawa decomposition method (ID) [21], Ying’s Lie group optimization method (LGO) [18]. At last, conclusions are presented and possible improvements of our method for future work are discussed.

2. Geometry of Special Euclidean Groups

In this section, we review some basics about the special Euclidean group , which is the underlying geometric space for rigid registration. The special Euclidean group consists of rotations and translations in the Euclidean space. It is the semi-direct product of the special orthogonal group and ,

Being a Lie group, is equipped with the smooth structure and group structure simultaneously.

Any element g in can be represented as a pair , where and stands for the rotation and translation respectively. In the matrix form, g can be written as which is an one-to-one correspondence with the pair . The action of the special group on the Euclidean space is defined as

The Lie algebra of , denoted by , contains infinitesimal rotations and translations. A general element h in is a pair , where J is a skew-symmetric matrix and v is a vector in . In matrix form, h can be rewritten as

The exponential mapping and the logarithm mapping provide methods to transfer information between the nonlinear manifold and the linear space . Concretely, the exponential map of a general Lie algebra element is defined by

In low dimensions the exponential map can be written in a compact form, as presented by J.M. Selig [34]. Denote , where

For any , direct computation shows that it satisfies the quartic equation . Thus the higher order terms in the taylor expansion (5) can be simplified. Explicitly, we have

Especially, for pure rotations, where in the expression (4) and h satisfies the cubic equation , the exponential map can be written as while for pure translations, where in the expression (4) and , the exponential map satisfies

The exponential map for pure rotation degenerates on to be as the following,

An inner product on the Lie algebra can be extended to a Riemannian metric on the Lie group G by left translation. Specifically, the inner product at is defined by where is a fixed parameter, . The left-invariant metric, for two tangent vectors and at an arbitrary group element is defined as

The right-invariant metric is defined similarly. In physics, the left-invariance and right-invariance correspond to the independence of the choices of the inertial frame and the body-fixed frame respectively [35]. A bi-invariant metric means it is both left-invariant and right-invariant. However, there is no bi-invariant metric on [36]. Thus, we adopt the left-invariant metric on throughout this paper.

Let f be a function defined on Riemannian manifold M, and be Euclidean gradient and Riemannian gradient, respectively. The relation between and is governed by [32] where and denote Euclidean inner product and Riemannian metric respectively. For the Riemannian metric defined in (12) on , we have

3. Extended Hamiltonian Learning on SE(n)

We first introduce the extended Hamilton learning on general Riemannian manifolds [32,33]. Then we specify on special Euclidean group .

On a Riemannian manifold M, equipped with a metric at , the extended Hamiltonian principle is

In Equation (15), stands for the trajectory of a particle moving along M and is the corresponding instantaneous velocity. The function denotes the kinetic energy and V represents associated potential energy. The symbol represents the variation of the action in the dynamical system, namely, the particle slides from a point to an infinitely close point. An element indicates the dissipation force at the point x. This system degenerates to conservative if the dissipation force disappears everywhere. Following [32,33], the kinetic energy is adopted as the symmetric bilinear form under the assumption of unit mass. From viscosity theory the dissipation term is assumed to be with in the paper. Then, the equations for such dynamical system read where is the Christoffel symmetric form, and denotes the Riemannian gradient of V.

In order to implement simulations, we turn the continuous version (16) into a discrete one, which yields where is the exponential map and is the selected step size.

When it comes to the special Euclidean group , we need to compute the Christoffel symmetric form. However, in the discrete case, it suffices to know the form in a subspace of . Specifically, the Christoffel symmetric form for , reads

Then Equation (17) can be rewritten as

Note that from the iteration (19) one cannot ensure that is skew-symmetric. Consequently will not be an element in . Thus, should be modified to keep the validity of the iteration. To do this we make an orthogonal projection from the Euclidean space to the subspace consisting of skew-symmetric matrices. i.e., For any matrix , it can be decomposed as the sum of a skew-symmetric part and a symmetric part . From (19), the symmetric part would be a negligible error when the step is small enough. Thus we apply the skew-symmetric part in the tangent space . See Figure 1. We remark that this step is nothing but the to compute the covariate derivative on Euclidean submanifolds by definition [37]. Hence, the iteration becomes

Figure 1

Iteration with decomposition of .

For details about the convergence of extended Hamiltonian learning, refer to [32]. The necessary condition for convergence is that satisfies , where is the maximum eigenvalue of the Hessian matrix of V.

4. 2D/3D Rigid Shape Registration Based on Extended Hamiltonian Learning

Given two n-D (n = 2,3) data sets and and a correspondence , denoting , find an element such that the cost function achieves it minimum.

The fundamental steps based on ICP for rigid registration are: First, for the current fixed transformation , find with elements such that the subset minimizes

Second, once the correspondent data set is obtained, update the transformation as for minimizing

In fact, the translation t can eliminated by coordinating the centers of the two data sets. Let and be centers of X and Y, respectively. For the obtained set , the least squares solution of (23) indicates

Thus, with the centralized data , optimization problems (22) and (23) can be simplified as

Here, we regard the registration error as potential of extended Hamiltonian system on the special Euclidean group . The Euclidean gradient of is given as from which we can compute the Riemannian gradient of by (14).

Therefore, for two given data sets and , we summarize the method based on extended Hamiltonian learning as Algorithm 1.

Initialize parameters , , ;

Set , ;

Centralize as , , and set ;

fordo

Search for by minimizing ;

Calculate the registration error ;

Search for by minimizing by Equation (19);

;

if then

return

else

;

end if

end for

5. Numerical Results

All data samples used in this section is from the MPEG-7 shape B database and all programs are written in Matlab 2018a and run by PC with AMD Athlon II P340 Dual-Core processor, 2.20 GHz CPU, and 2 GB RAM.

5.1. 2D Rigid Shape Registration

In 2D case, our method appears to be robust and insensitive to initial values. Moreover, a near-optimal registration can be obtained within a few steps. To give a visualization of our method, we test chicken-2 and chicken-3 in MPEG-7 shape B database as model data set and test data set. The initial rotation is set to be the identity and the initial translation is set to be the difference of two means of data. The numerical results are shown in Figure 2.

Figure 2

The (top left) figure is the model data, the (top right) figure is the test data, the (bottom left) is the figure after five iterations and the (bottom right) is the figure after final registration.

We select nine groups of rigid data in MPEG-7 shape B database to run the experiments and compare with Du’s SVD method, Ying’s Iwasawa decomposition method and Ying’s Lie group optimization method. To give a quantitive comparison we define the root mean square (RMS) error to be for model data and test data . Though the geodesic distance on is more reasonable in theory, practically it is difficult to compute since we do not know the true rigid transformation.

The precision is set to be . Comparison results are displayed in Figure 3. The resulted RMS errors are displayed in Table 1. Other methods require the careful choice of initial values, whereas for our method we simply choose identity as the initial rotation and difference of means as initial translation. We can find that EHL-ICP algorithm is more robust to the size and shape of point cloud data. Other algorithms may have good performance on small point sets with simple shapes but lose precision when point sets are complex.

Figure 3

Experiments using different data sets. The test data are colored green; the model data are colored blue; the registration results are colored red.

Table 1

Comparison of the RMS errors for 2D rigid shape registration. The optimal results are in bold.

Group	Model	Test	SVD	ID	LGO	EHL-ICP
(1)	bird-3	bird-4	0.5841	0.9996	0.5690	0.4048
(2)	deer-1	deer-4	0.5263	2.7598	0.5272	2.8826
(3)	horse-3	horse-4	0.5107	1.4278	0.5112	0.3880
(4)	beetle-7	beetle-8	0.8749	0.8746	0.5242	0.4730
(5)	cattle-1	cattle-20	22.8580	28.3719	22.5155	1.1656
(6)	hammer-4	hammer-5	0.4846	0.8037	0.4232	0.3043
(7)	chicken-2	chicken-3	0.5484	2.7 843	0.5471	0.5202
(8)	butterfly-1	butterfly-2	18.1726	34.3588	6.9691	2.9062
(9)	horseshoe-9	horseshoe-17	0.5425	0.5690	0.5873	0.3577

5.2. 3D Shape Registration

Similar to the 2D case, we selected a group of 3D models including bunny, chair, cactus, dinosaur, elephant and block to verify the validity of our algorithm. The initial rotation is set to be identity and the initial translation is set to be the difference of means of model data and test data. The visualized results are represented in Figure 4. Note that we do not require a subtle choice of initial values and parameters. The results demonstrate the efficacy of our method.

Figure 4

Blue: Model data (fixed); Green: Test data (moving); Red: Final data (registration). Figures of Stanford Bunny, dinosaur, block, chair, cactus, elephant.

6. Conclusions and Future Works

Shape registration plays a significant role in computer vision, where the task is to transfer one set of points to the other. Since the iterative closest point method is widely used in registration problems yet having some shortcomings, this paper proposes the EHL-ICP method, which incorporating the extended Hamiltonian learning with the ICP algorithm, to deal with the 2D and 3D rigid shape registration problem. By regarding the registration error as the potential of the extended Hamiltonian system, we formulate rigid registration as an optimization problem on the special Euclidean group . Numerical results show that our method is more effective and accurate when compared with other methods. Moreover, our method is robust with respect to initial values in both dimensions, which provides a good choice for rough registration.

For future work, we may generalize the extended Hamiltonian learning method to different registration problems. There are two hot topics worth mentioning. The first is to use optimal transportation for data set registration [38,39]. The basic idea is that shape data can be viewed as a sum of Dirac measures in a given space and the difference of shapes is taken to be the Wasserstein distance. To find a (non)rigid transformation is to find a measure-preserving map. Another possible application is affine registration where shapes are distorted and there is no rigid transformation [40,41]. From this viewpoint, we should consider extended Hamiltonian learning on the general linear group , where more techniques should be developed.