1 Introduction
The alignment of a set of objects by means of transformations plays an important role in the field of computer vision and recognition. For instance, for the creation of statistical shape models (SSMs) [5] training shapes are initially aligned for removing pose differences in order to only model shape variability.
The most common way of shape representation is by encoding each shape as a pointcloud. In order to be able to process a set of shapes it is necessary that correspondences between all shapes are established. Whilst there is a vast amount of research in the field of shape correspondences (for an overview see [11, 18]), in this paper we focus on the alignment of shapes and we assume that correspondences have already been established.
The alignment of two objects by removing location, scale and rotation is known as Absolute Orientation Problem (AOP) [13] or Procrustes Analysis [8]
. For the AOP there are various closedform solutions, among them methods based on singular value decomposition (SVD)
[1, 16]; based on eigenvalue decomposition
[13]; based on unit quaternions [12] or based on dual quaternions [19]. A comparison of these methods [6] has revealed that the accuracy and the robustness of all methods are comparable.The alignment of more than two objects is known as Generalised Procrustes Analysis (GPA). Whilst a computationally expensive global solution for GPA in two and three dimensions has been presented in [15]
, the most common way for solving the GPA is to align the objects with a reference object. However, fixing any of the objects as reference induces a bias. An unbiased alternative is to align all objects with the adaptive mean object as reference. An iterative algorithm then alternatingly updates the reference object and estimates the transformations aligning the objects. The iterative nature of these methods constitutes a problem if the relative transformation between any pair of objects is noisy. This is for example the case if data is missing, correspondences are wrong or if the transformations are observed by independent sensors (
e.g. noncommunicating robots observe each other). Noisy relative transformations can be characterised by transitive inconsistency, i.e. transforming to and to might lead to a different result than transforming directly to .This paper presents a novel method for synchronising the set of all pairwise transformations in such a way that they globally exhibit transitive consistency. Experiments demonstrate the effectiveness of this method in denoising noisy pairwise transformations. Furthermore, using this novel method the GPA is solved in an unbiased manner in closedform, i.e. noniterative. Transformation synchronisation is applied to solve the GPA with missing data as well as with wrong correspondence assignments and results in superior performance compared to existing methods.
Our main contribution is a generalisation of the techniques presented by Singer et al. [3, 9, 10, 17], who have introduced a method for minimising global selfconsistency errors between pairwise orthogonal transformations based on eigenvalue decomposition and semidefinite programming. With permutation transformations being a subset of orthogonal transformations, in [14] the authors demonstrate that the method by Singer et al. is also able to effectively synchronise permutation transformations for globally consistent matchings.
In our case, rather than considering the special case of orthogonal matrices, we present a synchronisation method for invertible linear transformations. Furthermore, it is demonstrated how this method can be applied for the synchronisation of similarity, euclidean and rigid transformations, which are of special interest for the groupwise alignment of shapes.
Whilst the proposed synchronisation method is applicable in many other fields where noisy pairwise transformations are to be denoised (e.g. groupwise image registration or multiview registration), in this paper GPA is used as illustrating example.
2 Methods
For the presentation of our novel transformation synchronisation method the notation and some foundations are introduced first. Subsequently, a formulation for the case of perfect information is given. Motivated by these elaborations, a straightforward extension to handle noisy pairwise transformations is presented. Finally, various types of transformations are discussed.
2.1 Notation and Foundations
are matrices representing pointclouds with points in dimensions where in the following all are simply referred to as pointclouds. Let
be the identity matrix and
be the vector containing only zeros, both having appropriate dimensions according to their context. The Frobenius norm is denoted by
. Let be an invertible transformation matrix aligning pointcloud with (for all ), where . Furthermore, is the set of all pairwise transformations.A desirable property of the set of transformations is that it complies with the following transitive consistency condition:
Definition 1.
The set of relative transformations is said to be transitively consistent if
Definition 1 states that the transformation from to followed by the transformation from to must be the same as directly transforming from to .
Lemma 1.
The set of relative transformations is transitively consistent if and only if there is a set of invertible transformations such that
Proof.
For the sake of completeness a proof is provided.
“”: Transitive consistency of follows directly from the definition , since for all it holds that
(1)  
(2)  
(3)  
(4) 
“”: In the rest of the proof we direct our attention towards the necessity of the existence of the transformations.
First of all, if the transformations in are transitively consistent
(5) 
This follows by the fact that the needs to be invertible while satisfying, by Definition 1, that .
Let for all . Now we show that are such matrices we seek. Since, by using (5), , we can write .
Now for any , we can use that
(6) 
Thus,
(7) 
∎
2.2 Perfect Information
Due to Lemma 1, there is a reference coordinate frame, denoted by , from which there are transformations such that for all . Note that the reference coordinate frame is merely used as a tool for deriving our method and it is irrelevant what the actual reference frame is. Let us introduce
(8)  
(9)  
(10)  
(11) 
where
Using this notation, finding either or (up to an invertible linear transformation) gives the transitively consistent transformations.
Definition 2.
Let . The set
is the column space of and the set
is the null space of .
Note that due to the invertability of (for all ) it holds that the matrix has rank and thus the dimensionality of the column space of is exactly .
Proposition 1.
Let . The linear subspace has dimension and is equal to .
Proof.
First it is shown that the columns of are contained in the null space of , i.e. , and then it is shown that the null space of has exactly dimension .
Note that , which we will make use of shortly. Multiplication of on the right to gives
(12)  
(13)  
(14)  
(15)  
(16) 
From (16) it can be seen that all the columns of are contained in the null space of , so .
However, it still remains to be shown that the dimensionality of is exactly , i.e. spans the entire null space of and not just a part of it. This is done by showing that there are no nonzero vectors that are not contained in but are contained in .
Formally this is expressed by the requirement that the set
is empty. Suppose now that is not empty, so it contains the element . Using the orthogonal decomposition theorem, the vector can be rewritten as , where and . The definition of states that , which implies that . Further, the definition of states that
(17)  
(18)  
(19)  
(20)  
Per definition , so it follows that , leading to  
(21)  
Multiplication of on the left gives  
(22)  
(23)  
(24)  
(25)  
Per definition , so , leading to  
(26)  
(27) 
Equation (27) is a contradiction to , thus, the set is empty. ∎
Proposition 1 states that in (16) can, up to an invertible linear transformation, be retrieved by finding the dimensional null space of . Let be the singular value decomposition (SVD) of . The columns of corresponding to the zero singular values span and give a solution to (16).
As we are only able to retrieve the transformations in the blocks of up to invertible linear transformations, w.l.o.g. we create a new version of , call it , with the first block being equal to the identity, as
(28) 
2.3 Noisy Pairwise Transformations
Up until this point, the matrix is obtained under perfect information, i.e. the transitivity condition in Definition 1 holds for all transformations contained in the blocks of . However, we are interested in the case when the transitivity condition does not hold due to measurement noise. Assume now that we have a noisy observation of , denoted as . Also, let the noisy version of be . Now, in general it is not the case that the null space of is dimensional. Instead, the leastsquares approximation of the dimensional null space is considered, which leads to the following optimisation problem:
Problem 1.
Leastsquares Transformation Synchronisation
subject to  
The rank approximation of the null space of can be retrieved using the SVD of . In this case the columns of corresponding to the smallest singular values span the rank approximation of the null space of , giving , the estimate for . By using (28), can be retrieved from .
2.4 Affine Transformations in Homogeneous Coordinates
In this section it is shown that the method is also applicable for invertible affine transformations, rather than invertible linear transformations. This is done by representing the dimensional affine transformations by using homogeneous matrices.
Each affine transformation can be written as
(29) 
where is the (invertible) linear transformation matrix and is the dimensional row vector representing the translation. The inverse of is given by
(30) 
Similar to the linear case described in (8), the matrix is constructed from all . It is assumed that the matrix , corresponding to the noisy observation of , contains blocks that are proper affine transformations, i.e. the last column of each block is .
A simple way to ensure that the synchronised transformations are affine transformations in homogeneous coordinates is to add the row vector , with , to the matrix . By adding the vector to , the vector is removed from the null space of . Using this approach, a solution is then found by solving Problem (1) with the updated . Then the resulting gives an estimate of the first columns of and these are the columns we seek.
2.5 Similarity Transformations
Similarity transformations are transformations that allow for translations, isotropic scaling, rotations and reflections. To retrieve similarity transformations, the estimates of the synchronised affine transformations are determined first. The translation component of can directly be extracted from since it has the structure presented in (29). To obtain the scaling factor and the orthogonal transformation, the linear component is factorised using SVD, resulting in . The orthogonal component is then given by
(31) 
and the isotropic scaling factor is given by
(32) 
where is the th element on the diagonal of .
Remark 1.
It can be shown that retrieving the orthogonal component as presented in eq. (31) is the leastsquares solution to the projection onto the set of orthogonal matrices. However, in eq. (32
) the isotropic scaling factor is retrieved as the geometric mean of the individual axisaligned scaling factors. The leastsquares solution to the projection onto the set of similarity transformations is given by the arithmetic mean,
i.e. .2.6 Euclidean Transformations
Similarity transformations without isotropic scaling are called euclidean transformations. To obtain euclidean transformations, the similarity transformations are extracted and the scaling factors (for all ) are set to .
2.7 Rigid Transformations
Euclidean transformations without reflections are called rigid transformations. Rigid transformations can be obtained by ensuring that the determinant of the rotational component described in (31) equals . This can be achieved by setting
(33)  
(34) 
3 Experiments
By generating ground truth data and adding Gaussian noise to it, we first compare the error of the synchronised transformations using our method to the error of the unsynchronised transformations. Furthermore, the transformation synchronisation method is applied for solving the Generalised Procrustes Problem with missing points and with wrong correspondence assignments.
3.1 Noisy Transformations
In this section it is described how the ground truth transformations are generated, how noisy versions thereof are generated and eventually results of the transformation synchronisation method are presented.
3.1.1 Ground Truth Transformations
For the analysis of the performance of our method we generate a set of random transformations , that are used in turn to generate the transitively consistent set of pairwise transformations , serving as ground truth for the evaluation. The generation of is described in the following.
The dotnotation is used to illustrate that
is a random variable with a particular probability distribution. For generating the set
, we assume that the pointclouds that lead to the transformations have some structural similarity, i.e. the transformations are not entirely random. In particular, the scaling factors, the translation components and the linear part of the transformation are restricted in the sense that they cannot be arbitrary. However, arbitrary orientations in dimensional space are allowed for.The set contains the elements , which are samples of
(35) 
where is a scaling factor and is a translation, with denoting the
dimensional uniform distribution having the open interval
as support. Samples of the random rotation matrixare drawn by extracting the rotational component of a nonsingular random matrix as described in (
31). The random noise matrix is given by , where is arandom matrix with each element having univariate normal distribution
. The purpose of creating the noise in the way using the random matrix is to restrict the linear component in the transformation and thus to avoid illconditionedness with very high probability.Depending on the type of transformation that is evaluated, the parameters of have different properties, which are summarised in Table 1.
linear  
affine  
similarity  
euclidean  =1  
rigid  =1 
Once the ground truth set of transitively consistent transformations has been established, a noisy version thereof is synthetically created, as described in the next section.
The error between two sets of pairwise transformations and is defined as
(36) 
3.1.2 Additive Gaussian Noise

The set of noisy pairwise transformations is created by adding to each element of the matrix a sample from , which is conducted for all matrices with . In the case of homogeneous transformation matrices , no noise is added to the last column, which shall always be .
Results of the simulations are shown in Fig. 1. The first row of graphs show that for all types of transformations the error of synchronised transformations is smaller than the error of the unsynchronised transformations and that the slope of the error in the synchronised case is smaller than in the unsynchronised case. In the second row it can be seen that, even with a high amount of noise (), the error of the synchronised transformations decreases with an increasing number of objects . As anticipated, with increasing there is more information available, directly resulting in a lower error. The last row of graphs shows that increasing the dimensionality results in an increasing error; however, the error of the synchronised transformations increases slower than for the unsynchronised ones.
3.2 Generalised Procrustes Analysis

In addition to evaluating the synchronisation of noisy pairwise transformations we have applied our method for solving the Generalised Procrustes Problem (GPP), which is done on the one hand with missing data and on the other hand with wrong correspondence assignments. For both simulations the 2D fish shapes from the ChuiRangarajan data set [4] with different levels of deformation and noise have been used (refer [4] for more details). For each level of deformation and noise the data set contains shapes, each comprising points in dimensions.
Finding the similarity transformation that best aligns two shapes, which is a subroutine for the evaluated referencebased, the iterative mean shapebased and the synchronisationbased method, is performed by an AOP implementation with symmetric scaling factors [13]. In the referencebased solution of GPP one shape is randomly selected as reference and all other shapes are aligned with the reference. For the iterative mean shapebased method the initial reference shape is selected randomly and then the mean shape is iteratively updated. In the synchronisationbased solution of GPP all pairwise AOPs are solved first, followed by the synchronisation of the resulting transformations in order to aggregate all information contained in the pairwise transformations. Additionally, the stratified GPA method proposed in [2]
is evaluated for solving the GPP. In our experiments we have observed that by using the stratified GPA method the linear part of the resulting transformations may collapse to the zero matrix; in order to enable a comparison with the other methods in these cases the linear part of the transformation has simply been set to the identity matrix.
In the missing data experiments as well as the wrong correspondence experiments for each single run out of shapes are randomly selected. For the experiments in the missing points case the missing points are simulated by discarding points according to a given probability. For the experiments with wrong correspondences the correct correspondences are randomly disturbed in order to simulate wrong correspondences.
In contrast to solving the AOPs, in both experiments the computation of the error is performed using the original shape (i.e. with all points and with perfect correspondences). With that we investigate up to which amount recovering the original shapes from corrupt shape data is possible. The average shape error of a set of shapes is defined as .
3.2.1 Missing Points
In every run, additionally to randomly selecting out of shapes, each data point of a shape is considered to be missing with probability . As the implemented methods solve the AOP only for common points in each pair of shapes, values for larger than have not been investigated because with the cases that the number of common points in a pair of shapes is less than occur too frequently (for dimensional data, there must be at least points in each shape in order to result in a system that is not underdetermined). Also, for it is possible that the number of common points in a pair of shapes is less than ; in these cases the draw of missing data is simply repeated.
In Fig. 2 the resulting error of the referencebased, the iterative mean shapebased, the synchronisationbased and the stratified solution of the GPP with missing data are shown for different levels of deformation and noise. It can be seen that even with an increasing amount of missing data, when using the synchronisationbased method the error increases only slightly, whilst the error of the referencebased method increases significantly with a larger amount of missing points. With respect to the error, the transformation synchronisation method performs only marginally better than the iterative mean shapebased method and the stratified method. However, the average runtimes for solving a single GPP instance was s for the referencebased method, s for the synchronisationbased method, s for the iterative mean shapebased method and s for the stratified method, illustrating that our method performs significantly better than all other methods when taking runtime and error into account at the same time.
3.2.2 Wrong Correspondence Assignments
Additionally to the case of missing points, we have applied our method to solve the GPP with wrong correspondence assignments between shapes. In order to mimic practical applications, where it is frequently the case that the true correspondences are unknown and thus it must be assumed that wrong correspondences are present, we do not make any efforts to correct these wrong correspondences (such as using RANSAC [7] or permutation synchronisation [14]). Instead, for each pair of shapes the AOP is solved whilst being aware that some of the points in the one shape have wrong counterparts in the other shape. Of course this will have influence on the resulting transformations. Thus, the objective of the simulations described in this section is to assess to what extent the transformations from shapes with wrong correspondences can be reconstructed using transformation synchronisation.
In every run, additionally to randomly selecting out of shapes, the correspondences between the points in each shape are disturbed. For disturbing the correspondence assignments each pair of shapes that is to be aligned is considered independently. For that, a proportion of points from the total number of points is selected. Then, as correspondences between the pair of pointclouds are implicitly given by the ordering of the rows, the rows corresponding to the previously selected points are reordered randomly in one of the pointclouds, directly resulting in disturbed correspondence assignments between the pair of pointclouds .
In Fig. 3 the referencebased, the iterative mean shapebased, the synchronisationbased and the stratified solution of the GPP with wrong correspondences are shown for different levels of deformation and noise. On the right of Fig. 3 examples of the correspondences between pairs of shapes are depicted for different values of .
It can be seen that for different levels of deformation and different levels of noise with of wrong correspondences the outcome is only marginally affected when using our proposed method. In contrast, all other evaluated methods result in significantly larger errors, which can be explained by the fact that our method is the only one that is able to make use of the information that is contained in all pairwise transformations.
4 Conclusion
The alignment of multiple (corresponding) pointclouds simultaneously is generally tackled by iteratively aligning all pointclouds to some reference. Whereas this approach is biased (selecting a fixed reference) or initialisationdependent (using the adaptive mean as reference) we have presented a method that is completely unbiased and does not depend on initialisation.
Our key observation is that the underlying noisefree transformations can be retrieved from the null space of a matrix that can directly be obtained from pairwise alignments. Whilst related approaches for rotation matrices [9, 10, 17] or permutation matrices [14] have been proposed, we have generalised the synchronisation method to handle general linear and affine transformations as well as similarity, euclidean and rigid transformations. Experimentally we were able to demonstrate that the proposed method is able to effectively reduce noise from the set of pairwise transformations and to solve the Generalised Procrustes Problem at least as good as existing approaches for the missing data case whilst significantly outperforming other methods for the presented wrong correspondence case.
Acknowledgements
Supported by the Fonds National de la Recherche, Luxembourg (5748689, 6538106, 8864515).
References
 [1] K. S. Arun, T. S. Huang, and S. D. Blostein. Leastsquares fitting of two 3D point sets. Pattern Analysis and Machine Intelligence, IEEE Transactions on, (5):698–700, 1987.
 [2] A. Bartoli, D. Pizarro, and M. Loog. Stratified generalized procrustes analysis. International Journal of Computer Vision, 101(2):227–253, 2013.
 [3] K. N. Chaudhury, Y. Khoo, and A. Singer. Global registration of multiple point clouds using semidefinite programming. arXiv.org, June 2013.
 [4] H. Chui and A. Rangarajan. A new point matching algorithm for nonrigid registration. Computer Vision and Image Understanding, 89(2):114–141, 2003.
 [5] T. F. Cootes and C. J. Taylor. Active Shape Models  Smart Snakes. In In British Machine Vision Conference, pages 266–275. SpringerVerlag, 1992.
 [6] D. W. Eggert, A. Lorusso, and R. B. Fisher. Estimating 3D rigid body transformations: a comparison of four major algorithms. Machine Vision and Applications, 9(56):272–290, Mar. 1997.
 [7] M. A. Fischler and R. C. Bolles. Random sample consensus: a paradigm for model fitting with applications to image analysis and automated cartography. Communications of the ACM, 24(6), June 1981.
 [8] J. C. Gower and G. B. Dijksterhuis. Procrustes problems, volume 3. Oxford University Press Oxford, 2004.
 [9] R. Hadani and A. Singer. Representation theoretic patterns in three dimensional CryoElectron Microscopy I: The intrinsic reconstitution algorithm. Annals of mathematics, 174(2):1219, 2011.
 [10] R. Hadani and A. Singer. Representation Theoretic Patterns in ThreeDimensional CryoElectron Microscopy II—The Class Averaging Problem. Foundations of computational mathematics (New York, N.Y.), 11(5):589–616–616, 2011.
 [11] T. Heimann and H.P. Meinzer. Statistical shape models for 3D medical image segmentation: A review. Medical Image Analysis, 13(4):543–563, 2009.
 [12] B. K. P. Horn. Closedform solution of absolute orientation using unit quaternions. Journal of the Optical Society of America A, 4(4):629–642, 1987.
 [13] B. K. P. Horn, H. M. Hilden, and S. Negahdaripour. Closedform solution of absolute orientation using orthonormal matrices. Journal of the Optical Society of America A, 5(7):1127, 1988.
 [14] D. Pachauri, R. Kondor, and V. Singh. Solving the multiway matching problem by permutation synchronization. In Advances in neural information processing systems, pages 1860–1868, 2013.

[15]
D. Pizarro and A. Bartoli.
Global optimization for optimal generalized procrustes analysis.
In
Computer Vision and Pattern Recognition (CVPR), 2011 IEEE Conference on
, pages 2409–2415. IEEE, 2011.  [16] P. H. Schönemann. A generalized solution of the orthogonal procrustes problem. Psychometrika, 31(1):1–10, Mar. 1966.

[17]
A. Singer and Y. Shkolnisky.
ThreeDimensional Structure Determination from Common Lines in CryoEM by Eigenvectors and Semidefinite Programming.
SIAM journal on imaging sciences, 4(2):543–572, June 2011.  [18] O. Van Kaick, H. Zhang, G. Hamarneh, and D. Cohen Or. A survey on shape correspondence. In Computer Graphics Forum, pages 1681–1707. Wiley Online Library, 2011.
 [19] M. W. Walker, L. Shao, and R. A. Volz. Estimating 3D location parameters using dual number quaternions. CVGIP: Image Understanding, 54(3):358–367, Nov. 1991.
Comments
There are no comments yet.