Loading [MathJax]/extensions/tex2jax.js
JSBSim Flight Dynamics Model 1.2.2 (22 Mar 2025)
An Open Source Flight Dynamics and Control Software Library in C++
All Classes Functions Variables Enumerations Enumerator Friends Pages
FGQuaternion Class Reference

Detailed Description

Models the Quaternion representation of rotations.

FGQuaternion is a representation of an arbitrary rotation through a quaternion. It has vector properties. This class also contains access functions to the euler angle representation of rotations and access to transformation matrices for 3D vectors. Transformations and euler angles are therefore computed once they are requested for the first time. Then they are cached for later usage as long as the class is not accessed trough a nonconst member function.

Note: The order of rotations used in this class corresponds to a 3-2-1 sequence, or Y-P-R, or Z-Y-X, if you prefer.

See also
Cooke, Zyda, Pratt, and McGhee, "NPSNET: Flight Simulation Dynamic Modeling Using Quaternions", Presence, Vol. 1, No. 4, pp. 404-420 Naval Postgraduate School, January 1994
D. M. Henderson, "Euler Angles, Quaternions, and Transformation Matrices", JSC 12960, July 1977
Richard E. McFarland, "A Standard Kinematic Model for Flight Simulation at NASA-Ames", NASA CR-2497, January 1975
Barnes W. McCormick, "Aerodynamics, Aeronautics, and Flight Mechanics", Wiley & Sons, 1979 ISBN 0-471-03032-5
Bernard Etkin, "Dynamics of Flight, Stability and Control", Wiley & Sons, 1982 ISBN 0-471-08936-2
Author
Mathias Froehlich, extended FGColumnVector4 originally by Tony Peden and Jon Berndt

Definition at line 86 of file FGQuaternion.h.

#include <FGQuaternion.h>

+ Inheritance diagram for FGQuaternion:
+ Collaboration diagram for FGQuaternion:

Public Member Functions

 FGQuaternion ()
 Default initializer.
 
 FGQuaternion (const FGMatrix33 &m)
 Initializer by matrix.
 
 FGQuaternion (const FGQuaternion &q)
 Copy constructor.
 
 FGQuaternion (double angle, const FGColumnVector3 &axis)
 Initializer by a rotation axis and an angle.
 
 FGQuaternion (double phi, double tht, double psi)
 Initializer by euler angles.
 
 FGQuaternion (FGColumnVector3 vOrient)
 Initializer by euler angle vector.
 
 FGQuaternion (int idx, double angle)
 Initializer by one euler angle.
 
 ~FGQuaternion ()
 Destructor.
 
FGQuaternion Conjugate (void) const
 Conjugate of the quaternion.
 
std::string Dump (const std::string &delimiter) const
 
double & Entry (unsigned int idx)
 Write access the entries of the vector.
 
double Entry (unsigned int idx) const
 Read access the entries of the vector.
 
double GetCosEuler (int i) const
 Retrieves cosine of the given euler angle.
 
double GetEuler (int i) const
 Retrieves the Euler angles.
 
const FGColumnVector3GetEuler (void) const
 Retrieves the Euler angles.
 
double GetEulerDeg (int i) const
 Retrieves the Euler angles.
 
FGColumnVector3 const GetEulerDeg (void) const
 Retrieves the Euler angle vector.
 
FGQuaternion GetQDot (const FGColumnVector3 &PQR) const
 Quaternion derivative for given angular rates.
 
double GetSinEuler (int i) const
 Retrieves sine of the given euler angle.
 
const FGMatrix33GetT (void) const
 Transformation matrix.
 
const FGMatrix33GetTInv (void) const
 Backward transformation matrix.
 
FGQuaternion Inverse (void) const
 Inverse of the quaternion.
 
double Magnitude (void) const
 Length of the vector.
 
void Normalize (void)
 Normalize.
 
 operator FGMatrix33 () const
 Conversion from Quat to Matrix.
 
bool operator!= (const FGQuaternion &q) const
 Comparison operator "!=".
 
double & operator() (unsigned int idx)
 Write access the entries of the vector.
 
double operator() (unsigned int idx) const
 Read access the entries of the vector.
 
FGQuaternion operator* (const FGQuaternion &q) const
 Arithmetic operator "*".
 
const FGQuaternionoperator*= (const FGQuaternion &q)
 Arithmetic operator "*=".
 
const FGQuaternionoperator*= (double scalar)
 Arithmetic operator "*=".
 
FGQuaternion operator+ (const FGQuaternion &q) const
 Arithmetic operator "+".
 
const FGQuaternionoperator+= (const FGQuaternion &q)
 
FGQuaternion operator- (const FGQuaternion &q) const
 Arithmetic operator "-".
 
const FGQuaternionoperator-= (const FGQuaternion &q)
 Arithmetic operator "-=".
 
const FGQuaternionoperator/= (double scalar)
 Arithmetic operator "/=".
 
const FGQuaternionoperator= (const FGQuaternion &q)
 Assignment operator "=".
 
bool operator== (const FGQuaternion &q) const
 Comparison operator "==".
 
double SqrMagnitude (void) const
 Square of the length of the vector.
 
- Public Member Functions inherited from FGJSBBase
 FGJSBBase ()
 Constructor for FGJSBBase.
 
virtual ~FGJSBBase ()
 Destructor for FGJSBBase.
 
void disableHighLighting (void)
 Disables highlighting in the console output.
 

Static Public Member Functions

static FGQuaternion zero (void)
 Zero quaternion vector.
 
- Static Public Member Functions inherited from FGJSBBase
static const std::string & GetVersion (void)
 Returns the version number of JSBSim.
 
static constexpr double KelvinToFahrenheit (double kelvin)
 Converts from degrees Kelvin to degrees Fahrenheit.
 
static constexpr double CelsiusToRankine (double celsius)
 Converts from degrees Celsius to degrees Rankine.
 
static constexpr double RankineToCelsius (double rankine)
 Converts from degrees Rankine to degrees Celsius.
 
static constexpr double KelvinToRankine (double kelvin)
 Converts from degrees Kelvin to degrees Rankine.
 
static constexpr double RankineToKelvin (double rankine)
 Converts from degrees Rankine to degrees Kelvin.
 
static constexpr double FahrenheitToCelsius (double fahrenheit)
 Converts from degrees Fahrenheit to degrees Celsius.
 
static constexpr double CelsiusToFahrenheit (double celsius)
 Converts from degrees Celsius to degrees Fahrenheit.
 
static constexpr double CelsiusToKelvin (double celsius)
 Converts from degrees Celsius to degrees Kelvin.
 
static constexpr double KelvinToCelsius (double kelvin)
 Converts from degrees Kelvin to degrees Celsius.
 
static constexpr double FeetToMeters (double measure)
 Converts from feet to meters.
 
static bool EqualToRoundoff (double a, double b)
 Finite precision comparison.
 
static bool EqualToRoundoff (float a, float b)
 Finite precision comparison.
 
static bool EqualToRoundoff (float a, double b)
 Finite precision comparison.
 
static bool EqualToRoundoff (double a, float b)
 Finite precision comparison.
 
static constexpr double Constrain (double min, double value, double max)
 Constrain a value between a minimum and a maximum value.
 
static constexpr double sign (double num)
 

Friends

FGQuaternion operator* (double, const FGQuaternion &)
 Scalar multiplication.
 
FGQuaternion QExp (const FGColumnVector3 &omega)
 Quaternion exponential.
 

Additional Inherited Members

- Public Types inherited from FGJSBBase
enum  { eL = 1 , eM , eN }
 Moments L, M, N. More...
 
enum  { eP = 1 , eQ , eR }
 Rates P, Q, R. More...
 
enum  { eU = 1 , eV , eW }
 Velocities U, V, W. More...
 
enum  { eX = 1 , eY , eZ }
 Positions X, Y, Z. More...
 
enum  { ePhi = 1 , eTht , ePsi }
 Euler angles Phi, Theta, Psi. More...
 
enum  { eDrag = 1 , eSide , eLift }
 Stability axis forces, Drag, Side force, Lift. More...
 
enum  { eRoll = 1 , ePitch , eYaw }
 Local frame orientation Roll, Pitch, Yaw. More...
 
enum  { eNorth = 1 , eEast , eDown }
 Local frame position North, East, Down. More...
 
enum  { eLat = 1 , eLong , eRad }
 Locations Radius, Latitude, Longitude. More...
 
enum  {
  inNone = 0 , inDegrees , inRadians , inMeters ,
  inFeet
}
 Conversion specifiers. More...
 
- Static Public Attributes inherited from FGJSBBase
static char highint [5] = {27, '[', '1', 'm', '\0' }
 highlights text
 
static char halfint [5] = {27, '[', '2', 'm', '\0' }
 low intensity text
 
static char normint [6] = {27, '[', '2', '2', 'm', '\0' }
 normal intensity text
 
static char reset [5] = {27, '[', '0', 'm', '\0' }
 resets text properties
 
static char underon [5] = {27, '[', '4', 'm', '\0' }
 underlines text
 
static char underoff [6] = {27, '[', '2', '4', 'm', '\0' }
 underline off
 
static char fgblue [6] = {27, '[', '3', '4', 'm', '\0' }
 blue text
 
static char fgcyan [6] = {27, '[', '3', '6', 'm', '\0' }
 cyan text
 
static char fgred [6] = {27, '[', '3', '1', 'm', '\0' }
 red text
 
static char fggreen [6] = {27, '[', '3', '2', 'm', '\0' }
 green text
 
static char fgdef [6] = {27, '[', '3', '9', 'm', '\0' }
 default text
 
static short debug_lvl = 1
 
- Static Protected Member Functions inherited from FGJSBBase
static std::string CreateIndexedPropertyName (const std::string &Property, int index)
 
- Static Protected Attributes inherited from FGJSBBase
static constexpr double radtodeg = 180. / M_PI
 
static constexpr double degtorad = M_PI / 180.
 
static constexpr double hptoftlbssec = 550.0
 
static constexpr double psftoinhg = 0.014138
 
static constexpr double psftopa = 47.88
 
static constexpr double fttom = 0.3048
 
static constexpr double ktstofps = 1852./(3600*fttom)
 
static constexpr double fpstokts = 1.0 / ktstofps
 
static constexpr double inchtoft = 1.0/12.0
 
static constexpr double m3toft3 = 1.0/(fttom*fttom*fttom)
 
static constexpr double in3tom3 = inchtoft*inchtoft*inchtoft/m3toft3
 
static constexpr double inhgtopa = 3386.38
 
static constexpr double slugtolb = 32.174049
 Note that definition of lbtoslug by the inverse of slugtolb and not to a different constant you can also get from some tables will make lbtoslug*slugtolb == 1 up to the magnitude of roundoff.
 
static constexpr double lbtoslug = 1.0/slugtolb
 
static constexpr double kgtolb = 2.20462
 
static constexpr double kgtoslug = 0.06852168
 
static const std::string needed_cfg_version = "2.0"
 
static const std::string JSBSim_version = JSBSIM_VERSION " " __DATE__ " " __TIME__
 

Constructor & Destructor Documentation

◆ FGQuaternion() [1/7]

FGQuaternion ( )
inline

Default initializer.

Default initializer, initializes the class with the identity rotation.

Definition at line 90 of file FGQuaternion.h.

90 : mCacheValid(false) {
91 data[0] = 1.0;
92 data[1] = data[2] = data[3] = 0.0;
93 }
+ Here is the caller graph for this function:

◆ FGQuaternion() [2/7]

FGQuaternion ( const FGQuaternion q)

Copy constructor.

Copy constructor, initializes the quaternion.

Parameters
qa constant reference to another FGQuaternion instance

Definition at line 64 of file FGQuaternion.cpp.

64 : mCacheValid(q.mCacheValid)
65{
66 data[0] = q(1);
67 data[1] = q(2);
68 data[2] = q(3);
69 data[3] = q(4);
70 if (mCacheValid) {
71 mT = q.mT;
72 mTInv = q.mTInv;
73 mEulerAngles = q.mEulerAngles;
74 mEulerSines = q.mEulerSines;
75 mEulerCosines = q.mEulerCosines;
76 }
77}

◆ FGQuaternion() [3/7]

FGQuaternion ( double  phi,
double  tht,
double  psi 
)

Initializer by euler angles.

Initialize the quaternion with the euler angles.

Parameters
phiThe euler X axis (roll) angle in radians
thtThe euler Y axis (attitude) angle in radians
psiThe euler Z axis (heading) angle in radians

Definition at line 82 of file FGQuaternion.cpp.

82 : mCacheValid(false)
83{
84 InitializeFromEulerAngles(phi, tht, psi);
85}

◆ FGQuaternion() [4/7]

Initializer by euler angle vector.

Initialize the quaternion with the euler angle vector.

Parameters
vOrientThe euler axis angle vector in radians (phi, tht, psi)

Definition at line 89 of file FGQuaternion.cpp.

89 : mCacheValid(false)
90{
91 double phi = vOrient(ePhi);
92 double tht = vOrient(eTht);
93 double psi = vOrient(ePsi);
94
95 InitializeFromEulerAngles(phi, tht, psi);
96}

◆ FGQuaternion() [5/7]

FGQuaternion ( int  idx,
double  angle 
)
inline

Initializer by one euler angle.

Initialize the quaternion with the single euler angle where its index is given in the first argument.

Parameters
idxIndex of the euler angle to initialize
angleThe euler angle in radians

Definition at line 117 of file FGQuaternion.h.

118 : mCacheValid(false) {
119
120 double angle2 = 0.5*angle;
121
122 double Sangle2 = sin(angle2);
123 double Cangle2 = cos(angle2);
124
125 if (idx == ePhi) {
126 data[0] = Cangle2;
127 data[1] = Sangle2;
128 data[2] = 0.0;
129 data[3] = 0.0;
130
131 } else if (idx == eTht) {
132 data[0] = Cangle2;
133 data[1] = 0.0;
134 data[2] = Sangle2;
135 data[3] = 0.0;
136
137 } else {
138 data[0] = Cangle2;
139 data[1] = 0.0;
140 data[2] = 0.0;
141 data[3] = Sangle2;
142
143 }
144 }

◆ FGQuaternion() [6/7]

FGQuaternion ( double  angle,
const FGColumnVector3 axis 
)
inline

Initializer by a rotation axis and an angle.

Initialize the quaternion to represent the rotation around a given angle and an arbitrary axis.

Parameters
angleThe angle in radians
axisThe rotation axis

Definition at line 152 of file FGQuaternion.h.

153 : mCacheValid(false) {
154
155 double angle2 = 0.5 * angle;
156
157 double length = axis.Magnitude();
158 double Sangle2 = sin(angle2) / length;
159 double Cangle2 = cos(angle2);
160
161 data[0] = Cangle2;
162 data[1] = Sangle2 * axis(1);
163 data[2] = Sangle2 * axis(2);
164 data[3] = Sangle2 * axis(3);
165 }
+ Here is the call graph for this function:

◆ FGQuaternion() [7/7]

FGQuaternion ( const FGMatrix33 m)

Initializer by matrix.

Initialize the quaternion with the matrix representing a transform from one frame to another using the standard aerospace sequence, Yaw-Pitch-Roll (3-2-1).

Parameters
mthe rotation matrix

Definition at line 139 of file FGQuaternion.cpp.

139 : mCacheValid(false)
140{
141 data[0] = 0.50*sqrt(1.0 + m(1,1) + m(2,2) + m(3,3));
142 double t = 0.25/data[0];
143 data[1] = t*(m(2,3) - m(3,2));
144 data[2] = t*(m(3,1) - m(1,3));
145 data[3] = t*(m(1,2) - m(2,1));
146
147 Normalize();
148}
void Normalize(void)
Normalize.
+ Here is the call graph for this function:

◆ ~FGQuaternion()

~FGQuaternion ( )
inline

Destructor.

Definition at line 174 of file FGQuaternion.h.

174{}

Member Function Documentation

◆ Conjugate()

FGQuaternion Conjugate ( void  ) const
inline

Conjugate of the quaternion.

Compute and return the conjugate of the quaternion. This one is equal to the inverse iff the quaternion is normalized.

Definition at line 450 of file FGQuaternion.h.

450 {
451 return FGQuaternion( data[0], -data[1], -data[2], -data[3] );
452 }
FGQuaternion()
Default initializer.

◆ Dump()

std::string Dump ( const std::string &  delimiter) const

Definition at line 239 of file FGQuaternion.cpp.

240{
241 std::ostringstream buffer;
242 buffer << std::setprecision(16) << data[0] << delimiter;
243 buffer << std::setprecision(16) << data[1] << delimiter;
244 buffer << std::setprecision(16) << data[2] << delimiter;
245 buffer << std::setprecision(16) << data[3];
246 return buffer.str();
247}

◆ Entry() [1/2]

double & Entry ( unsigned int  idx)
inline

Write access the entries of the vector.

Parameters
idxthe component index.

Return a reference to the vector entry at the given index. Indices are counted starting with 1.

This function is just a shortcut for the double& operator()(unsigned int idx) function. It is used internally to access the elements in a more convenient way.

Note that the index given in the argument is unchecked.

Definition at line 300 of file FGQuaternion.h.

300 {
301 mCacheValid = false;
302 return data[idx-1];
303 }

◆ Entry() [2/2]

double Entry ( unsigned int  idx) const
inline

Read access the entries of the vector.

Parameters
idxthe component index.

Return the value of the matrix entry at the given index. Indices are counted starting with 1.

This function is just a shortcut for the double operator()(unsigned int idx) const function. It is used internally to access the elements in a more convenient way.

Note that the index given in the argument is unchecked.

Definition at line 285 of file FGQuaternion.h.

285{ return data[idx-1]; }

◆ GetCosEuler()

double GetCosEuler ( int  i) const
inline

Retrieves cosine of the given euler angle.

Returns
the sine of the Euler angle theta (pitch attitude) corresponding to this quaternion rotation.

Definition at line 245 of file FGQuaternion.h.

245 {
246 ComputeDerived();
247 return mEulerCosines(i);
248 }
+ Here is the caller graph for this function:

◆ GetEuler() [1/2]

double GetEuler ( int  i) const
inline

Retrieves the Euler angles.

Parameters
ithe Euler angle index. units radians.
Returns
a reference to the i-th euler angles corresponding to this quaternion rotation.

Definition at line 210 of file FGQuaternion.h.

210 {
211 ComputeDerived();
212 return mEulerAngles(i);
213 }

◆ GetEuler() [2/2]

const FGColumnVector3 & GetEuler ( void  ) const
inline

Retrieves the Euler angles.

Returns
a reference to the triad of Euler angles corresponding to this quaternion rotation. units radians

Definition at line 199 of file FGQuaternion.h.

199 {
200 ComputeDerived();
201 return mEulerAngles;
202 }
+ Here is the caller graph for this function:

◆ GetEulerDeg() [1/2]

double GetEulerDeg ( int  i) const
inline

Retrieves the Euler angles.

Parameters
ithe Euler angle index.
Returns
a reference to the i-th euler angles corresponding to this quaternion rotation. units degrees

Definition at line 220 of file FGQuaternion.h.

220 {
221 ComputeDerived();
222 return radtodeg*mEulerAngles(i);
223 }

◆ GetEulerDeg() [2/2]

FGColumnVector3 const GetEulerDeg ( void  ) const
inline

Retrieves the Euler angle vector.

Returns
an Euler angle column vector corresponding to this quaternion rotation. units degrees

Definition at line 229 of file FGQuaternion.h.

229 {
230 ComputeDerived();
231 return radtodeg*mEulerAngles;
232 }

◆ GetQDot()

FGQuaternion GetQDot ( const FGColumnVector3 PQR) const

Quaternion derivative for given angular rates.

Returns the derivative of the quaternion corresponding to the angular velocities PQR.

Computes the quaternion derivative which results from the given angular velocities

Parameters
PQRa constant reference to a rotation rate vector
Returns
the quaternion derivative
See also
Stevens and Lewis, "Aircraft Control and Simulation", Second Edition, Equation 1.3-36.

See Stevens and Lewis, "Aircraft Control and Simulation", Second Edition, Equation 1.3-36. Also see Jack Kuipers, "Quaternions and Rotation Sequences", Equation 11.12.

Definition at line 158 of file FGQuaternion.cpp.

159{
160 return FGQuaternion(
161 -0.5*( data[1]*PQR(eP) + data[2]*PQR(eQ) + data[3]*PQR(eR)),
162 0.5*( data[0]*PQR(eP) - data[3]*PQR(eQ) + data[2]*PQR(eR)),
163 0.5*( data[3]*PQR(eP) + data[0]*PQR(eQ) - data[1]*PQR(eR)),
164 0.5*(-data[2]*PQR(eP) + data[1]*PQR(eQ) + data[0]*PQR(eR))
165 );
166}
+ Here is the call graph for this function:

◆ GetSinEuler()

double GetSinEuler ( int  i) const
inline

Retrieves sine of the given euler angle.

Returns
the sine of the Euler angle theta (pitch attitude) corresponding to this quaternion rotation.

Definition at line 237 of file FGQuaternion.h.

237 {
238 ComputeDerived();
239 return mEulerSines(i);
240 }
+ Here is the caller graph for this function:

◆ GetT()

const FGMatrix33 & GetT ( void  ) const
inline

Transformation matrix.

Returns
a reference to the transformation/rotation matrix corresponding to this quaternion rotation.

Definition at line 188 of file FGQuaternion.h.

188{ ComputeDerived(); return mT; }
+ Here is the caller graph for this function:

◆ GetTInv()

const FGMatrix33 & GetTInv ( void  ) const
inline

Backward transformation matrix.

Returns
a reference to the inverse transformation/rotation matrix corresponding to this quaternion rotation.

Definition at line 193 of file FGQuaternion.h.

193{ ComputeDerived(); return mTInv; }
+ Here is the caller graph for this function:

◆ Inverse()

FGQuaternion Inverse ( void  ) const
inline

Inverse of the quaternion.

Compute and return the inverse of the quaternion so that the orientation represented with *this multiplied with the returned value is equal to the identity orientation.

Definition at line 436 of file FGQuaternion.h.

436 {
437 double norm = SqrMagnitude();
438 if (norm == 0.0)
439 return *this;
440 double rNorm = 1.0/norm;
441 return FGQuaternion( data[0]*rNorm, -data[1]*rNorm,
442 -data[2]*rNorm, -data[3]*rNorm );
443 }
double SqrMagnitude(void) const
Square of the length of the vector.

◆ Magnitude()

double Magnitude ( void  ) const
inline

Length of the vector.

Compute and return the euclidean norm of this vector.

Definition at line 460 of file FGQuaternion.h.

460{ return sqrt(SqrMagnitude()); }
+ Here is the caller graph for this function:

◆ Normalize()

void Normalize ( void  )

Normalize.

Normalize the vector to have the Magnitude() == 1.0. If the vector is equal to zero it is left untouched.

Definition at line 170 of file FGQuaternion.cpp.

171{
172 // Note: this does not touch the cache since it does not change the orientation
173 double norm = Magnitude();
174 if (norm == 0.0 || fabs(norm - 1.000) < 1e-10) return;
175
176 double rnorm = 1.0/norm;
177
178 data[0] *= rnorm;
179 data[1] *= rnorm;
180 data[2] *= rnorm;
181 data[3] *= rnorm;
182}
double Magnitude(void) const
Length of the vector.
+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ operator FGMatrix33()

operator FGMatrix33 ( ) const
inline

Conversion from Quat to Matrix.

Definition at line 330 of file FGQuaternion.h.

330{ return GetT(); }
const FGMatrix33 & GetT(void) const
Transformation matrix.

◆ operator!=()

bool operator!= ( const FGQuaternion q) const
inline

Comparison operator "!=".

Parameters
qa quaternion reference
Returns
true if both quaternions do not represent the same rotation.

Definition at line 343 of file FGQuaternion.h.

343{ return ! operator==(q); }
bool operator==(const FGQuaternion &q) const
Comparison operator "==".

◆ operator()() [1/2]

double & operator() ( unsigned int  idx)
inline

Write access the entries of the vector.

Parameters
idxthe component index.

Return a reference to the vector entry at the given index. Indices are counted starting with 1.

Note that the index given in the argument is unchecked.

Definition at line 270 of file FGQuaternion.h.

270{ mCacheValid = false; return data[idx-1]; }

◆ operator()() [2/2]

double operator() ( unsigned int  idx) const
inline

Read access the entries of the vector.

Parameters
idxthe component index.

Return the value of the matrix entry at the given index. Indices are counted starting with 1.

Note that the index given in the argument is unchecked.

Definition at line 259 of file FGQuaternion.h.

259{ return data[idx-1]; }

◆ operator*()

FGQuaternion operator* ( const FGQuaternion q) const
inline

Arithmetic operator "*".

Multiplication of two quaternions is like performing successive rotations.

Parameters
qa quaternion to be multiplied.
Returns
a quaternion representing Q, where Q = Q * q.

Definition at line 406 of file FGQuaternion.h.

406 {
407 return FGQuaternion(data[0]*q.data[0]-data[1]*q.data[1]-data[2]*q.data[2]-data[3]*q.data[3],
408 data[0]*q.data[1]+data[1]*q.data[0]+data[2]*q.data[3]-data[3]*q.data[2],
409 data[0]*q.data[2]-data[1]*q.data[3]+data[2]*q.data[0]+data[3]*q.data[1],
410 data[0]*q.data[3]+data[1]*q.data[2]-data[2]*q.data[1]+data[3]*q.data[0]);
411 }

◆ operator*=() [1/2]

const FGQuaternion & operator*= ( const FGQuaternion q)
inline

Arithmetic operator "*=".

Multiplication of two quaternions is like performing successive rotations.

Parameters
qa quaternion to be multiplied.
Returns
a quaternion reference representing Q, where Q = Q * q.

Definition at line 417 of file FGQuaternion.h.

417 {
418 double q0 = data[0]*q.data[0]-data[1]*q.data[1]-data[2]*q.data[2]-data[3]*q.data[3];
419 double q1 = data[0]*q.data[1]+data[1]*q.data[0]+data[2]*q.data[3]-data[3]*q.data[2];
420 double q2 = data[0]*q.data[2]-data[1]*q.data[3]+data[2]*q.data[0]+data[3]*q.data[1];
421 double q3 = data[0]*q.data[3]+data[1]*q.data[2]-data[2]*q.data[1]+data[3]*q.data[0];
422 data[0] = q0;
423 data[1] = q1;
424 data[2] = q2;
425 data[3] = q3;
426 mCacheValid = false;
427 return *this;
428 }

◆ operator*=() [2/2]

const FGQuaternion & operator*= ( double  scalar)
inline

Arithmetic operator "*=".

Parameters
scalara multiplicative value.
Returns
a quaternion reference representing Q, where Q = Q * scalar.

Definition at line 370 of file FGQuaternion.h.

370 {
371 data[0] *= scalar;
372 data[1] *= scalar;
373 data[2] *= scalar;
374 data[3] *= scalar;
375 mCacheValid = false;
376 return *this;
377 }

◆ operator+()

FGQuaternion operator+ ( const FGQuaternion q) const
inline

Arithmetic operator "+".

Parameters
qa quaternion to be summed.
Returns
a quaternion representing Q, where Q = Q + q.

Definition at line 389 of file FGQuaternion.h.

389 {
390 return FGQuaternion(data[0]+q.data[0], data[1]+q.data[1],
391 data[2]+q.data[2], data[3]+q.data[3]);
392 }

◆ operator+=()

const FGQuaternion & operator+= ( const FGQuaternion q)
inline

Definition at line 344 of file FGQuaternion.h.

344 {
345 // Copy the master values ...
346 data[0] += q.data[0];
347 data[1] += q.data[1];
348 data[2] += q.data[2];
349 data[3] += q.data[3];
350 mCacheValid = false;
351 return *this;
352 }

◆ operator-()

FGQuaternion operator- ( const FGQuaternion q) const
inline

Arithmetic operator "-".

Parameters
qa quaternion to be subtracted.
Returns
a quaternion representing Q, where Q = Q - q.

Definition at line 397 of file FGQuaternion.h.

397 {
398 return FGQuaternion(data[0]-q.data[0], data[1]-q.data[1],
399 data[2]-q.data[2], data[3]-q.data[3]);
400 }

◆ operator-=()

const FGQuaternion & operator-= ( const FGQuaternion q)
inline

Arithmetic operator "-=".

Parameters
qa quaternion reference.
Returns
a quaternion reference representing Q, where Q = Q - q.

Definition at line 357 of file FGQuaternion.h.

357 {
358 // Copy the master values ...
359 data[0] -= q.data[0];
360 data[1] -= q.data[1];
361 data[2] -= q.data[2];
362 data[3] -= q.data[3];
363 mCacheValid = false;
364 return *this;
365 }

◆ operator/=()

const FGQuaternion & operator/= ( double  scalar)
inline

Arithmetic operator "/=".

Parameters
scalara divisor value.
Returns
a quaternion reference representing Q, where Q = Q / scalar.

Definition at line 382 of file FGQuaternion.h.

382 {
383 return operator*=(1.0/scalar);
384 }
const FGQuaternion & operator*=(double scalar)
Arithmetic operator "*=".

◆ operator=()

const FGQuaternion & operator= ( const FGQuaternion q)
inline

Assignment operator "=".

Assign the value of q to the current object. Cached values are conserved.

Parameters
qreference to an FGQuaternion instance
Returns
reference to a quaternion object

Definition at line 310 of file FGQuaternion.h.

310 {
311 // Copy the master values ...
312 data[0] = q.data[0];
313 data[1] = q.data[1];
314 data[2] = q.data[2];
315 data[3] = q.data[3];
316 ComputeDerived();
317 // .. and copy the derived values if they are valid
318 mCacheValid = q.mCacheValid;
319 if (mCacheValid) {
320 mT = q.mT;
321 mTInv = q.mTInv;
322 mEulerAngles = q.mEulerAngles;
323 mEulerSines = q.mEulerSines;
324 mEulerCosines = q.mEulerCosines;
325 }
326 return *this;
327 }

◆ operator==()

bool operator== ( const FGQuaternion q) const
inline

Comparison operator "==".

Parameters
qa quaternion reference
Returns
true if both quaternions represent the same rotation.

Definition at line 335 of file FGQuaternion.h.

335 {
336 return data[0] == q.data[0] && data[1] == q.data[1]
337 && data[2] == q.data[2] && data[3] == q.data[3];
338 }

◆ SqrMagnitude()

double SqrMagnitude ( void  ) const
inline

Square of the length of the vector.

Compute and return the square of the euclidean norm of this vector.

Definition at line 466 of file FGQuaternion.h.

466 {
467 return data[0]*data[0] + data[1]*data[1]
468 + data[2]*data[2] + data[3]*data[3];
469 }

◆ zero()

static FGQuaternion zero ( void  )
inlinestatic

Zero quaternion vector.

Does not represent any orientation. Useful for initialization of increments

Definition at line 480 of file FGQuaternion.h.

480{ return FGQuaternion( 0.0, 0.0, 0.0, 0.0 ); }

Friends And Related Symbol Documentation

◆ operator*

FGQuaternion operator* ( double  scalar,
const FGQuaternion q 
)
friend

Scalar multiplication.

Parameters
scalarscalar value to multiply with.
qVector to multiply.

Multiply the Vector with a scalar value.

Definition at line 539 of file FGQuaternion.h.

539 {
540 return FGQuaternion(scalar*q.data[0], scalar*q.data[1], scalar*q.data[2], scalar*q.data[3]);
541}

◆ QExp

FGQuaternion QExp ( const FGColumnVector3 omega)
friend

Quaternion exponential.

Parameters
omegarotation velocity Calculate the unit quaternion which is the result of the exponentiation of the vector 'omega'.

Definition at line 548 of file FGQuaternion.h.

548 {
549 FGQuaternion qexp;
550 double angle = omega.Magnitude();
551 double sina_a = angle > 0.0 ? sin(angle)/angle : 1.0;
552
553 qexp.data[0] = cos(angle);
554 qexp.data[1] = omega(1) * sina_a;
555 qexp.data[2] = omega(2) * sina_a;
556 qexp.data[3] = omega(3) * sina_a;
557
558 return qexp;
559}

The documentation for this class was generated from the following files: