NASA Technical Memorandum 101089

An Algorithm for the

Systematic Disturbance of

Optimal Rotational Solutions

Arthur J. Grunwald and Mary K. Kaiser

(NA_A-TM-!OIOSV) AN ALGORITHM FOR TflP

Sy3TEMATIC DIqTUR,_ANCF nF OPTIPIAL ROTATIONAL

SF)LUT!QNS (NASA) 29 p C_'_CL i2A

NQ0-i057_

August 1989

NASA

National Aeronautics and

Space Administration
" 'IlF- !.

llE_-

i
NASA Technical Memorandum 101089

An Algorithm for the

Systematic Disturbance of

Optimal Rotational Solutions

Arthur J. Grunwald, Technion, Haifa, Israel

Mary K. Kaiser, Ames Research Center, Moffett Field, California

August 1989

\

N/L_A

National Aeronautics and

Space Administration

Ames Research Center

Moffett Field, Califomia 94035

NOTATION

b

bo

bf

bs

Du-->v

e -= {ex,ey,ez }T

Mu----)v

M U'-_V

qu-->v =- {ql,q2,q3,q4} T

Angular orientation at time t, not necessarily along the most direct rotational path; also used to indicate the body system

Initial angular orientation

Final angular orientation

Angular orientation, for following the most direct rotational path

3x3 element direction cosine matrix (DCM), for transforming a coordinate

of the u-system into the v-system.

Exponent of disturbance function

Unity vector, indicating the direction of the Euler axis for rotation from an

initial system to a final one, expressed in the initial system

4x4 element quatemion matrix of the quaternion for rotation from u to v

Transmutated quaternion matrix

Quaternion for rotation from u to v, uniquely describing the orientation of

v into u

T Transpose of a matrix; also used to indicate a column vector

Time

tO

tf

x b - {xb, yb, zb}T

x i = {xi,yi,zi} T

Initial time

Final time

Orthogonal right-hand object-body coordinate system; in aeronautical

applications the xb-axis is pointing forward along the aircraft body main

axis, the yb-axis is pointing outward through the right wing, and the

zb-axis is pointing downward

Orthogonal right-hand Earth-based inertial reference system, with the

xi-axis pointing to the north, the yi-axis to the east, and the zi-axis downward, toward the center of Earth; in the experiments the xi-axis coincides

with the direction of the viewing axis

iii PRECEDING PAGE BLANK NOT FILMED
_xqo

xql ;xq2;xq3

O_

ctO

7,8

Aq

V

Like the initial system bo, fixed with respect to the inertial system, with

the xqO-axis aligned with the Euler axis and the yqO-axis in the xbO-o-ybO

plane

Intermediate systems to obtain the plane in which to rotate o_b in order to

be disturbed

"Angular distance" between orientation b and bf, i.e., the angle of rotation

for the remaining most direct path rotation from b to bf

Total angle of rotation for the most direct path rotation from initial orientation bo to final orientation bf

"Rotational deviation" from the most direct path; it is the amount of rotation of the quatenlion Aq

Averaged rotational deviation from the most direct path over the interval

tO _<t_<tf

Azimuth angle and elevation angle, respectively, which specify the orientation of _3 in the bo-system

Quaternion for the rotational deviation, i.e., for rotation from bs to b

Angle of disturbance; it is the angle over which cob is rotated in order to

be disturbed

v0

_d

coO_lC_l

b

0_d

Angle of disturbance at time tO

Azimuth angle and elevation angle, respectively, which specify the

orientation of c0q° in the q0-system

Euler angles for yaw, pitch, and roll, respectively, specifying the orientation of body with respect to the inertial system

Angle of plane in which _b is rotated to be disturbed

Vector of angular rotation, for the most direct rotational path between orientation b0 and bf, expressed in the b0-system

Constant rotational speed

Disturbed vector of angular rotation; it is the vector of rotation at which

the object proceeds rotating at time t, expressed in the b-system

iv
(.0$

b cob o)i___>b=

_r_ S

Vector of angular rotation, for the most direct rotational path between orientation b0 and bf, rotating at reduced angular rate, expressed in the

b0-system

Vector of angular rotation, for the remaining most direct rotational path

between orientation b and bf, expressed in the b-system

The vector _b expressed in the q_0-system

4×4 element matrix, containing the components of COdb, used for computing the rotational path

4×4 element matrix, containing the components of D.s, used for

computing the most direct rotational path

V

AN ALGORITHM FOR THE SYSTEMATIC DISTURBANCE

OF OPTIMAL ROTATIONAL SOLUTIONS

Arthur J. Grunwald I and Mary K. Kaiser

Ames Research Center

SUMMARY

An algorithm for introducing a systematic rotational disturbance into an optimal (i.e., single axis)

rotational trajectory is described. This disturbance introduces a motion vector orthogonal to the

quatemion-defined optimal rotation axis. By altering the magnitude of this vector, the degree of nonoptimality can be controlled. The metric properties of the distortion parameter are described, with

analogies to two-dimensional translational motion.

This algorithm has been implemented in a motion-control program on a three-dimensional graphic

workstation. It supports a series of human performance studies on the detectability of rotational trajectory optimality by naive observers.

1. INTRODUCTION

This paper describes an algorithm for generating kinematically suboptimal ("warped") rotational

trajectories. First, the basic idea of creating a suboptimal trajectory is demonstrated with a twodimensional translatory analogy. Next, the mathematical formulation for describing the rotation of an

object is given for three methods: (1) Euler angles, (2) direction cosine matrix, and (3) quaternions, and

the relation between these methods is discussed. Then, the rotational equivalent of the straight path for

translatory motion is defined and ways for computing the rotational path and the deviation from the

direct path are delineated. Finally, the method used for creating a reproducible "smoothly warped" suboptirnal trajectory is outlined. This method allows both the shape of the trajectory and the magnitude of

the deviation from the optimal path to be defined by two independent parameters. A flow chart is presented which summarizes the computations performed to create displays employing this algorithm.

lTcchnion, Haifa, Israel.
2. DEVIATION FROM THE OPTIMAL PATH TRAJECTORY FOR

TRANSLATORY MOTION

2.1 Method for Creating a Suboptimal Trajectory

The basic idea of creating a suboptimal trajectory can best be illustrated with a simple, twodimensional translational analogy. Consider two points on a plane as illustrated in figure 1. An object

moves from the initial location b0, at time to, to the final location bf, at time tf, at a constant velocity,

V0. The kinematically optimal trajectory between b0 and bf is a straight line; the traveled distance

along this line is R0. At time t, the object is at location b, which is not necessarily along this most

direct path. However, the remaining most direct path from b to bf is again given by a straight line

which connects these two locations. If the velocity vector, V, points from b to location bf, the object

would proceed from t until it reaches tf along this remaining most direct path. However, in order to

produce a remaining suboptimal trajectory, a disturbance is introduced by rotating V over the angle of

disturbance v, to obtain the disturbed vector, Vd. The angle v is chosen to be a simple exponential

function of the range R between b and bf, according to

tan v(t)/R(t)t e

=/R 0 ] tan v 0

(1)

where the exponent e determines the characteristics of the disturbance and therefore the trajectory

shape. The parameter vO, which is the angle of disturbance at to, determines the magnitude of the

disturbance.

For e > 1, v(t) will reduce to zero quickly, such that the trajectory will be curved mainly at the

beginning and straight toward the end, while for 0 < e _< 1 the opposite is the case. Typical trajectory

shapes for these two ranges of e are shown in figure 2.

2.2 Measures of"Suboplimality" of Motion

Various measures can be considered for specifying the "degree of suboptimality" of motion. The

first one is the difference between the traveled distance along the nondirect path, and the shortest distance, R0. For constant velocity V0, this is equivalent to the difference between the actual travel time

(tf- tO) and the most direct path travel time, R0/V0. Since this difference will be small relative to the

total travel distance (or travel time), this measure will be fairly insensitive for trajectories with small

deviations.

2
A secondcandidatefor ameasureof indirectness is the averaged angle of disturbance, computed by

Lf

- _ 1 v(t)dt (2)

v (tf- t0) to

This method has the disadvantage that only the disturbance and not the actual deviation from the straight

path is considered. Thus, for a straight trajectory section (like the final section of a trajectory with

e > 1), v is zero and does not add to the measure of indirectness, although the deviation from the most

direct path definitely exists (fig. 2a).

A third measure of indirectness is the averaged deviation from the straight path. The way in which

this score is computed is shown in figure 3. It is assumed that the object travels along the nondirect path

with a constant velocity V0. At each instance of time t, the distance d(t) between a point on the trajectory b and the equivalent position on the optimal path bs is computed; bs is the position on the

optimal path at which the object would arrive at time t when traveling with a reduced speed of

Vs = R0/(tf- tO) < V0. The positions b and bs are computed by solving the differential vector equations

R0

X-'b(t) = V--d(t) ; X-'bs(t) = (tf- t0)V 0 V0 (3)

with initial conditions

where

and the distance d(t) is

Xb(tO) = _Xbs(tO)= Xbo

xb = ix]; =

Y b Xbs Y bs

d(t) = IXb(t) - Xbs(t)l (4)

The averaged deviation from the optimal path over the interval tO < t < tf is then computed by

1 d(t)dt

a - (tf- to) to

(5)

3
3. DEVIATION FROM THE OPTIMAL TRAJECTORY FOR ROTATIONS

3.1 Definition of Coordinate Systems and Object Angular Rotation

The description of the angular rotation of an object in space involves rotational transformations

between coordinate systems. Two basic systems are defined: an Earth-based inertial reference system,

i, and an object body coordinate system, b. The Earth-based system, defined according to the aeronautical convention, is a fight-hand system with its xi-axis toward the north, its yi-axis to the east, and its

zi-axis pointing downward toward the center of Earth. In creating displays on the workstation, the

inertial system is chosen to coincide with the observer's viewing system, with the xi-axis being the

viewing direction. The object is described in the body system. For aerospace applications, for example,

this system is attached to the aircraft body, with the xb-axis pointing forward along the aircraft body

main axis, the yb-axis pointing outward through the right wing, and the zb-axis pointing downward.

(NOTE: Coordinate systems other than the North-East-Down convention are frequently used by nonaerospace disciplines. However, by means of simple geometrical computations, coordinates expressed

in one system can be transformed to another; the generality of the algorithm presented in this paper is

not affected.)

The angular orientation of the body system with respect to the Earth system can be described in a

variety of ways. The first one is by a set of three Euler angles, which define three successive rotations.

In aerospace applications these angles are the yaw angle _, which sets the heading plane; the pitch

angle 0, which sets the aircraft elevation angle with respect to the horizontal plane; and the roI1 angle 0,

in this order (fig. 4). A coordinate defined in the inertial system, x i - {xi,yi, zi} T, can now be expressed

in the rotated body system, x b - {xb,yb,zb} T, by

x b = [@] [0] [_] x i = D{_t,o,o }x i - Di___>b3. i (6)

where

[_] -

cgt sx¢ 0

-s_ cgt 0

0 0 1

; [0] -

cO 0 -sO

0 1 0

sO 0 cO I 1 0 0 ]

; [0]-- 0 c0 sO

0 -s0 c0

(7)

where c and s denote sine and cosine, respectively, and D is the nine-element direction cosine matrix

(DCM). The superscript T denotes transpose, since a column vector is specified. Likewise, by the

inverse transformation the coordinates of the b-system can be expressed in the i-system by

xi = [_]T [0IT [0I T x b = DT{_,0,q_ }x._b= Db__)ixb (8)

4
wherethe superscript T denotes the transposed matrix and, since D is an orthonormal matrix, transpose and inverse are identical.

Thus, the second way of expressing the angular orientation of the b-system in the i-system is by

means of the nine-element DCM. A third way is through the use of quaternions. Euler's theorem states

that, regardless of the initial orientation (defined by the system u) and the final orientation (defined by

the system v), it is always possible to find one axis about which the object can be rotated to bring it

from orientation u to v. The orientation of this axis of rotation, expressed in the initial system u, is

given by the unity vector e - {ex,ey,ez} T, and the amount of rotation by the angle t_. The quaternion

for rotation from u to v, which uniquely describes the orientation of v with respect to u, is given by

qu--_v _-

qt

q2

q3

q4

cos o_/2

e x sin _2

ey sin ct/2

e z sin _2

(9)

Since e is unity, it follows from equation (9) that

Iql= q +q2+q2+ =1 (10)

which is a unique property of the quaternion, and

1/2

Io0 = 2 cos-lql = 2 sinqlq2 + q2 + q_! ;

0<o_< 180 ° (11)

It also follows from equation (9) that the inverse quaternion, for rotation from v to u (i.e., the quaternion

which describes the orientation of the u-system in the v-system), is given by

-I

qv-+u = q-u-_v=

-ql

q2

q3

q4

(12)

which simply means that the direction of rotation about the Euler axis is inverted.
It isvery usefulto find thequatemionfor asequenceof rotations.Supposethattheangularorientationof systemv with respectto u is givenby quay and the orientation of a third system, w, with

respect to v by q__v_w. Then, the orientation of w with respect to u is given by

qu-,w = {Mv--)w}qu_v = (Mu_v} q_v_w (13)

with

ME

ql 42 43 44

q2 ' ql q4 -q3

I

q3 ,_14 ql q2

!

q4 , q3 -q2 ql

• M* -

q l -q2 -q3 -q4

q2 ' ql -q4 q3

!

q3 ' q4 ql -q2

!

q4 ,-q3 q2 ql

(14)

where M is the quatemion matrix, composed of the elements of q, and M* is the transmutated quaternion matrix, obtained by transposing the vector kernel (or minor) of the first element, indicated by the

dotted partition.

By using equation (13) in a sequence of rotations, the quaternion qi_b for rotation from inertial

system i to body system b can be found as a function of the Euler angles. Following equation (9) the

quaternions for the yaw, pitch, and roll rotations are given by

qi--')U =

cW2

0

0

sW2

c0/2

; q U----)V =

0

s0/2

0

s_12

; q__v_)b= (15)

0

0

respectively, where u and v denote the intermediate stages. Then, q._i_b is given by

, q i_b = {Mv_b} {Mu-)v} qi--)u (16)

and after evaluating equation (16) with equations (14) and (15)

qi---_b =

c¢/2 c0/2 c_/2 + sqb/2 sO/2

s¢/2 c0/2 cW2 - c0/2 sO/2

c¢/2 sO/2 c_/2 + s¢,/2 c0/2

-s¢/2 sO/2 c_1//2 + c0/2 c0/2

s_/2 ]

s_/2

s_/2

s_/2 ]

(17)
TheDCM Di-+b canbeexpressedin termsof qi_b as follows:

2 2 2

ql + q2 - q3 - q_

Di___)b= 2(q2q3 - qlq4)

2(q2q4 + qlq3)

2(q2q3 + qlq4)

ql - q2 + -

2(q4q3 - qlq2)

2(q2q4 - q 1q3)

2(q4q3 + qlq2)

ql- q -q +4

(18)

3.2 Rotational Equivalent of a Translatory Straight Path

It is clear from the definition of the Euler axis that the rotational analogue of a straight path between

two points in the translatory case is a rotation about the Euler axis. This rotation brings the object from

an initial orientation to a final one along the most direct path. The total angle of rotation about this axis

is equivalent to the distance between two points in the translatory case. Equivalently, this means that

the instantaneous axis of rotation (or vector of angular velocity) is fixed in space. This also means that

each coordinate of the object will describe a path, which will be the great-arc of a circle, located in a

plane perpendicular to this fixed axis of rotation. For points of the object located on this axis, the radius

of the circle will be zero. Mathematically, this means that the coordinates of these points are not

affected by the rotational transformation.

In contrast, when the rotation is not proceeding along the optimal trajectory, the instantaneous axis

of rotation will not be fixed in space and will perform a "wobbling" motion, somewhat like the nutation

of a gyroscope.

3.3 Computation of the Rotational Path

Similar to translatory motion, in which the path is obtained by integrating the instantaneous velocity

vector, the rotational path is obtained by integrating the components of the vector of angular velocity.

For the rotational motion of the body system b, this vector is given by

_b.__)b(t ) = {C0x,_y,_z} T (19)

where the subscript i--+b indicates the rotation of the b-system with respect to the i-system, and the

superscript b indicates that this vector is expressed in the rotating b-system. The rotational path is then

obtained by solving the differential equation

q_'i__)b(t)= {_(t) }q_i__b(t) ; qi_b(tO) = qi_b 0

(20)

where bo indicates the orientation of b at time tO and

7
f2(t) = 1/2- 1

0 -tox -toy -toz

tox 0 toz -t0y

COy -toz 0 COx

toz COy -fO x 0

(21)

3.4 Computation of the "Deviation" from the Optimal Rotational Trajectory

For the optimal rotational trajectory, the vector c0b_b(t) = tob_b(t 0) = 030 = constant will be fixed in

space and will coincide with the Euler axis. coo -10301 is the constant rotational speed and c_3 is the

total angle of rotation, equivalent to the shortest distance, R0, in the translatory case. The computation

of the "deviation" from the optimal trajectory is analogous to the translatory case (fig. 3). At each

instant of time the quaternion for the orientation b, qi_b(t), and the quatemion for the corresponding

orientation along the optimal trajectory bs, qi_bs(t), are computed, bs is the orientation on the optimal

trajectory which would be reached at time t- when rotating at a reduced angular rate of COs= o_0/(tf- tO)

< COo. The quatemion qi_bs(t) is obtained by solving the differential equation

cl_i_obs(t) = {f2s} q_i_bs (t) ; q_i-+bs(t0) = qi_+b 0

(22)

where f2s is fixed and generated according to equation (21) with

o_0

cos - coo (23)

(tf- t0)_ 0

The "rotational deviation" from the optimal trajectory (which is equivalent with the deviation d

from the straight path in the translatory case) is given by the amount of rotation [3 of the quaternion

Aq(t) = qbs_+b(t) = {Mi_b(t) }q_-lbs(t) (24)

-1

where qiob s is computed from qi--_bs by equation (12) and 13 is computed with equation (11). The

averaged deviation angle over the rotation interval tO < t < tf is then computed by

1 _3(t)dt (25)

_- (tf- to) to
3.5 Introduction of the Disturbance

Consider that the b-system rotates from initial orientation b0 at time to to final orientation bf at

time tf, not necessarily along the optimal trajectory. At time t the orientation of the b-system is b.

However, analogous to the translational case (in which the direction of the remaining most direct path is

V_.), in the rotational case the remaining most direct path from time t to tf onward is given by the Euleraxis rotation from b to bf

q_b_.)bf(t) = {Mi._)bf}qb_i(t) = {Mi_bf}qi-lb(t) (26)

-1

where qi_b is computed from q_i_b by equation (12), qi_b is computed by solving the differential

equation (eq. (20)), and Mi__.)bl is given by the final orientation of b in i. For 0 < ot < 180 °, the orientation of the Euler axis for rotation from b to bf (expressed in the b-system) is specified by

e(t)= q3/s

q4/s b----)bf

s = q2 + q3 + = sin(o_/2) (27)

and theb"angular distance" from b to bf is c_ = 2 sin-Is (eq. (11)). If the instantaneous vector of rotation _i__b(t) = o_b(t) (the subscript i---)b is henceforth omitted for clarity) would be chosen along e(t),

the rotation would proceed from t onward until it reaches tf along a remaining most direct path. However, in order to produce a remaining suboptimal trajectory, a disturbance is introduced to deviate o_b(t)

from e(t). As in the translatory case, the deviation of o)b(t) from e(t) should be a function of the angular distance o_(t). The vector 0_b(t) is deviated from e(t) by rotating it over the angle v(t), which is

chosen to be a simple exponential function of or(t) according to

tan v(t)=lOt(t)] e[ / tanv0 (28)

where v0 is the disturbance at time to, and e is an exponent. The parameter v0 determines the averaged deviation from the optimal trajectory and thus the magnitude of the disturbance, and the exponent

e determines the shape of the disturbance. Since b globally approaches bf, c_(t) will decrease monotonically such that the time-history of v(t) (and thus the shape of the trajectory) will be determined by

the exponent e. For e > 1, v(t) will reduce to zero quickly such that the trajectory will be "warped"

mainly at the beginning. For 0 < e < 1, the greatest warping occurs toward the end of the trajectory.

It should be noted that, although the angle of disturbance v(t) is defined, the plane in which Nb is

rotated is still undefined. A sequence of rotational transformations is required to choose a plane which

ensures a reproducible, smoothly warped rotational path. The method by which this plane is chosen is

described below.
The first rotationaltransformationthatis requiredto introducethedisturbanceinvolvesatransformationto the q0-system.Like the b0-system,the q0-systemis alsofixedwith respectto theinertial

system,but it hasthe xq0-axisalignedwith the _0-axis,andthe yq0-axisin the xb0-o-yb° plane

(fig. 5). Sincetheorientationof _ in the b0-systemis givenbyits azimuthangley andelevation

angle _5,thefixedDCM for rotationfrom b0to q0 is obtainedby asuccessiveyawrotationby the

angle Y andapitch rotationby theangle 5, following equations(6) and(7),accordingto

Db0___q 0 =

0 01 0 -sY c"[ 0

s_5 0 c_5 0 0 1

c8c3' cSsy -s5

-sy cy 0

sSc_' sSs_' c5

(29)

First, _b(t) = o_e(t) for the remaining most direct path from time t until tf, as computed by equation (27), is transformed from the b-system into the q0-system by

__qo(t) = [Db__x]ol_b(t) (30)

where

[Db__xt0] = [Db0_q0][Di___)b0][Db_i] (31)

The first two matrices on the right-hand side of equation (31) can be precomputed before starting the

motion-control program, and the third matrix is computed from qb_i - qiTb by solving equation (20)

and using equation (18).

The orientation of ¢oq0 in the q0-system is shown in figure 6. If the rotation were totally along the

optimal trajectory, coq0 would coincide with ._ and thus with the xq°-axis. However, its deviation

from this axis can be expressed in the q0-system by the azimuth angle Z and the elevation angle _. A

system ql is now defined, with the xql-axis along ¢.oq0 and the yql-axis in the xq0-o-yq 0 horizontal

plane. Similar to equation (29), the DCM for rotation from q0 to ql is given by

Dq0__q 1 =

-I

c cz c sz -st [

J

-sz cz 0

s cz sgsz

(32)

Since

10
COq0-{COx, COy, COz}T; I_q01= COO (33)

it follows from thegeometryof figure 6 that

cosX= C0x/d;

sinX = coy/d ;

(34)

where

d = COx2 + COy) 211/2 (35)

Substituting these expressions in equation (32) yields

cox/co0

Dq0___xt1 = -coy/d

-coxcoz/(dco0)

coy/C°0 coz/co0

(36)

Equations (35) and (36) show that Deo__Xl1 can be computed simply from the components of _q0 in the

q0-system.

Next, the q2-system is defined now rotated about the xql-axis by the angle Cd (fig. 7). This angle

is randomly chosen from a look-up table, but remains fixed throughout a trajectory calculation. The

DCM for rotation from ql to q2 is then given by

t 1 0 0

Dql__Xt2 = 0 C(_d SCd

0 -SCd CCd

(37)

The plane xq2-o-yq 2 is the required plane in which the vector co is rotated to be disturbed (fig. 7).

Last, the system q3 is defined, rotated about the zq2-axis by the angle of disturbance v(t), and computed by equation (28) (fig. 7). The DCM for rotation from q2 to q3 is then given by

11
Dq2_q3 =

-q

cv sv 0 |

-sv cv 0

0 0 1

(38)

With the DCMs previously computed with equations (31) and (36)-(38), the DCM for rotation from

b to q3 can be computed according to

[Db_xt3] = [Oq2__3][Dql__,Zl2][Dq0__,Xtl][Db_xt0] (39)

The xq3-axis is now aligned with the direction of the "disturbed" vector cob. This vector, expressed in

the b-system, is then computed by

cob = [Db__q3ITcoq3 (40)

where

¢0q3 = It00,0,0 ] T (41)

2.5 Computer Implementation

The motion-control program involves (1) precomputation of matrices, vectors, and quaternions,

which remain constant for all created trajectories, and (2) real-time computations for realizing the rotational path and computation of the deviations from the most direct path. This motion-control program is

incorporated into an experimental control program for psychophysical studies of human observers' sensitivity for trajectory optimality. This experimental control program presents the motion stimuli, records

observers' responses, and performs post-experiment data reduction.

A chart of the sequence of computations performed during the presentation of a rotational trajectory

stimulus in a typical psychophysical experiment is shown in figure 8. The prerun computations are

shown in bold blocks and the relevant input parameters for each stimulus are shown in the circles. For

each stimulus, the initial orientation b0 with respect to the inertial system i, the direction and magnitude of the most direct path vector of rotation COO,and the total angle of rotation v0 are specified.

The orientation of b0 in i is specified by the azimuth angle _, the elevation angle 0, and the roll

angle _ (fig. 4). The direction of ._ is specified with respect to the b0-system by the azimuth angle "/

and the elevation angle 5 (fig. 5). The magnitude of _ is the specified rotational speed, coo. For each

presentation the valuesof _,0,_ and 7,5 are picked at random and without replacement from a look-up

table.

12
First, thequaternionqi--)b 0 is precomputed with _, 0, and 0, using equation (17). The DCM Di--)b0

is computed simply with equation (18). Next, following the geometry of figure 5, from yand 8, the

components of o0 in the b0-system can be precomputed according to

c?c5]

co0 = cools7c8 t

L-sSj

(42)

Since cob - coO coincides with the Euler axis for rotation from b0 to bf and thus has the same direction,

the fixed quaternion for rotation from b0 to bf follows from equations (42) and (9):

q_-b0--)b f =

ca0/2-1

cycgso_0/2 l

s"/c_Ssa0/2 [

-s_Ssa0/2.J

(43)

and the corresponding quatemion matrix M b ---)b is computed with equation (14). The latter one, 0 f

together with the quaternion of the initial orientation qi_b 0 can be used to precompute the quaternion

of the final orientation

qi_bf = {Mb0---_f}qi_b 0

(44)

and the corresponding quatemion matrix Mi_bf is computed with equation (14).

Last, Db0_q 0 is precomputed with ',/and _5,using equation (29); D ,---xt., with _d, using equation (37) where 0d is picked at random from a table; and c0q 3 is defined _']th equation (41). The

matrix f2s is computed using equation (21), with cos computed with equation (23), and _o-)0is

computed with equation (42).

The on-line computations during stimulus presentation are aimed at computing the "disturbed" vector of angular rotation ob(t) with equation (40), which constitutes the vector of rotation at which the

object proceeds rotating at time t. This requires the computation of cob(t) - c00e(t) with equations (26)

and (27), which is the undisturbed vector of rotation for the remaining most direct path between orientation b and bf, and the computation of various DCMs for transffrmations between the b-, i-, q0-, ql-,

q2-, and q3-systems. It also requires the computation of the "at" gular distance" from b to bf, oc(t) with

equation (11). With c_(t), the angle of disturbance v(t) is corn 'uted with equation (28). The disturbed

vector of rotation is used to compute f_(t) with equation (21), _fhich in its turn, is used to compute the

rotational path, described by the quaternion qi_b(t) by solving the differential equation of equation

(2O).

13
It canbeshownthattheaveragedrotationaldeviationfrom themostdirectpath,13,andtheduration

of therotation,normalizedwith respectto thedurationof themostdirectpath,(tf - t0)co0/ct0, dependon

theparametersv0, e,ando_0only, andarenotaffectedby theinitial orientation,theorientationof the

Euleraxis,or theangleof theplaneof disturbance,00d. Hence these values can be precomputed and

tabulated. For each set of v0, e, and o% the rotational path is preexecuted twice: the first time to compute t_he time interval (tf- tO) and the second time, using this interval, to compute the averaged deviation, [3.

At run time the relevant parameter for specifying the magnitude of distortion of the rotational path is

13 rather thany0. Therefore, before starting the run, the value of v0 for realizing the trajectory with

the specified 13 is computed by linear interpolation from the tabulated values. In order to verify that the

specified trajectory was executed, the actual rotational deviation 13 is computed on line, and at the end

of the presentation of the rotational path, the averaged value is compared with the specified one.

4. CONCLUSION

This report has presented a systematic method for creating suboptimal rotational trajectories. Two

parameters of the distortion metric, v0 and e, independently control the magnitude of the distortion and

its "shape" (i.e., how the distortion is distributed over the time course of the trajectory). These systematic disturbances can be introduced into a motion-control program to examine human observers' sensitivity to rotational trajectory optimality. In addition, this algorithm can be used in other applications

requiring systematic and reproducible disturbances of rotational trajectories.

14
o

1 -Vd /

t v(t)

b

V

NON-DIRECT

PATH

FINAL

LOCATION

bf

/%

/ // // _Vd

/

\

\

\

\

(to) /

INITIAL

LOCATION

\

\

REMAINING

MOST DI RECT

PATH (t _ tf)

MOST DIRECT PATH (to -+ tf)

Y

Figure 1. Method for creating a suboptimal trajectory for translational motion.

15
(a)

0

e>l

•YO

b0 (t O)

INITIAL

LOCATION

J

J

J

J

FINAL

LOCATION

J bf (tf)

¥

]1

(b)

O<e _<1

b O (t0)

INITIAL

LOCATION y

Figure 2. The effect of the exponent e on the trajectory shape.

16
x

o

NON - DIRECT PATH,

TRAVELED AT

SPEED V 0

(t)

b

I

I

I

I

I

I

I

/

/

/

/ i

/ /

/ /

/ /

I /

/

bs (t) FINAL

1 I -V0 R0 LOCATION

0_/ _ DIRECT PATH TRAVELED

b0(t _ AT SPEEDR 0/(tf-t o )

INITIAL

LOCATION

Figure 3. Method for computing the averaged deviation from the most direct path.

17
\

\

\

\

\\

\

\

\\ I1I
xb0-0-y b0 PLANE

xq0

EULER AXIS b0-_ bf

xb0

yb0

yq0

zb0

Figure 5. Definition of the q0-system.

19
-(.d Z

\

xqO.o-y qO PLANE

q/

/

/ \

/,, _,/" _ zqI

zqO

Figure 6. Definition of the ql-system.

-->bf

yqO

2O
xql, xq2

qo

EULER AXIS b --* bf

• xq3

_ DIRECTION OF

DISTURBED

VECTOR OF

ROTATION

o

yql

xq2-o-y q2 PLANE

yq2

zq2

zql

Figure 7. Definition of the q2 and q3-systems.

21
I PRECOMPUTE qi_bo EQ. ] (17)

PRECOMPUTE H PRECOMPUTE I

qi_bf EQ. (44) Mi_"bf EQ. (14) I

Li P"EOOMPUTE I l

PRECOMPUTE

-qb0_bf EQ (43)

Mb0--bf EQ'(14)l I " I

PRECOMPUTE COMPUTE

Di-'b0 EQ. (18) _---_-J _ _ I Db-_i EQS. (18)(12)

PRECOMPUTE l ' _! C:M;:_E E

Db0_'q0 EQ. (29) Q (31)

SET INITIAL CONDITIONS

-qi-,b (to) = qi_b s (to) = qi_b 0

+qi-- b s (t)

_-- COMPUTE q_b_bf EQ. (26)

f

COMPUTE

e

- EQ. (27)

{

COMPUTE

__b = _0_e

COMPUTE __qo EQ. (30) 1

COMPUTE

(35)

Dq0_q I EQS. (36)

"= COMPUTE Db-_q3 EQ. (391

COMPUTE

EQ. (11)

COMPUTE u EQ. (28)

COMPUTE

Dq2-_q3 EQ. (38)

COMPUTE Dql_q2 EQ. (37) !_ _

(PART 1 OF 2)

Figure 8. Flow chart Of sequence of computations.

!:

22
?

I

I

?

PRECOMPUTE _b EQ. (42) {

f

PRECOMPUTE -_s EQ. (23) I

+

PRECOMPUTE _s EQ. 121) I

f

COMPUTE

q'-*bs EQ. (22)

RECORD AND

VERIFY

_, tf-t 0

{

J

qi -* b$ (t)_ qi -_ b (t)

COMPUTE Mi-*b EQ. (14)

+

COMPUTE

_q EO. (24)

COMPUTE ] /_ EQ. (11)

COMPUTE

I

I

EQ. (25)

NUMERICAL INTEGRATION

?

COMPUTE H

b . cuq3

-_d EQ. (40) -

f

COMPUTE

EQ. (21)

t

COMPUTE

-qi-'b EQ. (20)

i

1

GO TO NEXT ITERATION

WITH UPDATED VALUES

OF qi_b AND qi_,bs

PRECOMPUTE

EQ. (41)

(PART 2 OF 2)

Figure 8. Concluded.

23
,q ii i

NASA Ne tional Aeronautics end

Space Admin/$tatJ_

1. Report No.

NASA TM- 101089

Report Documentation Page

2. GoVernment' Accession No. r ,

4. Title and Subtitle

An Algorithm for the Systematic Disturbance of

Optimal Rotational Solutions

7. Author(s)

Arthur J. Grunwald (Technion, Haifa, israel) and

Mary K. Kaiser

9, Performing Organization Name and Address

Ames Research Center

Moffett Field, CA 94035

12. Sponsoring Agency Name and Address

National Aeronautics and Space Administration

Washington, DC 20546-0001

3. Recipient's Catalog No.

,5.Rel:_h_Date

August 1989

6. Performing Organ'ization Code

8. Performing Organization Report No.

A-89092

10. Work Unit NO. "

506-47-11

11'.' (_ntraet or'Grant No,

13. Type of Report and Period Covered

Technical Memorandum

14. Sponsoring Agency Code

15. Supplementary Notes

Point of Contact: Mary K. Kaiser, Ames Research Center, MS 239-3, Moffett Field, CA 94035

(415) 694-6243 or FTS 464-6243

16. Abstract

An algorithm for introducing a systematic rotational disturbance into an optimal (i.e., single

axis) rotational trajectory is described. This disturbance introduces a motion vector orthogonal to the

quaternion-defined optimal rotation axis. By altering the magnitude of this vector, the degree of nonoptimality can be controlled. The metric properties of the distortion parameter are described, with

analogies to two-dimensional translational motion.

This algorithm has been implemented in a motion-control program on a three-dimensional

graphic workstation. It supports a series of human performance studies on the detectability of rotational trajectory optimality by naive observers.

17. Key Words (Suggested by Author(s))

Quaternion

Systemic distortions

Rotational solutions

18. Distribution Statement

Unclassified-Unlimited

Subject Category - 59

19. Security Classif. (of this report)

Unclassified

20. Security Classif. (of this page)

Unclassified

21. No. of Pages

28

22. Price

A03

NASA FORM 1626 0CT86

For sale by the National Technical Information Service, Springfield, Virginia 22161