Tweet
Let's call R = [i j k] the 3x3 orientation matrix, where i, j, k are column
vectors. Notice that, according to the most widely used convention, i, j, k
are the positions of the tips of the versors of the local (technical or
anatomical) coordinate system L (versor = unit vector codirectional with a
directed axis), represented in the global coordinate system G.
DIFFERENTIATION
I would use smoothing quintic spline functions to compute the derivative of
each element of R. This is because spline functions can be differentiated
very easily (they are piecewise polynomials). With MATLAB, for instance, you
can use SPAPS to build the nine spline functions, and FNDER to obtain their
first and second derivatives. These derivatives are the linear velocities
and accelerations of the tips of the versors. In my opinion, quintic is
preferable to cubic spline (this is related to the fact that, in short,
quintic spline smooths by reducing jerk, while cubic spline does it by
reducing acceleration).
COMPUTING ANGULAR VELOCITY
1. Why should you first convert R into Cardan angles, if you can use
directly method 2?
2. [omega]=[M'(t)][M(t)]^-1 is called Poisson's equation, as far as I know.
And it does not yield an approximation of omega ([~omega]). The
approximation comes from the values of M' and M. Obviously,
[~omega]=[~M'(t)][~M(t)]^-1. This is the correct method, together with the
equivalent method which uses a simple linear combination of three cross
divisions:
2b. [omega]= 0.5 * (i'/i + j'/j + k'/k) = 0.5 * (i x i' + j x j' + k x k'),
where "/" stands for orthogonal cross division. But my article on
"Anticrossproducts and cross division" has not been published yet (in press,
coming soon on the Journal of Biomechanics; DOI:
10.1016/j.jbiomech.2007.09.030; when published on line, I invite you to read
sections see sections 5.3. and 5.4, and figures 3 and 8, which in my opinion
very nicely show the geometrical interpretation of the angular velocity
vector, and its relationship with the linear velocity of the tips of i, j,
k).
3. This PiG method can only give the cosines of three angles that are not
Euler or Cardan angles. How can you compute the angular velocity vector from
that? For instance, if by computing the arccosine of the thee dot products
you discover that i and j rotated by 0.05 rad in 0.05 s, and k rotated by 0
deg, then you can say that the angular velocity was 1 rad/s about the z
axis. But if all the versors rotated by a nonnull angle, then I guess this
method just doesn't work.
4. The cross products (and the ensuing normalization needed to obtain
versors i, j and k) amplify or attenuate the stereophotogrammetric error. In
the past, a lot of papers dealt with the problem of smoothing before or
after differentiation, and before or after center of mass and "joint center"
position estimation. Estimating R means estimating the position of three
points, at a distance of 1 mm or 1 m from the origin of L (depending on the
unit you are using for length). These points are not real markers, they are
"virtual" markers. Whatever you decide, you should consider that the noise
in the position of virtual markers is different from the
stereophotogrammetric error (error in the reconstructed position of real
markers).
With kind regards,
Paolo de Leva
-----Messaggio originale-----
Da: * Biomechanics and Movement Science listserver
[mailto:BIOMCH-L@NIC.SURFNET.NL] Per conto di Rettig, Oliver
Inviato: venerdì 9 novembre 2007 15.48
A: BIOMCH-L@NIC.SURFNET.NL
Oggetto: Estimation of angular velocity/acceleration from 3x3 matrices
Hallo,
I am interested in how to calculate angular velocity/acceleration from a
time serie of 3x3 rotation matrices which I get from marker positions
measured by a Vicon motion analysis system. Different ways are possible and
I want to learn which of the methods are used in practice and which
numerical advantages or disadvantages have the different methods.
The following possiblities are my startpoint:
1. First calculate time series of cardan angles with an arbitray rotation
order. Then you can estimate the derivations of these three time series of
double values by a) simply calculating differences or by b) low order
polynom derivations (a kind of filtering and differentiating in one step,
different orders and techniques are possible). One disadvantage of this
method can be gimbal lock.
2. Use the formula [~omega]=[M'(t)][M(t)]^-1 to calculate the Tensor which
includes the components of the angular velocity by time derivatives of each
component of the rotation Matrix M(t) multiplicated with the inverse of the
rotations matrix. For the estimation of the time derivative of the component
the same methods from 1. a),b) can be used. In Comparison with 1. gimbal
lock should be no problem.
3. You can simply estimate the velocity by calculating dot products of the
columns between two sequent matrices. The implementation of Vicon PiG seems
to do this by using matrices of frames with a time distance of 0.05s by
120Hz frame rate.
4. Because my rotations matrices are typically based on cross products of
vectors between markers I can analytical calculate formulars of the angular
velocity as function from the derivations of the marker positions.
Which of these methods do you use? Are there further methods? Which
filtering details you use? Are there further practical details I have to
look at? I am also interested in literature.
best regards
Oliver Rettig
---
Oliver Rettig
Stiftung Orthopädische Universitätsklinik Heidelberg
Ganglabor
Schlierbacher Landstr. 200a
69118 Heidelberg
Germany
Tel: +49 6221-96 6720
Fax: +49 6221-96 6725
vectors. Notice that, according to the most widely used convention, i, j, k
are the positions of the tips of the versors of the local (technical or
anatomical) coordinate system L (versor = unit vector codirectional with a
directed axis), represented in the global coordinate system G.
DIFFERENTIATION
I would use smoothing quintic spline functions to compute the derivative of
each element of R. This is because spline functions can be differentiated
very easily (they are piecewise polynomials). With MATLAB, for instance, you
can use SPAPS to build the nine spline functions, and FNDER to obtain their
first and second derivatives. These derivatives are the linear velocities
and accelerations of the tips of the versors. In my opinion, quintic is
preferable to cubic spline (this is related to the fact that, in short,
quintic spline smooths by reducing jerk, while cubic spline does it by
reducing acceleration).
COMPUTING ANGULAR VELOCITY
1. Why should you first convert R into Cardan angles, if you can use
directly method 2?
2. [omega]=[M'(t)][M(t)]^-1 is called Poisson's equation, as far as I know.
And it does not yield an approximation of omega ([~omega]). The
approximation comes from the values of M' and M. Obviously,
[~omega]=[~M'(t)][~M(t)]^-1. This is the correct method, together with the
equivalent method which uses a simple linear combination of three cross
divisions:
2b. [omega]= 0.5 * (i'/i + j'/j + k'/k) = 0.5 * (i x i' + j x j' + k x k'),
where "/" stands for orthogonal cross division. But my article on
"Anticrossproducts and cross division" has not been published yet (in press,
coming soon on the Journal of Biomechanics; DOI:
10.1016/j.jbiomech.2007.09.030; when published on line, I invite you to read
sections see sections 5.3. and 5.4, and figures 3 and 8, which in my opinion
very nicely show the geometrical interpretation of the angular velocity
vector, and its relationship with the linear velocity of the tips of i, j,
k).
3. This PiG method can only give the cosines of three angles that are not
Euler or Cardan angles. How can you compute the angular velocity vector from
that? For instance, if by computing the arccosine of the thee dot products
you discover that i and j rotated by 0.05 rad in 0.05 s, and k rotated by 0
deg, then you can say that the angular velocity was 1 rad/s about the z
axis. But if all the versors rotated by a nonnull angle, then I guess this
method just doesn't work.
4. The cross products (and the ensuing normalization needed to obtain
versors i, j and k) amplify or attenuate the stereophotogrammetric error. In
the past, a lot of papers dealt with the problem of smoothing before or
after differentiation, and before or after center of mass and "joint center"
position estimation. Estimating R means estimating the position of three
points, at a distance of 1 mm or 1 m from the origin of L (depending on the
unit you are using for length). These points are not real markers, they are
"virtual" markers. Whatever you decide, you should consider that the noise
in the position of virtual markers is different from the
stereophotogrammetric error (error in the reconstructed position of real
markers).
With kind regards,
Paolo de Leva
-----Messaggio originale-----
Da: * Biomechanics and Movement Science listserver
[mailto:BIOMCH-L@NIC.SURFNET.NL] Per conto di Rettig, Oliver
Inviato: venerdì 9 novembre 2007 15.48
A: BIOMCH-L@NIC.SURFNET.NL
Oggetto: Estimation of angular velocity/acceleration from 3x3 matrices
Hallo,
I am interested in how to calculate angular velocity/acceleration from a
time serie of 3x3 rotation matrices which I get from marker positions
measured by a Vicon motion analysis system. Different ways are possible and
I want to learn which of the methods are used in practice and which
numerical advantages or disadvantages have the different methods.
The following possiblities are my startpoint:
1. First calculate time series of cardan angles with an arbitray rotation
order. Then you can estimate the derivations of these three time series of
double values by a) simply calculating differences or by b) low order
polynom derivations (a kind of filtering and differentiating in one step,
different orders and techniques are possible). One disadvantage of this
method can be gimbal lock.
2. Use the formula [~omega]=[M'(t)][M(t)]^-1 to calculate the Tensor which
includes the components of the angular velocity by time derivatives of each
component of the rotation Matrix M(t) multiplicated with the inverse of the
rotations matrix. For the estimation of the time derivative of the component
the same methods from 1. a),b) can be used. In Comparison with 1. gimbal
lock should be no problem.
3. You can simply estimate the velocity by calculating dot products of the
columns between two sequent matrices. The implementation of Vicon PiG seems
to do this by using matrices of frames with a time distance of 0.05s by
120Hz frame rate.
4. Because my rotations matrices are typically based on cross products of
vectors between markers I can analytical calculate formulars of the angular
velocity as function from the derivations of the marker positions.
Which of these methods do you use? Are there further methods? Which
filtering details you use? Are there further practical details I have to
look at? I am also interested in literature.
best regards
Oliver Rettig
---
Oliver Rettig
Stiftung Orthopädische Universitätsklinik Heidelberg
Ganglabor
Schlierbacher Landstr. 200a
69118 Heidelberg
Germany
Tel: +49 6221-96 6720
Fax: +49 6221-96 6725