Means and Averaging in the Group of Rotations.
In this paper we give precise definitions of different, properly invariant notions of mean or average rotation. Each mean is associated with a metric in SO(3). The metric induced from the Frobenius inner product gives rise to a mean rotation that is given by the closest special orthogonal matrix to the usual arithmetic mean of the given rotation matrices. The mean rotation associated with the intrinsic metric on SO(3) is the Riemannian center of mass of the given rotation matrices. We show that the Riemannian mean rotation shares many common features with the geometric mean of positive numbers and the geometric mean of positive Hermitian operators. We give some examples with closed-form solutions of both notions of mean.
where is the angle of rotation of R. The kth root exp(1 Log R) is the one for which
the eigenvalues have the largest positive real part, and is the only one we denote by
R1=k. In the case it is the only square root with positive real part.
2.2. Metrics in SO(3). A straightforward way to dene a distance function in
SO(3) is to use the Euclidean distance of the ambient space M(3), i.e., if R1 and
are two rotation matrices then
where k kF is the Frobenius norm which is induced by the Euclidean inner product,
known as the Frobenius inner product, dened by hR1; It is
easy to see that this distance is bi-invariant in SO(3), i.e., dF (P R1Q; P
dF (R1; R2) for all P ; Q in SO(3).
Another way to dene a distance function in SO(3) is to use its Riemannian
structure. The Riemannian distance between two rotations
is the length of the shortest geodesic curve that connects R1 and given by
4 M. MOAKHER
Note that the geodesic curve of minimal length may not be unique. If RT1 R2 is an
involution, in other words if (RT1 rotation through an angle , then
and can be connected by two curves of equal length. In such a case, the
rotations and are said to be antipodal points in SO(3) and is said to be
the cut point of R1 and vice versa.
The Riemannian distance (2.6) is also bi-invariant in SO(3). Indeed, using the
fact [3], we can show that dR(P R1Q; P
Remark 2.1. The Euclidean distance (2.5) represents the chordal distance between
and R2, i.e., the length of the Euclidean line segment in the space of M(3)
(except for the end points R1 and R2, this line segment does not lie in SO(3)), whereas
the Riemannian distance (2.6) represents the arc-length of the shortest geodesic curve
(great-circle arc), which lies entirely in SO(3), passing through R1 and R2.
Remark 2.2. If denotes the angle of rotation of RT1
Therefore, when theprotations R1 and are
suciently close, i.e., is small, we have dF (R1; R2) 2 dR(R1; R2).
2.3. Covariant derivative and Hessian. We recall that the tangent space at
a point R of SO(3) is the space of all matrices such that RT is skew symmetric
and that the normal space (associated with the Frobenius inner product) at R consists
of all matrices N such that RT N is symmetric [5].
For a real-valued function f(R) dened on SO(3), the covariant derivative rf is
the unique tangent vector at R such that

d
dt
where Q(t) is a geodesic emanating from R in the direction of , i.e.,
R exp(tA) and
The Hessian of f(R) is given by the quadratic form
d2
dt2
where Q(t) is a geodesic and is in the tangent space at R as above.
2.4. Geodesic convexity. We recall that a subset A of a Riemannian manifold
M is said to be convex if the shortest geodesic curve between any two points x and y
in A is unique in M and lies in A. A real-valued function dened on a convex subset
A of M is said to be convex if its restriction to any geodesic path is convex, i.e., if
its domain for all x 2 M and u 2 Tx(M),where exp is the exponential map at x.
x
With these denitions, one can readily see that any geodesic ball Br(Q) in SO(3)
of radius r less than around Q is convex and that the real-valued function f denedon Br(Q) by when r is less than 2 . Geodesic balls
with radius greater or equal than are not convex.3. Mean rotation. For a given set of N points xn; in IRd the
mean x is given by the barycenter of the N points. The
mean also has a variational property; it minimizes the sum of the squared
distances to the given points xn;
x2IRd n=1
where here de(; ) represents the usual Euclidean distance in IRd.
One can also use the arithmetic mean to average N positive real numbers xn >
and the mean is itself a positive number. In many applications,
however, it is more appropriate to use the geometric mean to average positive numbers,
which is possible because positive numbers form a multiplicative group. The geometric
1=N 1=N
mean also has a variational property; it minimizes the sum of the
squared hyperbolic distances to the given data
where dh(x; log yj is the hyperbolic distance1 between x and y.
As we have seen, for the set of positive real numbers dierent notions of mean can
be associated with dierent metrics. In what follows, we will extend these notions of
mean to the group of proper orthogonal matrices.
By analogy with IRd, a plausible denition of the mean of N rotation matrices
is that it is the minimizer in SO(3) of the sum of the squared
distances from that rotation matrix to the given rotation matrices
represents a distance in
SO(3). Now the two distance functions (2.5) and (2.6) dene the two dierent means.
Definition 3.1. The mean rotation in the Euclidean sense, i.e., associated with
the metric (2.5), of N given rotation matrices is dened as
Definition 3.2. The mean rotation in the Riemannian sense, i.e., associated
with the metric (2.6), of N given rotation matrices is dened as
The minimum here is understood to be the global minimum. We remark that
in IRd, or in the set of positive numbers, the objective functions to be minimized
are convex over their domains, and therefore the means are well dened and unique.
However, in SO(3), as we shall see, the objective functions in (3.3) and (3.4) are not
(geodesically) convex, and therefore the means may not be unique.
Before we proceed to study these two means, we note that both satisfy the following
desirable properties that one would expect from a mean in SO(3), and that
are counterparts of properties of means of numbers, namely,
1We borrow this terminology from the hyperbolic geometry of the Poincare upper half-plane. In
fact, the hyperbolic length of the geodesic segment joining the points P (a; y1) and Q(a; y2), y1; y2 > 0
is j log y1 j, (see [26]).
6 M. MOAKHER
1. Invariance under permutation: For any permutation of the numbers 1
through N, we have
2. Bi-invariance: If R is the mean rotation of fRng;
is the mean rotation of fP RnQg; every P and Q in SO(3). This
property follows immediately from the bi-invariance of the two metrics dened above.
3. Invariance under transposition: If R is the mean rotation of fRng;
then RT is the mean rotation of fRT
We remark that the bi-invariance property is in some sense the counterpart of the
homogeneity property of means of positive numbers (but here left and right multiplication
are both needed because the rotation group is not commutative).
3.1. Characterization of the Euclidean mean. The following proposition
gives a relation between the Euclidean mean and the usual arithmetic mean.
Proposition 3.3. The mean rotation
is the orthogonal projection of R = onto the special orthogonal group SO(3).
In other words, the mean rotation in the Euclidean sense is the projection of the
arithmetic mean R of in the linear space M(3) onto SO(3).
Proof. As are all orthogonal, it follows that
On the other hand, the orthogonal projection of R onto SO(3) is given by
"XN XN Rn RTm XN RTn #
Because of Proposition 3.3, the mean in the Euclidean sense will be termed the
projected arithmetic mean to reect the fact that it is the orthogonal projection of the
usual arithmetic mean in M(3) onto SO(3).
Remark 3.4. The projected arithmetic mean can now be seen to be related to the
classical orthogonal Procrustes problem [10], which seeks the orthogonal matrix that
most closely transforms a given matrix into a second one.
Proposition 3.5. If det R is positive, then the mean rotation in the Euclidean
sense given by the unique polar factor in
the polar decomposition [10] of R.
Proof. Critical points of the objective function
dened on SO(3) and corresponding to the minimization problem (3.3) are those
elements of SO(3) for which the covariant derivative of (3.5) vanishes. Using (2.8)

we get Therefore, critical points of (3.5) are the

rotation matrices R such that PNn=1 R RTn R RT equivalently, for
which the matrix S dened by
is symmetric.
Since R is orthogonal, and both S and are symmetric, it follows
that N2M. Therefore, there exists an orthogonal matrix U such that
of M. The eight possible square roots of M are UT diag(
To determine the square root of N2M that corresponds
to the minimum of (3.5) we require that the Hessian of the objective function (3.5)
at R given by (3.6) be positive for all tangent vectors at R. From (2.9) we obtain

Hess therefore at R given by (3.6) we have
Hess
where a; b; c are such that = UT RBU and
b a 0
As we are looking for a proper rotation matrix, i.e., an orthogonal matrix with
determinant one, it follows from (3.6) that det
that Hess F(; ) is positive for all tangent vectors at R if and only
positive
and In fact, (3.5) has four cpriticapl poinpts belonging to
which consist of apminimpum [(1p; 2;
Hence, the projected arithmetic mean is given by
which, when det R > 0, coincides with the polar factor of the polar decomposition of
R. Of course uniqueness fails when the smallest eigenvalue of M is not simple.
Remark 3.6. The case where det is a degenerate case. However, if R
has rank 2, i.e., when nd a unique closest proper
orthogonal matrix to R (see [6] for details), and hence can dene the mean rotation
in the Euclidean sense.
3.2. Characterization of the Riemannian mean. First, we compute the
derivative of the real-valued function H(P
to t where P is the geodesic emanating from R in the direction of
RA. As is in the tangent space at R, we have
be the angle of rotation of QT P (t), i.e., such that
8 M. MOAKHER
Dierentiate (3.8) to get d H(P is the
dt t=0 sin
angle of rotation of QT R and we have used the fact that H(P
Recall that, since A is skew symmetric, symmetric matrix S.
It follows that tr(QT
Then, with the help of (2.4) we obtain d H(P
dt t=0
fore, the covariant derivative of H is given by
The second derivative of (3.8) gives
d2 sin cos
Let U be an orthogonal matrix and B the skew-symmetric matrix such that
Then, as tr(QT it is easy to see that
d2 sin
The RHS of (3.10) is always positive for arbitrary a, b, c in IR and 2 (; ). It
follows that Hess H(; ) is positive for all tangent vectors .
denote the objective function of the minimization problem (3.4), i.e.,
Using the above, the covariant derivative of G is found to be
Therefore, a necessary condition for regular extrema of (3.11) is
By (3.10) we conclude that the Hessian Hess G(; ) of the objective function (3.11)
is positive for all tangent vectors . Therefore, equation (3.12) characterizes local
minima of (3.11) only. As a matter of fact, local maxima are not regular points, i.e.,
they are points where (3.11) is not dierentiable.
It is worth noting that, as R = R , the characterization for the Riemannian
mean given in (3.12) is similar to the characterization
of the geometric mean (3.2) of positive numbers. However, while in the scalar case
the characterization (3.13) has the geometric mean as unique solution, the characterization
(3.12) has multiple solutions, and hence is a necessary but not a sucient
condition to determine the Riemannian mean. The lack of uniqueness of solutions of
(3.12) is akin to the fact that, due to the existence of a cut point for each element of
SO(3), the objective function (3.11) is not convex over its domain.
In general, closed-form solutions to (3.12) cannot be found. However, for some
special cases solutions can be given explicitly. In the following subsections, we will
present some of these special cases.
Remark 3.7. The Riemannian mean of may also be called the Riemannian
barycenter of which is a notion introduced by Grove, Karcher
and Ruh [11]. In [17] it was proven that for manifolds with negative sectional curva-
ture, the Riemannian barycenter is unique.
3.2.1. Riemannian mean of two rotations. Intuitively, in the case
the mean rotation in the Riemannian sense should lie midway between R1 and
along the shortest geodesic curve connecting them, i.e., it should be the rotation
R2)1=2. Indeed, straightforward computation shows that R1(RT1 R2)1=2 does
satisfy condition (3.12). Alternatively, equation (3.12) can be solved analytically as
follows. First, we rewrite it as
then we take the exponential of both sides to obtain multiplying
both sides with RT1 R we get (RT1 Such an equation has two
solutions in SO(3) that correspond to local minima of (3.11). However, the global
minimum is the one that corresponds to taking the square root of the above equation
that has eigenvalues with positive real part, i.e., (RT1 R2)1=2. Therefore, for two non-
antipodal rotation matrices R1 and R2, the mean in the Riemannian sense is given
explicitly by
The second equality can be easily veried by pre-multiplying R1(RT1 R2)1=2 by R2RT2
which is equal to I. This makes it clear that G is symmetric with respect to R1 and
R2, i.e., G(R1;
3.2.2. Riemannian mean of rotations in a one-parameter subgroup. In
the case where all matrices Rn; belong to a one-parameter subgroup
of SO(3), i.e., they represent rotations about a common axis, we expect that their
mean is also in the same subgroup. Further, one can easily show that equation (3.12)
Y
reduces to saying that R is an Nth root of Rn. Therefore, the Riemannian mean
is the Nth root that yields the minimum value of the objective function (3.11).
In this case, all rotations lie on a single geodesic curve. One can show that
the geometric mean G(R1; R2; R3) of three rotations R1, and R3 such that
3, is the rotation that is located at 32 of the length
of the shortest geodesic segment connecting R1 and G(R2; R3), i.e., the rotation
have
M. MOAKHER
This explicit formula does not hold in the general case due to the inherent curvature
of SO(3), see the discussion at the end of Example 2 below.
When the rotations belong to a geodesic segment of length less than
and centered at the identity, the above formula reduces to
1=N 1=N
Once again we see the close similarity between the geometric mean of positive numbers
and the Riemannian mean of rotations. This is to be expected since both the set of
positive numbers and SO(3) are multiplicative groups, and we have used their intrinsic
metrics to dene the mean. For this reason, we will call the mean in the Riemannian
sense the geometric mean.
3.3. Equivalence of both notions of mean of two rotations. In the follow-
ing, we show that for two rotations the projected arithmetic mean and the geometric
mean coincide. First, we prove the following lemma.
Lemma 3.8. Let R1 and R2 be two elements of SO(3), then det(R1
Proof. Consider the real-valued function dened on [0; 1] by
We see that this function is continuous with
Assume that f(1) < 0, i.e., det(R1 there exists in [0; 1] such that
Hence, must be in the spectrum of RT2 R1 which is a proper orthogonal matrix. But
this cannot happen, which contradicts the assumption that det(R1
In general, the result of the above lemma does not hold for more than two rotations
matrices. We will see examples of three rotation matrices for which the determinant
of their sum can be negative.
Proposition 3.9. The polar factor of the polar decomposition of
and are two rotation matrices, is given by R1(RT1 R2)1=2.
Proof. Let Q be the proper orthogonal matrix and S be the positive-denite
matrix such that QS is the unique polar decomposition of
R1. One can easily verify that (RT1 R2)1=2
(RT1 R2)1=2 is the positive-denite square root of 2I and that the
inverse of this square root is given by Hence, the
polar factor is
Since the polar decomposition is unique, the result of this proposition together
with the previous lemma shows that both notions of mean agree for the case of
two rotation matrices. For more than two rotations, however, both notions of mean
coincide only in special cases that present certain symmetries. In Example 2 of x 4
below, we shall consider a two-parameter family of cases illustrating this coincidence.
4. Analytically solvable examples. In this section we present two cases in
which we can solve for both the projected arithmetic mean and the geometric mean
explicitly. These examples help us gain a deeper and concrete insight to both notions
of mean. Furthermore, Example 2 conrms our intuitive idea that for \symmetric"
cases, both notions of mean agree.
4.1. Example 1. We begin with a simple example where all rotation matrices
for which we want to nd the mean lie in a one-parameter subgroup of SO(3). Using
the bi-invariance property we can reduce the problem to that of nding the mean of
MEAN ROTATION 11
Projected arithmetic mean: The arithmetic sum of these matrices has a positive
determinant n)2. Hence, the projected arithmetic mean
of the given matrices is given by the polar factor of the polar decomposition of their
sum. After performing such a decomposition we nd thatN
2cos a sin a >< cos
> sin a = sin n:
Such a mean is well dened as long as This mean agrees with the notionof directional mean used in the statistics literature for circular and spherical data
[20, 7, 9, 8]. The quantity 1 r=N, which is called the circular variance, is a measure
of dispersion of the circular data . The direction dened by the angle a
is called the mean direction of the directions dened by
Geometric mean: Solutions of (3.12) are given by
cos l sin l 0 N
4sin l cos l 05; where
The geometric mean of these rotation matrices is therefore the solution that yields
the minimum value of the objective function (3.11). Of course, as we have seen in x 3,
the geometric mean is given explicitly by (3.15).
Note that, even though elements of a one-parameter subgroup commute, the two
rotations (3.15) and (3.16) are dierent. This is due to choice of the kth root of a
rotation matrix to be the one with eigenvalues that have the largest positive real parts.
To see this, consider the case
where P is a rotation of an angle about the z-axis while R
If the rotation matrices Rn are such that n <
certain number 2 IR, then their geometric mean is a rotation about the z-axis of
an angle
The geometric mean rotation of the rotations given by (4.1) coincides with the
concept of median direction of circular data [20, 7].
Remark 4.1. When neither the
projected arithmetic mean nor the geometric mean is well dened. On the one hand
so the projected arithmetic mean is not dened, while on the other hand
the objective function (3.11) for the geometric mean has two local minima with the
same value, namely, R1(RT1 R2)1=2 and its cut value andtherefore the global minimum is not unique.
Let F~ and G~ be the functions dened on [; ] such that
any rotation R about the z-axis through an angle , i.e., F~ and
G~ are the restrictions of the objective functions (3.5) and (3.11) to the subgroup
considered in this example. In Fig. 4.1 we give the plots of F~ and G~ for the sets of
data takes several dierent
values. It is clear that neither (3.5) nor (3.11) is convex. While the function (3.5)
M. MOAKHER
Projected Arithmetic Mean Geometric Mean10
~
a=p/4
~
a=p/4
a=p
Fig. 4.1. Plots of the objective functions F~() and G~() for dierent values of . Note that
is constant and G~ has four local minima with an equal value. Consequently, neither
the projected arithmetic mean nor the geometric mean is well dened.
is smooth the function (3.11) has cusp points but only at local maxima. However,
if the given rotations are located in a geodesic ball of radius less than =2, i.e., in
this example have angles i such that ji jj < then the objective
functions restricted to this geodesic ball are convex and hence the means are well
dened. Such case is illustrated in Fig. 4.2 which shows plots of F~ and G~ for the
following sets of data takes several
dierent values.
Projected Arithmetic Mean Geometric Mean8
~
a=p/4
~
a=p/4
Fig. 4.2. Plots of the objective functions F~() and G~() for dierent values of . Restricted to
[=4; 3=4], i.e., between the dashed lines, the objective functions are indeed convex.
4.2. Example 2. In the second example we consider N elements of SO(3) that
represent rotations through an angle about the axes dened by the unit vectors
sin sin n; cos ]T , where
Projected arithmetic mean: Straightforward computations show that the projected
arithmetic mean is given by
cos a sin a 0 ><cos
a cos a 05;
By using half-angle tangent formulas in the above we obtain the following simple
relation between a and
a
Geometric mean: Since the rotation axes are symmetric about the z-axis, and
the rotations share the same angle, we expect that their geometric mean is a rotation
about the z-axis through a certain angle g. Furthermore, because of this symmetry
we also expect that the mean in the Euclidean sense agrees with the one in the
Riemannian sense.
From the Campbell-Baker-Hausdor formula for elements of SO(3) [23] we have
where the coecients a; b; c and are given by
a
Therefore, the characterization (3.12) of the geometric mean reduces to
a Log Rn bN Log R
This is a matrix equation in so(3), which is equivalent to a system of three nonlinear
equations. Because the axes of rotation of Rn are symmetric about the z-axis
we have cos It follows that [Log Rn; Log
R. Therefore, this system reduces to the following
single equation for the angle g
which when compared with (4.2) indeed shows that a = g and therefore the projected
arithmetic mean and the geometric mean coincide.
This example provides a family of mean problems parameterized by and
where the projected arithmetic and geometric mean coincide. We now further examine
the problem of nding the mean of three rotations about the three coordinate axes
14 M. MOAKHER
through the same angle , which, by the bi-invariance property of both means, can be
considered as a special case of this two-parameter family with .Therefore the mean of these three rotations is a rotatiopn through an angle about
the axis generated by the vector [1; 1; 1]T with tan . The rotations R1,
and R3 form a geodesic equilateral triangle in SO(3). By symmetry arguments
the geometric mean should be the intersection of the three geodesic medians, i.e., the
geodesic segments joining the vertices of the geodesic triangle to the midpoints of the
opposite sides. In at geometry, this intersection is located at two-thirds from the
vertices of the triangle. However, in the case of SO(3), due to its intrinsic curvature,
this is not true. The ratio of the length of the geodesic segment joining one rotation
and the geometric mean, to the length of the geodesic median joining this rotation
and the midpoint of the geodesic curve joining the two other rotations is plotted as a
function of the angles in Fig. 4.3.0.65g0.620 p/4 p/2 3p/4 p
Fig. 4.3. Plot of the ratio of the geodesic distance from one vertex to the barycenter over
the geodesic distance from this vertex to the midpoint of the opposed edge in the geodesic equilateral
triangle in SO(3). The departure of from 2/3, which is due to the curvature of SO(3), increases
with the length, , of the sides of the triangle.
5. Weighted means and power means. Our motivation of this work was to
construct a lter that smooths the rotation data giving the relative orientations of
successive base pairs in a DNA fragment, see [19] for details. Such a lter can be a
generalization of moving window lters, which are based on weighted averages, used
in linear spaces to smooth noisy data. The construction of such lters and the direct
analogy we have found between the arithmetic and geometric means in the group of
positive numbers, and the projected arithmetic and geometric means in the group
of rotations, have led us to the introduction of weighted means and power means of
rotations that we discuss next.
Definition 5.1. The weighted projected arithmetic mean of N given rotations
dened as
This mean satises the bi-invariance property. Using similar arguments as for the
projected arithmetic mean one can show that the weighted projected arithmetic mean
is given by the polar factor of the polar decomposition of the matrix
MEAN ROTATION 15
provided that det A is positive.
Definition 5.2. The weighted geometric mean of N rotations
weights dened as
This mean also satises the bi-invariance property. Using arguments similar to
those used for the geometric mean, we can show that the weighted geometric mean is
characterized by n=1 wn Log(Rn
Definition 5.3. For a real number s such that 0 < jsj 1, we dene the
weighted s-th power mean rotation of N rotations
h is
We note that Mw
1 this is the weighted projected arithmetic mean. Because elements of SO(3) are
orthogonal, and the trace operation is invariant under transposition, the weighted
s-th power mean is the same as the weighted (s)-th power mean. Therefore, it is
immediate that the weighted projected harmonic mean, dened by
coincides with the weighted projected arithmetic mean.
This is a natural generalization of the s-th power mean of positive numbers and
it is in line with the fact that for positive numbers the s-th power mean
is given by the s-th root of the arithmetic mean of One has to
note, however, that for s such that 0 < jsj < 1 this mean is not invariant under the
action of elements of SO(3). This is not a surprise as the power mean of positive
numbers also does not satisfy the homogeneity property.
For the set of positive numbers [12] and similarly for the set of Hermitian denite
positive operators [25], there is a natural ordering of elements and the classical
mean inequalities holds. Furthermore, it is well known
[12, 25] that the s-th power mean converges to the geometric mean as s goes to 0.
However, for the group of rotations such a natural ordering does not exists. Nonethe-
less, one can show that if all rotations belong to a geodesic ball of radius
less than centered at the identity, then the projected power mean indeed convergesto the geometric mean as s tends to 0.
Analysis of numerical algorithms for computing the geometric mean rotation and
the use of the dierent notions of mean rotation for smoothing three-dimensional
orientation data will be published elsewhere.

Acknowledgment

. The author is grateful to Professor J. H. Maddocks for
suggesting this problem and for his valuable comments on this paper. He also thanks
the anonymous referee for his helpful comments.


--R

Dierential Geometry: Manifolds
Matrix Computations
Jacobi elds and Finsler metrics on compact Lie groups with an application to dierentiable pinching problem

optimization and Dynamical Systems

Maximum likelihood estimation for the matrix von Mises-Fisher and Bingham distributions
Fitting smooth paths to spherical data
Riemannian center of mass and mollier smoothing
The von Mises-Fisher matrix distribution in orientation statis- tics
A continuum rod model of sequence-dependent DNA structure
Statistics of directional data
A Mathematical Introduction to Robotic Manipula- tion
Fitting smooth paths to rotation data
Geometrical Methods in Robotics
optimization techniques on Riemannian manifolds
Hermitian semidenite matrix means and related matrix inequalities-an intro- duction
Convex functions and optimization methods on Riemannian manifolds
Equatorial distributions on a sphere
--TR

--CTR

Doug L. James , Christopher D. Twigg, Skinning mesh animations, ACM Transactions on Graphics (TOG), v.24 n.3, July 2005
Kavan , Steven Collins , Ji ra , Carol O'Sullivan, Skinning with dual quaternions, Proceedings of the 2007 symposium on Interactive 3D graphics and games, April 30-May 02, 2007, Seattle, Washington
Xavier Pennec, Intrinsic Statistics on Riemannian Manifolds: Basic Tools for Geometric Measurements, Journal of Mathematical Imaging and Vision, v.25 n.1, p.127-154, July      2006
Christophe Lenglet , Mikal Rousson , Rachid Deriche , Olivier Faugeras, Statistics on the Manifold of Multivariate Normal Distributions: Theory and Application to Diffusion Tensor MRI Processing, Journal of Mathematical Imaging and Vision, v.25 n.3, p.423-444, October   2006
Peter J. Basser , Sinisa Pajevic, Spectral decomposition of a 4th-order covariance tensor: Applications to diffusion tensor MRI, Signal Processing, v.87 n.2, p.220-236, February, 2007
special orthogonal group;rotation;geodesics;averaging;operator means