pracma (version 1.9.9)

procrustes: Solving the Procrustes Problem

Description

procrustes solves for two matrices A and B the `Procrustes Problem' of finding an orthogonal matrix Q such that A-B*Q has the minimal Frobenius norm.

kabsch determines a best rotation of a given vector set into a second vector set by minimizing the weighted sum of squared deviations. The order of vectors is assumed fixed.

Usage

procrustes(A, B)
kabsch(A, B, w = NULL)

Arguments

A, B

two numeric matrices of the same size.

w

weights , influence the distance of points

Value

procrustes returns a list with components P, which is B*Q, then Q, the orthogonal matrix, and d, the Frobenius norm of A-B*Q.kabsch returns a list with U the orthogonal matrix applied, R the translation vector, and d the least root mean square between U*A+R and B.

Details

The function procrustes(A,B) uses the svd decomposition to find an orthogonal matrix Q such that A-B*Q has a minimal Frobenius norm, where this norm for a matrix C is defined as sqrt(Trace(t(C)*C)), or norm(C,'F') in R.

Solving it with B=I means finding a nearest orthogonal matrix.

kabsch solves a similar problem and uses the Procrustes procedure for its purpose. Given two sets of points, represented as columns of the matrices A and B, it determines an orthogonal matrix U and a translation vector R such that U*A+R-B is minimal.

References

Golub, G. H., and Ch. F. van Loan (1996). Matrix Computations. 3rd Edition, The John Hopkins University Press, Baltimore London. [Sect. 12.4, p. 601]

Kabsch, W. (1976). A solution for the best rotation to relate two sets of vectors. Acta Cryst A, Vol. 32, p. 9223.

See Also

svd

Examples

##  Procrustes
U <- randortho(5)               # random orthogonal matrix
P <- procrustes(U, eye(5))

##  Kabsch
P <- matrix(c(0, 1, 0, 0, 1, 1, 0, 1,
              0, 0, 1, 0, 1, 0, 1, 1,
              0, 0, 0, 1, 0, 1, 1, 1), nrow = 3, ncol = 8, byrow = TRUE)
R <- c(1, 1, 1)
phi <- pi/4
U <- matrix(c(1, 0, 0,
              0, cos(phi), -sin(phi),
              0, sin(phi),  cos(phi)), nrow = 3, ncol = 3, byrow = TRUE)

Q <- U %*% P + R
K <- kabsch(P, Q)
# K$R == R  and  K$U %*% P + c(K$R) == Q