PDA

View Full Version : joint attitude math & s/w



Herman J. Woltring
02-21-1990, 10:40 PM
Dear Biomch-L readers,

Some readers on the list have shown interest in the mathematics and programming
details underlying the current joint attitude parametrization debate. For this
purpose, I'd like to refer them to a test package PRP FORTRAN which is available
to the readership by sending the request

SEND PRP FORTRAN BIOMCH-L

to LISTSERV@HEARN.EARN (or .BITNET); the request should be sent from the address
at which the Biomch-L listserver knows the requester. A FORTRAN-77 package is
then returned which accomodates all possible Cardanic conventions and the
proposed helical convention. While the PRP package was written for use under a
VAX/VMS system, only a few I/O statements in the package's testprogramme may
have to be modified for use in different environments.

Basically, the relation between a 3*3 orthonormal matrix Rijk (with |Rijk| = +1)
and the selected Cardanic angles can be described as follows:

Let Rijk = Ri(PHIi) * Rj(PHIj) * Rk(PHIk), where PHIi, PHIj, and PHIk denote
`planar' rotations, in the right-handed sense, about selected axes of a Carte-
sian co-ordinate system (1: x-axis, 2: y-axis, 3: z-axis). Defining ijk as
cyclic (i.e., 123, 231, or 312), the elements of Ri(PHIi) are defined as

rii = 1; rij = rji = rik = rki = 0
rjj = rkk = cos(PHIi); rkj = -rjk = sin(PHIi)

whence the elements r.. of Rijk can be expressed as:

rkj + rji = {1 + sin(PHIj)} sin(PHIi + PHIk)
rjj - rki = {1 + sin(PHIj)} cos(PHIi + PHIk)
rkj - rji = {1 - sin(PHIj)} sin(PHIi - PHIk)
rjj + rki = {1 - sin(PHIj)} cos(PHIi - PHIk)

This yields two unique pairs -- as per Codman's Paradox -- for PHIi and PHIk
unless |PHIj| = PI/2 -- gimbal-lock; in the latter case, two of the four rela-
tions above are zero. Note that all angles are modulo 2*PI. Given the choosen
PHIi and PHIk, PHIj can be derived from

rik = sin(PHIj)
{rii cos(PHIk) + rkk cos(PHIi)} - {rij sin(PHIk) + rjk sin(PHIi)} = cos(PHIj)

If ijk is anti-cyclic, i.e., 321, 213 or 132, the cyclic case can be obtained
via the intermediate transformation

Rijk(PHIi,PHIj,PHIk) = Rkji(-PHIk,-PHIj,-PHIi)'

Note that for other ijk, i.e., i=k, the n e u t r a l attitude corresponds to
gimbal-lock, with either PHIi + PHIk or PHIi - PHIk undefined. Interestingly,
Leonhard Euler's original publication `De Immutatione Coordinatarum, Caput IV,
Appendix de Superficiebus' in his `Introductio in Analysin Infinitorum' (Lausan-
ne, France 1748) provides such a case.

For the helical convention, the following formulae apply:

(Rijk - Rijk')/2 = sin(theta) A(N)
(Rijk + Rijk')/2 = cos(theta) I + {1-cos(theta)} N N'
THETA = theta N, N = (n1,n2,n3)', |N| = 1.

where A(N) is a skew-symmetric matrix with elements aji = -aij = nk and akk = 0,
for cyclic ijk. This system has a unique solution for 0 < theta < PI. If theta
= 0, no solution exists for N, but THETA is zero; for theta = PI, there are two
solutions, one THETA, the other -THETA, as one would expect for any periodic
solution modulo 2*PI.

Those who are interested in written material may get in touch with me directly,
and I'll be glad to send a mathematically oriented paper presented during a
1989 CAMARC Workshop on Clinical Protocols in Ancona/Italy.

Regards -- Herman J. Woltring, tel & fax +31.40.41 37 44