Quaternion


1   package JSci.maths;
2   
3   import JSci.GlobalSettings;
4   import JSci.maths.matrices.AbstractDoubleSquareMatrix;
5   import JSci.maths.matrices.DoubleSquareMatrix;
6   import JSci.maths.matrices.MatrixDimensionException;
7   import JSci.maths.vectors.Double3Vector;
8   import JSci.maths.groups.AbelianGroup;
9   import JSci.maths.fields.*;
10  import JSci.maths.algebras.*;
11  
12  /**
13  * The Quaternion class encapsulates quaternions.
14  * @jsci.planetmath Quaternions
15  * @version 1.1
16  * @author Mark Hale
17  *
18  * Quaternions are a noncommutative C<sup>*</sup>-algebra over the reals.
19  */
20  public final class Quaternion implements Field.Member, CStarAlgebra.Member {
21          private static final long serialVersionUID = 1605315490425547301L;
22  
23          private double re;
24          private double imi, imj, imk;
25  
26          public static final Quaternion ONE=new Quaternion(1.0, 0.0, 0.0, 0.0);
27          public static final Quaternion I=new Quaternion(0.0, 1.0, 0.0, 0.0);
28          public static final Quaternion J=new Quaternion(0.0, 0.0, 1.0, 0.0);
29          public static final Quaternion K=new Quaternion(0.0, 0.0, 0.0, 1.0);
30          /**
31          * Constructs a quaternion.
32          */
33          public Quaternion(final double real,final Double3Vector imag) {
34                  re=real;
35                  imi=imag.getComponent(0);
36                  imj=imag.getComponent(1);
37                  imk=imag.getComponent(2);
38          }
39          /**
40          * Constructs the quaternion q<sub>0</sub>+iq<sub>1</sub>+jq<sub>2</sub>+kq<sub>3</sub>.
41          */
42          public Quaternion(final double q0,final double q1,final double q2,final double q3) {
43                  re=q0;
44                  imi=q1;
45                  imj=q2;
46                  imk=q3;
47          }
48          /**
49          * Constructs a quaternion representing a rotation matrix.
50          * Note: if the matrix is not orthogonal then the quaternion will not have unit norm.
51          * Note: unit quaternions are a double cover of SO(3).
52          * @param m a rotation matrix
53          * @author Steve Lamont
54          */
55          public static Quaternion rotation(AbstractDoubleSquareMatrix m) {
56                  if(m.rows() != 3 && m.columns() != 3)
57                          throw new MatrixDimensionException("The matrix is not 3-dimensional.");
58                  double re, imi, imj, imk;
59                  double wSqr = ( 1.0 + m.trace() ) / 4.0;
60                  if(wSqr > GlobalSettings.ZERO_TOL) {
61                          re = Math.sqrt( wSqr );
62                          imi = ( m.getElement(2,1) - m.getElement(1,2) ) / ( re * 4.0 );
63                          imj = ( m.getElement(0,2) - m.getElement(2,0) ) / ( re * 4.0 );
64                          imk = ( m.getElement(1,0) - m.getElement(0,1) ) / ( re * 4.0 );
65              } else {
66                          double xSqr = -( m.getElement(1,1) + m.getElement(2,2) ) / 2.0;
67                          re = 0.0;
68                          if(xSqr > GlobalSettings.ZERO_TOL) {
69                          imi = Math.sqrt( xSqr );
70                          imj = m.getElement(1,0) / ( 2.0 * imi );
71                          imk = m.getElement(2,0) / ( 2.0 * imi );
72                          } else {
73                          double ySqr = ( 1.0 - m.getElement(2,2) ) / 2.0;
74                          imi = 0.0;
75                                  if(ySqr > GlobalSettings.ZERO_TOL) {
76                                          imj = Math.sqrt( ySqr );
77                                          imk = m.getElement(2,1) / ( 2.0 * imj );
78                          } else {
79                                          imj = 0.0;
80                                          imk = 1.0;
81                          }
82                          }
83              }
84                  return new Quaternion(re, imi, imj, imk);
85          }
86          /**
87          * Returns a 3D rotation matrix representing this quaternion.
88          * Note: if this quaternion does not have unit norm then the matrix will not be orthogonal.
89          * @author Steve Lamont
90          */
91          public AbstractDoubleSquareMatrix toRotationMatrix() {
92                  final double[][] array = new double[3][3];
93  
94                  array[0][0] = 1.0 - 2.0 * ( imj * imj + imk * imk );
95                  array[0][1] = 2.0 * ( imi * imj - re * imk );
96              array[0][2] = 2.0 * ( imi * imk + re * imj );
97  
98              array[1][0] = 2.0 * ( imi * imj + re * imk );
99              array[1][1] = 1.0 - 2.0 * ( imi * imi + imk * imk );
100             array[1][2] = 2.0 * ( imj * imk - re * imi );
101 
102                 array[2][0] = 2.0 * ( imi * imk - re * imj );
103                 array[2][1] = 2.0 * ( imj * imk + re * imi );
104                 array[2][2] = 1.0 - 2.0 * ( imi * imi + imj * imj );
105 
106                 return new DoubleSquareMatrix( array );
107         }
108         /**
109         * Compares two quaternions for equality.
110         * @param obj a quaternion
111         */
112         public boolean equals(Object   obj) {
113                 if(obj instanceof Quaternion) {
114                         final Quaternion q = (Quaternion)obj;
115                         return this.subtract(q).norm() <= GlobalSettings.ZERO_TOL;
116                 } else
117                         return false;
118         }
119         /**
120         * Returns a string representing the value of this quaternion.
121         */
122         public String   toString() {
123                 final StringBuffer   buf=new StringBuffer  (40);
124                 buf.append(re);
125                 if(imi>=0.0)
126                         buf.append("+");
127                 buf.append(imi);
128                 buf.append("i");
129                 if(imj>=0.0)
130                         buf.append("+");
131                 buf.append(imj);
132                 buf.append("j");
133                 if(imk>=0.0)
134                         buf.append("+");
135                 buf.append(imk);
136                 buf.append("k");
137                 return buf.toString();
138         }
139         /**
140         * Returns a hashcode for this quaternion.
141         */
142         public int hashCode() {
143                 return (int)(Math.exp(norm()));
144         }
145         /**
146         * Returns true if either the real or imaginary part is NaN.
147         */
148         public boolean isNaN() {
149                 return (re==Double.NaN) || (imi==Double.NaN) ||
150                 (imj==Double.NaN) || (imk==Double.NaN);
151         }
152         /**
153         * Returns true if either the real or imaginary part is infinite.
154         */
155         public boolean isInfinite() {
156                 return (re==Double.POSITIVE_INFINITY) || (re==Double.NEGATIVE_INFINITY)
157                         || (imi==Double.POSITIVE_INFINITY) || (imi==Double.NEGATIVE_INFINITY)
158                         || (imj==Double.POSITIVE_INFINITY) || (imj==Double.NEGATIVE_INFINITY)
159                         || (imk==Double.POSITIVE_INFINITY) || (imk==Double.NEGATIVE_INFINITY);
160         }
161         /**
162         * Returns the real part of this quaternion.
163         */
164         public double real() {
165                 return re;
166         }
167         /**
168         * Returns the imaginary part of this quaternion.
169         */
170         public Double3Vector imag() {
171                 return new Double3Vector(imi, imj, imk);
172         }
173     public Object   getSet() {
174         throw new RuntimeException  ("Not yet implemented: please file bug report");
175     }
176         /**
177         * Returns the l<sup>2</sup>-norm (magnitude),
178         * which is also the C<sup>*</sup> norm.
179         */
180         public double norm() {
181                 return Math.sqrt(sumSquares());
182         }
183         /**
184         * Returns the sum of the squares of the components.
185         */
186         public double sumSquares() {
187                 return re*re+imi*imi+imj*imj+imk*imk;
188         }
189 
190 //============
191 // OPERATIONS
192 //============
193 
194         /**
195         * Returns the negative of this quaternion.
196         */
197         public AbelianGroup.Member negate() {
198                 return new Quaternion(-re, -imi, -imj, -imk);
199         }
200         /**
201         * Returns the inverse of this quaternion.
202         */
203         public Field.Member inverse() {
204                 final double sumSqr=sumSquares();
205                 return new Quaternion(re/sumSqr, -imi/sumSqr, -imj/sumSqr, -imk/sumSqr);
206         }
207         /**
208         * Returns the involution of this quaternion.
209         */
210         public CStarAlgebra.Member involution() {
211                 return conjugate();
212         }
213         /**
214         * Returns the conjugate of this quaternion.
215         */
216         public Quaternion conjugate() {
217                 return new Quaternion(re, -imi, -imj, -imk);
218         }
219 
220 // ADDITION
221 
222         /**
223         * Returns the addition of this number and another.
224         */
225         public AbelianGroup.Member add(final AbelianGroup.Member x) {
226                 if(x instanceof Quaternion)
227                         return add((Quaternion)x);
228                 else if(x instanceof MathDouble)
229                         return addReal(((MathDouble)x).value());
230                 else if(x instanceof MathInteger)
231                         return addReal(((MathInteger)x).value());
232                 else
233                         throw new IllegalArgumentException  ("Member class not recognised by this method.");
234         }
235         /**
236         * Returns the addition of this quaternion and another.
237         * @param q a quaternion
238         */
239         public Quaternion add(final Quaternion q) {
240                 return new Quaternion(re+q.re, imi+q.imi, imj+q.imj, imk+q.imk);
241         }
242         /**
243         * Returns the addition of this quaternion with a real part.
244         * @param real a real part
245         */
246         public Quaternion addReal(final double real) {
247                 return new Quaternion(re+real, imi, imj, imk);
248         }
249         /**
250         * Returns the addition of this quaternion with an imaginary part.
251         * @param imag an imaginary part
252         */
253         public Quaternion addImag(final Double3Vector imag) {
254                 return new Quaternion(re, imi+imag.getComponent(0), imj+imag.getComponent(1), imk+imag.getComponent(2));
255         }
256 
257 // SUBTRACTION
258 
259         /**
260         * Returns the subtraction of this number and another.
261         */
262         public AbelianGroup.Member subtract(final AbelianGroup.Member x) {
263                 if(x instanceof Quaternion)
264                         return subtract((Quaternion)x);
265                 else if(x instanceof MathDouble)
266                         return subtractReal(((MathDouble)x).value());
267                 else if(x instanceof MathInteger)
268                         return subtractReal(((MathInteger)x).value());
269                 else
270                         throw new IllegalArgumentException  ("Member class not recognised by this method.");
271         }
272         /**
273         * Returns the subtraction of this quaternion by another.
274         * @param q a quaternion
275         */
276         public Quaternion subtract(final Quaternion q) {
277                 return new Quaternion(re-q.re, imi-q.imi, imj-q.imj, imk-q.imk);
278         }
279         /**
280         * Returns the subtraction of this quaternion by a real part.
281         * @param real a real part
282         */
283         public Quaternion subtractReal(final double real) {
284                 return new Quaternion(re-real, imi, imj, imk);
285         }
286         /**
287         * Returns the subtraction of this quaternion by an imaginary part.
288         * @param imag an imaginary part
289         */
290         public Quaternion subtractImag(final Double3Vector imag) {
291                 return new Quaternion(re, imi-imag.getComponent(0), imj-imag.getComponent(1), imk-imag.getComponent(2));
292         }
293 
294 // MULTIPLICATION
295 
296         /**
297         * Returns the multiplication of this number by a real scalar.
298         */
299         public Module.Member scalarMultiply(final Ring.Member x) {
300                 if(x instanceof MathDouble)
301                         return multiply(((MathDouble)x).value());
302                 else if(x instanceof MathInteger)
303                         return multiply(((MathInteger)x).value());
304                 else
305                         throw new IllegalArgumentException  ("Member class not recognised by this method.");
306         }
307         /**
308         * Returns the multiplication of this number and another.
309         */
310         public Ring.Member multiply(final Ring.Member x) {
311                 if(x instanceof Quaternion)
312                         return multiply((Quaternion)x);
313                 else if(x instanceof MathDouble)
314                         return multiply(((MathDouble)x).value());
315                 else if(x instanceof MathInteger)
316                         return multiply(((MathInteger)x).value());
317                 else
318                         throw new IllegalArgumentException  ("Member class not recognised by this method.");
319         }
320         /**
321         * Returns the multiplication of this quaternion and another.
322         * @param q a quaternion
323         */
324         public Quaternion multiply(final Quaternion q) {
325                 return new Quaternion(
326                         re*q.re-imi*q.imi-imj*q.imj-imk*q.imk,
327                         re*q.imi+q.re*imi+(imj*q.imk-q.imj*imk),
328                         re*q.imj+q.re*imj+(imk*q.imi-q.imk*imi),
329                         re*q.imk+q.re*imk+(imi*q.imj-q.imi*imj)
330                 );
331         }
332         /**
333         * Returns the multiplication of this quaternion by a scalar.
334         * @param x a real number
335         */
336         public Quaternion multiply(final double x) {
337                 return new Quaternion(x*re, x*imi, x*imj, x*imk);
338         }
339 
340 // DIVISION
341 
342         /**
343         * Returns the division of this number by a real scalar.
344         */
345         public VectorSpace.Member scalarDivide(final Field.Member x) {
346                 if(x instanceof MathDouble)
347                         return divide(((MathDouble)x).value());
348                 else
349                         throw new IllegalArgumentException  ("Member class not recognised by this method.");
350         }
351         /**
352         * Returns the division of this number and another.
353         */
354         public Field.Member divide(final Field.Member x) {
355                 if(x instanceof Quaternion)
356                         return divide((Quaternion)x);
357                 else if(x instanceof MathDouble)
358                         return divide(((MathDouble)x).value());
359                 else
360                         throw new IllegalArgumentException  ("Member class not recognised by this method.");
361         }
362         /**
363         * Returns the division of this quaternion by another.
364         * @param q a quaternion
365         * @exception ArithmeticException If divide by zero.
366         */
367         public Quaternion divide(final Quaternion q) {
368                 final double qSumSqr=q.sumSquares();
369                 return new Quaternion(
370                         (re*q.re+imi*q.imi+imj*q.imj+imk*q.imk)/qSumSqr,
371                         (q.re*imi-re*q.imi-(imj*q.imk-q.imj*imk))/qSumSqr,
372                         (q.re*imj-re*q.imj-(imk*q.imi-q.imk*imi))/qSumSqr,
373                         (q.re*imk-re*q.imk-(imi*q.imj-q.imi*imj))/qSumSqr
374                 );
375         }
376         /**
377         * Returns the division of this quaternion by a scalar.
378         * @param x a real number
379         * @exception ArithmeticException If divide by zero.
380         */
381         public Quaternion divide(final double x) {
382                 return new Quaternion(re/x, imi/x, imj/x, imk/x);
383         }
384         public Quaternion normalize() {
385                 return this.divide(norm());
386         }
387 
388 //============
389 // FUNCTIONS
390 //============
391 
392         public static Quaternion exp(Quaternion q) {
393                 final double k = Math.exp(q.re);
394                 Double3Vector imag = q.imag();
395                 final double imagNorm = imag.norm();
396                 return new Quaternion(k*Math.cos(imagNorm), (Double3Vector) imag.normalize().scalarMultiply(k*Math.sin(imagNorm)));
397         }
398         public static Quaternion log(Quaternion q) {
399                 final double norm = q.norm();
400                 return new Quaternion(Math.log(norm), (Double3Vector) q.imag().normalize().scalarMultiply(Math.acos(q.re/norm)));
401         }
402         public static Quaternion sin(Quaternion q) {
403                 Double3Vector imag = q.imag();
404                 final double imagNorm = imag.norm();
405                 return new Quaternion(Math.sin(q.re)*ExtraMath.cosh(imagNorm), (Double3Vector) imag.normalize().scalarMultiply(Math.cos(q.re)*ExtraMath.sinh(imagNorm)));
406         }
407         public static Quaternion cos(Quaternion q) {
408                 Double3Vector imag = q.imag();
409                 final double imagNorm = imag.norm();
410                 return new Quaternion(Math.cos(q.re)*ExtraMath.cosh(imagNorm), (Double3Vector) imag.normalize().scalarMultiply(-Math.sin(q.re)*ExtraMath.sinh(imagNorm)));
411         }
412         public static Quaternion tan(Quaternion q) {
413                 return sin(q).divide(cos(q));
414         }
415         public static Quaternion sinh(Quaternion q) {
416                 Double3Vector imag = q.imag();
417                 final double imagNorm = imag.norm();
418                 return new Quaternion(ExtraMath.sinh(q.re)*Math.cos(imagNorm), (Double3Vector) imag.normalize().scalarMultiply(ExtraMath.cosh(q.re)*Math.sin(imagNorm)));
419         }
420         public static Quaternion cosh(Quaternion q) {
421                 Double3Vector imag = q.imag();
422                 final double imagNorm = imag.norm();
423                 return new Quaternion(ExtraMath.cosh(q.re)*Math.cos(imagNorm), (Double3Vector) imag.normalize().scalarMultiply(ExtraMath.sinh(q.re)*Math.sin(imagNorm)));
424         }
425         public static Quaternion tanh(Quaternion q) {
426                 return sinh(q).divide(cosh(q));
427         }
428 }
429 
430
A to Z: JavaDoc & Examples Daily Java News & Articles Open Source Projects Open Source Codes Free Computer Books Remove Frame
Popular Tags