Graphics Reference
In-Depth Information
S
3
that's in the
q
1
-
q
2
plane and is perpendicular to
q
1
.
We subtract the projection
(
q
2
·
1. Find a vector
v
∈
q
1
)
q
1
of
q
2
onto
q
1
from
q
2
to get a vector
perpendicular to
q
1
, and normalize:
q
2
−
(
q
2
·
q
1
)
q
1
v
=
.
(11.41)
q
2
−
(
q
2
·
q
1
)
q
1
2. Find a unit-speed path along the great circle from
q
1
through
v
. That's
just
γ
(
t
)=cos(
t
)
q
1
+sin(
t
)
v
. This path reaches
q
2
at
t
=
θ
, where
θ
=
cos
−
1
(
q
1
·
q
2
)
is the angle between the two vectors.
3. Adjust the path so that it reaches
q
2
at time 1 rather than at time
θ
,by
multiplying
t
by
θ
.
The resultant code is shown in Listing 11.4.
Listing 11.4: Code for spherical linear interpolation between two quaternions.
1
2
3
4
5
6
7
8
9
10
double[4] slerp(double[4] q1, double[4] q2, double t)
{
assert(dot(q1, q1) == 1);
assert(dot(q2, q2) == 1);
// build a vector in q1-q2 plane that's perp. to q1
double[4] u = q2 - dot(q1, q2)
*
q2;
u = u / length(u); // ...and make it a unit vector.
double angle = acos(dot(q1, q2));
return cos(t
*
angle)
*
q1 + sin(t
*
angle)
*
u;
}
As the argument to the cosine ranges from
0
to
angle
, the result varies from
q1
to
q2
.
11.2.6.2 Interpolating between Rotations
We now have all the tools we need to interpolate between rotations. Suppose
that
M
1
and
M
2
are rotation matrices, with
M
1
corresponding to the quaternions
±
q
2
, as indicated schematically in
Figure 11.8. Starting with
q
1
, we determine which of
q
2
and
q
1
and
M
2
corresponding to the quaternions
±
q
2
is closer, and
then interpolate along an arc between
q
1
and this point; to find interpolating
rotations, we project to
SO
(
3
)
via
K
.
Listing 11.5 shows the pseudocode.
−
Listing 11.5: Code to interpolate between two rotations expressed as matrices.
1
2
3
4
5
6
7
8
9
10
11
12
13
14
Mat33 RotInterp(Mat33 m1, Mat33 m2, double t)
// find a rotation that's t of the way from m1 to m2 in SO(3).
// m1 and m2 must be rotation matrices.
{
if (
m
1
m
2
=
−
I
){
Report error; can't interpolate between opposite rotations.
}
Quaternion q1, q1p, q2, q2p;
RotationToQuaternion(m1, q1, q1p);
RotationToQuaternion(m2, q2, q2p);
if (Dot(q1, q2) < 0) q2 = q2p;
Quaternion qi = Quaternion.slerp(q1, q2, t);
return K(qi);
// K is the projection from S3 to SO(3)
}