ROOT  6.07/01
Reference Guide
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
3DConversions.cxx
Go to the documentation of this file.
1 // @(#)root/mathcore:$Id$
2 // Authors: W. Brown, M. Fischler, L. Moneta 2005
3 
4  /**********************************************************************
5  * *
6  * Copyright (c) 2005, LCG ROOT FNAL MathLib Team *
7  * *
8  * *
9  **********************************************************************/
10 
11 // Source file for something else
12 //
13 // Created by: Mark Fischler and Walter Brown Thurs July 7, 2005
14 //
15 // Last update: $Id$
16 //
17 
18 // TODO - For now, all conversions are grouped in this one compilation unit.
19 // The intention is to seraparte them into a few .cpp files instead,
20 // so that users needing one form need not incorporate code for them all.
21 
22 #include <cmath>
23 #include <limits>
24 
25 
26 #include "Math/Math.h"
27 
29 
38 
39 
40 namespace ROOT {
41 namespace Math {
42 namespace gv_detail {
43 
48 };
49 
50 
51 // ----------------------------------------------------------------------
52 void convert( Rotation3D const & from, AxisAngle & to)
53 {
54  // conversions from Rotation3D
55  double m[9];
56  from.GetComponents(m, m+9);
57 
58  const double uZ = m[kYX] - m[kXY];
59  const double uY = m[kXZ] - m[kZX];
60  const double uX = m[kZY] - m[kYZ];
61 
62 
63  // in case of rotaiton of an angle PI, the rotation matrix is symmetric and
64  // uX = uY = uZ = 0. Use then conversion through the quaternion
68  Quaternion tmp;
69  convert (from,tmp);
70  convert (tmp,to);
71  return;
72  }
73 
75 
76  u.SetCoordinates( uX, uY, uZ );
77 
78  static const double pi = M_PI;
79 
80  double angle;
81  const double cosdelta = (m[kXX] + m[kYY] + m[kZZ] - 1.0) / 2.0;
82  if (cosdelta > 1.0) {
83  angle = 0;
84  } else if (cosdelta < -1.0) {
85  angle = pi;
86  } else {
87  angle = std::acos( cosdelta );
88  }
89 
90 
91  //to.SetAngle(angle);
92  to.SetComponents(u, angle);
93  to.Rectify();
94 
95 } // convert to AxisAngle
96 
97 static void correctByPi ( double& psi, double& phi ) {
98  static const double pi = M_PI;
99  if (psi > 0) {
100  psi -= pi;
101  } else {
102  psi += pi;
103  }
104  if (phi > 0) {
105  phi -= pi;
106  } else {
107  phi += pi;
108  }
109 }
110 
111 void convert( Rotation3D const & from, EulerAngles & to)
112 {
113  // conversion from Rotation3D to Euler Angles
114  // Mathematical justification appears in
115  // http://www.cern.ch/mathlibs/documents/eulerAngleComputation.pdf
116 
117  double r[9];
118  from.GetComponents(r,r+9);
119 
120  double phi, theta, psi;
121  double psiPlusPhi, psiMinusPhi;
122  static const double pi = M_PI;
123  static const double pi_2 = M_PI_2;
124 
125  theta = (std::fabs(r[kZZ]) <= 1.0) ? std::acos(r[kZZ]) :
126  (r[kZZ] > 0.0) ? 0 : pi;
127 
128  double cosTheta = r[kZZ];
129  if (cosTheta > 1) cosTheta = 1;
130  if (cosTheta < -1) cosTheta = -1;
131 
132  // Compute psi +/- phi:
133  // Depending on whether cosTheta is positive or negative and whether it
134  // is less than 1 in absolute value, different mathematically equivalent
135  // expressions are numerically stable.
136  if (cosTheta == 1) {
137  psiPlusPhi = atan2 ( r[kXY] - r[kYX], r[kXX] + r[kYY] );
138  psiMinusPhi = 0;
139  } else if (cosTheta >= 0) {
140  psiPlusPhi = atan2 ( r[kXY] - r[kYX], r[kXX] + r[kYY] );
141  double s = -r[kXY] - r[kYX]; // sin (psi-phi) * (1 - cos theta)
142  double c = r[kXX] - r[kYY]; // cos (psi-phi) * (1 - cos theta)
143  psiMinusPhi = atan2 ( s, c );
144  } else if (cosTheta > -1) {
145  psiMinusPhi = atan2 ( -r[kXY] - r[kYX], r[kXX] - r[kYY] );
146  double s = r[kXY] - r[kYX]; // sin (psi+phi) * (1 + cos theta)
147  double c = r[kXX] + r[kYY]; // cos (psi+phi) * (1 + cos theta)
148  psiPlusPhi = atan2 ( s, c );
149  } else { // cosTheta == -1
150  psiMinusPhi = atan2 ( -r[kXY] - r[kYX], r[kXX] - r[kYY] );
151  psiPlusPhi = 0;
152  }
153 
154  psi = .5 * (psiPlusPhi + psiMinusPhi);
155  phi = .5 * (psiPlusPhi - psiMinusPhi);
156 
157  // Now correct by pi if we have managed to get a value of psiPlusPhi
158  // or psiMinusPhi that was off by 2 pi:
159 
160  // set up w[i], all of which would be positive if sin and cosine of
161  // psi and phi were positive:
162  double w[4];
163  w[0] = r[kXZ]; w[1] = r[kZX]; w[2] = r[kYZ]; w[3] = -r[kZY];
164 
165  // find biggest relevant term, which is the best one to use in correcting.
166  double maxw = std::fabs(w[0]);
167  int imax = 0;
168  for (int i = 1; i < 4; ++i) {
169  if (std::fabs(w[i]) > maxw) {
170  maxw = std::fabs(w[i]);
171  imax = i;
172  }
173  }
174  // Determine if the correction needs to be applied: The criteria are
175  // different depending on whether a sine or cosine was the determinor:
176  switch (imax) {
177  case 0:
178  if (w[0] > 0 && psi < 0) correctByPi ( psi, phi );
179  if (w[0] < 0 && psi > 0) correctByPi ( psi, phi );
180  break;
181  case 1:
182  if (w[1] > 0 && phi < 0) correctByPi ( psi, phi );
183  if (w[1] < 0 && phi > 0) correctByPi ( psi, phi );
184  break;
185  case 2:
186  if (w[2] > 0 && std::fabs(psi) > pi_2) correctByPi ( psi, phi );
187  if (w[2] < 0 && std::fabs(psi) < pi_2) correctByPi ( psi, phi );
188  break;
189  case 3:
190  if (w[3] > 0 && std::fabs(phi) > pi_2) correctByPi ( psi, phi );
191  if (w[3] < 0 && std::fabs(phi) < pi_2) correctByPi ( psi, phi );
192  break;
193  }
194 
195  to.SetComponents( phi, theta, psi );
196 
197 } // convert to EulerAngles
198 
199 ////////////////////////////////////////////////////////////////////////////////
200 /// conversion from Rotation3D to Quaternion
201 
202 void convert( Rotation3D const & from, Quaternion & to)
203 {
204  double m[9];
205  from.GetComponents(m, m+9);
206 
207  const double d0 = m[kXX] + m[kYY] + m[kZZ];
208  const double d1 = + m[kXX] - m[kYY] - m[kZZ];
209  const double d2 = - m[kXX] + m[kYY] - m[kZZ];
210  const double d3 = - m[kXX] - m[kYY] + m[kZZ];
211 
212  // these are related to the various q^2 values;
213  // choose the largest to avoid dividing two small numbers and losing accuracy.
214 
215  if ( d0 >= d1 && d0 >= d2 && d0 >= d3 ) {
216  const double q0 = .5*std::sqrt(1+d0);
217  const double f = .25/q0;
218  const double q1 = f*(m[kZY]-m[kYZ]);
219  const double q2 = f*(m[kXZ]-m[kZX]);
220  const double q3 = f*(m[kYX]-m[kXY]);
221  to.SetComponents(q0,q1,q2,q3);
222  to.Rectify();
223  return;
224  } else if ( d1 >= d2 && d1 >= d3 ) {
225  const double q1 = .5*std::sqrt(1+d1);
226  const double f = .25/q1;
227  const double q0 = f*(m[kZY]-m[kYZ]);
228  const double q2 = f*(m[kXY]+m[kYX]);
229  const double q3 = f*(m[kXZ]+m[kZX]);
230  to.SetComponents(q0,q1,q2,q3);
231  to.Rectify();
232  return;
233  } else if ( d2 >= d3 ) {
234  const double q2 = .5*std::sqrt(1+d2);
235  const double f = .25/q2;
236  const double q0 = f*(m[kXZ]-m[kZX]);
237  const double q1 = f*(m[kXY]+m[kYX]);
238  const double q3 = f*(m[kYZ]+m[kZY]);
239  to.SetComponents(q0,q1,q2,q3);
240  to.Rectify();
241  return;
242  } else {
243  const double q3 = .5*std::sqrt(1+d3);
244  const double f = .25/q3;
245  const double q0 = f*(m[kYX]-m[kXY]);
246  const double q1 = f*(m[kXZ]+m[kZX]);
247  const double q2 = f*(m[kYZ]+m[kZY]);
248  to.SetComponents(q0,q1,q2,q3);
249  to.Rectify();
250  return;
251  }
252 } // convert to Quaternion
253 
254 ////////////////////////////////////////////////////////////////////////////////
255 /// conversion from Rotation3D to RotationZYX
256 /// same Math used as for EulerAngles apart from some different meaning of angles and
257 /// matrix elements. But the basic algoprithms principles are the same described in
258 /// http://www.cern.ch/mathlibs/documents/eulerAngleComputation.pdf
259 
260 void convert( Rotation3D const & from, RotationZYX & to)
261 {
262  // theta is assumed to be in range [-PI/2,PI/2].
263  // this is guranteed by the Rectify function
264 
265  static const double pi_2 = M_PI_2;
266 
267  double r[9];
268  from.GetComponents(r,r+9);
269 
270  double phi,theta,psi = 0;
271 
272  // careful for numeical error make sin(theta) ourtside [-1,1]
273  double sinTheta = r[kXZ];
274  if ( sinTheta < -1.0) sinTheta = -1.0;
275  if ( sinTheta > 1.0) sinTheta = 1.0;
276  theta = std::asin( sinTheta );
277 
278  // compute psi +/- phi
279  // Depending on whether cosTheta is positive or negative and whether it
280  // is less than 1 in absolute value, different mathematically equivalent
281  // expressions are numerically stable.
282  // algorithm from
283  // adapted for the case 3-2-1
284 
285  double psiPlusPhi = 0;
286  double psiMinusPhi = 0;
287 
288  // valid if sinTheta not eq to -1 otherwise is zero
289  if (sinTheta > - 1.0)
290  psiPlusPhi = atan2 ( r[kYX] + r[kZY], r[kYY] - r[kZX] );
291 
292  // valid if sinTheta not eq. to 1
293  if (sinTheta < 1.0)
294  psiMinusPhi = atan2 ( r[kZY] - r[kYX] , r[kYY] + r[kZX] );
295 
296  psi = .5 * (psiPlusPhi + psiMinusPhi);
297  phi = .5 * (psiPlusPhi - psiMinusPhi);
298 
299  // correction is not necessary if sinTheta = +/- 1
300  //if (sinTheta == 1.0 || sinTheta == -1.0) return;
301 
302  // apply the corrections according to max of the other terms
303  // I think is assumed convention that theta is between -PI/2,PI/2.
304  // OTHERWISE RESULT MIGHT BE DIFFERENT ???
305 
306  //since we determine phi+psi or phi-psi phi and psi can be both have a shift of +/- PI.
307  // The shift must be applied on both (the sum (or difference) is knows to +/- 2PI )
308  //This can be fixed looking at the other 4 matrix terms, which have terms in sin and cos of psi
309  // and phi. sin(psi+/-PI) = -sin(psi) and cos(psi+/-PI) = -cos(psi).
310  //Use then the biggest term for making the correction to minimize possible numerical errors
311 
312  // set up w[i], all of which would be positive if sin and cosine of
313  // psi and phi were positive:
314  double w[4];
315  w[0] = -r[kYZ]; w[1] = -r[kXY]; w[2] = r[kZZ]; w[3] = r[kXX];
316 
317  // find biggest relevant term, which is the best one to use in correcting.
318  double maxw = std::fabs(w[0]);
319  int imax = 0;
320  for (int i = 1; i < 4; ++i) {
321  if (std::fabs(w[i]) > maxw) {
322  maxw = std::fabs(w[i]);
323  imax = i;
324  }
325  }
326 
327  // Determine if the correction needs to be applied: The criteria are
328  // different depending on whether a sine or cosine was the determinor:
329  switch (imax) {
330  case 0:
331  if (w[0] > 0 && psi < 0) correctByPi ( psi, phi );
332  if (w[0] < 0 && psi > 0) correctByPi ( psi, phi );
333  break;
334  case 1:
335  if (w[1] > 0 && phi < 0) correctByPi ( psi, phi );
336  if (w[1] < 0 && phi > 0) correctByPi ( psi, phi );
337  break;
338  case 2:
339  if (w[2] > 0 && std::fabs(psi) > pi_2) correctByPi ( psi, phi );
340  if (w[2] < 0 && std::fabs(psi) < pi_2) correctByPi ( psi, phi );
341  break;
342  case 3:
343  if (w[3] > 0 && std::fabs(phi) > pi_2) correctByPi ( psi, phi );
344  if (w[3] < 0 && std::fabs(phi) < pi_2) correctByPi ( psi, phi );
345  break;
346  }
347 
348  to.SetComponents(phi, theta, psi);
349 
350 } // convert to RotationZYX
351 
352 // ----------------------------------------------------------------------
353 // conversions from AxisAngle
354 
355 void convert( AxisAngle const & from, Rotation3D & to)
356 {
357  // conversion from AxixAngle to Rotation3D
358 
359  const double sinDelta = std::sin( from.Angle() );
360  const double cosDelta = std::cos( from.Angle() );
361  const double oneMinusCosDelta = 1.0 - cosDelta;
362 
363  const AxisAngle::AxisVector & u = from.Axis();
364  const double uX = u.X();
365  const double uY = u.Y();
366  const double uZ = u.Z();
367 
368  double m[9];
369 
370  m[kXX] = oneMinusCosDelta * uX * uX + cosDelta;
371  m[kXY] = oneMinusCosDelta * uX * uY - sinDelta * uZ;
372  m[kXZ] = oneMinusCosDelta * uX * uZ + sinDelta * uY;
373 
374  m[kYX] = oneMinusCosDelta * uY * uX + sinDelta * uZ;
375  m[kYY] = oneMinusCosDelta * uY * uY + cosDelta;
376  m[kYZ] = oneMinusCosDelta * uY * uZ - sinDelta * uX;
377 
378  m[kZX] = oneMinusCosDelta * uZ * uX - sinDelta * uY;
379  m[kZY] = oneMinusCosDelta * uZ * uY + sinDelta * uX;
380  m[kZZ] = oneMinusCosDelta * uZ * uZ + cosDelta;
381 
382  to.SetComponents(m,m+9);
383 } // convert to Rotation3D
384 
385 void convert( AxisAngle const & from , EulerAngles & to )
386 {
387  // conversion from AxixAngle to EulerAngles
388  // TODO better : temporary make conversion using Rotation3D
389 
390  Rotation3D tmp;
391  convert(from,tmp);
392  convert(tmp,to);
393 }
394 
395 void convert( AxisAngle const & from, Quaternion & to)
396 {
397  // conversion from AxixAngle to Quaternion
398 
399  double s = std::sin (from.Angle()/2);
401 
402  to.SetComponents( std::cos(from.Angle()/2),
403  s*axis.X(),
404  s*axis.Y(),
405  s*axis.Z()
406  );
407 } // convert to Quaternion
408 
409 void convert( AxisAngle const & from , RotationZYX & to )
410 {
411  // conversion from AxisAngle to RotationZYX
412  // TODO better : temporary make conversion using Rotation3D
413  Rotation3D tmp;
414  convert(from,tmp);
415  convert(tmp,to);
416 }
417 
418 
419 // ----------------------------------------------------------------------
420 // conversions from EulerAngles
421 
422 void convert( EulerAngles const & from, Rotation3D & to)
423 {
424  // conversion from EulerAngles to Rotation3D
425 
426  typedef double Scalar;
427  const Scalar sPhi = std::sin( from.Phi() );
428  const Scalar cPhi = std::cos( from.Phi() );
429  const Scalar sTheta = std::sin( from.Theta() );
430  const Scalar cTheta = std::cos( from.Theta() );
431  const Scalar sPsi = std::sin( from.Psi() );
432  const Scalar cPsi = std::cos( from.Psi() );
433  to.SetComponents
434  ( cPsi*cPhi-sPsi*cTheta*sPhi, cPsi*sPhi+sPsi*cTheta*cPhi, sPsi*sTheta
435  , -sPsi*cPhi-cPsi*cTheta*sPhi, -sPsi*sPhi+cPsi*cTheta*cPhi, cPsi*sTheta
436  , sTheta*sPhi, -sTheta*cPhi, cTheta
437  );
438 }
439 
440 void convert( EulerAngles const & from, AxisAngle & to)
441 {
442  // conversion from EulerAngles to AxisAngle
443  // make converting first to quaternion
444  Quaternion q;
445  convert (from, q);
446  convert (q, to);
447 }
448 
449 void convert( EulerAngles const & from, Quaternion & to)
450 {
451  // conversion from EulerAngles to Quaternion
452 
453  typedef double Scalar;
454  const Scalar plus = (from.Phi()+from.Psi())/2;
455  const Scalar minus = (from.Phi()-from.Psi())/2;
456  const Scalar sPlus = std::sin( plus );
457  const Scalar cPlus = std::cos( plus );
458  const Scalar sMinus = std::sin( minus );
459  const Scalar cMinus = std::cos( minus );
460  const Scalar sTheta = std::sin( from.Theta()/2 );
461  const Scalar cTheta = std::cos( from.Theta()/2 );
462 
463  to.SetComponents ( cTheta*cPlus, -sTheta*cMinus, -sTheta*sMinus, -cTheta*sPlus );
464  // TODO -- carefully check that this is correct
465 }
466 
467 void convert( EulerAngles const & from , RotationZYX & to )
468 {
469  // conversion from EulerAngles to RotationZYX
470  // TODO better : temporary make conversion using Rotation3D
471  Rotation3D tmp;
472  convert(from,tmp);
473  convert(tmp,to);
474 }
475 
476 
477 // ----------------------------------------------------------------------
478 // conversions from Quaternion
479 
480 void convert( Quaternion const & from, Rotation3D & to)
481 {
482  // conversion from Quaternion to Rotation3D
483 
484  const double q0 = from.U();
485  const double q1 = from.I();
486  const double q2 = from.J();
487  const double q3 = from.K();
488  const double q00 = q0*q0;
489  const double q01 = q0*q1;
490  const double q02 = q0*q2;
491  const double q03 = q0*q3;
492  const double q11 = q1*q1;
493  const double q12 = q1*q2;
494  const double q13 = q1*q3;
495  const double q22 = q2*q2;
496  const double q23 = q2*q3;
497  const double q33 = q3*q3;
498 
499  to.SetComponents (
500  q00+q11-q22-q33 , 2*(q12-q03) , 2*(q02+q13),
501  2*(q12+q03) , q00-q11+q22-q33 , 2*(q23-q01),
502  2*(q13-q02) , 2*(q23+q01) , q00-q11-q22+q33 );
503 
504 } // conversion to Rotation3D
505 
506 void convert( Quaternion const & from, AxisAngle & to)
507 {
508  // conversion from Quaternion to AxisAngle
509 
510  double u = from.U();
511  if ( u >= 0 ) {
512  if ( u > 1 ) u = 1;
513  const double angle = 2.0 * std::acos ( from.U() );
515  axis (from.I(), from.J(), from.K());
516  to.SetComponents ( axis, angle );
517  } else {
518  if ( u < -1 ) u = -1;
519  const double angle = 2.0 * std::acos ( -from.U() );
521  axis (-from.I(), -from.J(), -from.K());
522  to.SetComponents ( axis, angle );
523  }
524 } // conversion to AxisAngle
525 
526 void convert( Quaternion const & from, EulerAngles & to )
527 {
528  // conversion from Quaternion to EulerAngles
529  // TODO better
530  // temporary make conversion using Rotation3D
531 
532  Rotation3D tmp;
533  convert(from,tmp);
534  convert(tmp,to);
535 }
536 
537 void convert( Quaternion const & from , RotationZYX & to )
538 {
539  // conversion from Quaternion to RotationZYX
540  // TODO better : temporary make conversion using Rotation3D
541  Rotation3D tmp;
542  convert(from,tmp);
543  convert(tmp,to);
544 }
545 
546 // ----------------------------------------------------------------------
547 // conversions from RotationZYX
548 void convert( RotationZYX const & from, Rotation3D & to) {
549  // conversion to Rotation3D (matrix)
550 
551  double phi,theta,psi = 0;
552  from.GetComponents(phi,theta,psi);
553  to.SetComponents( std::cos(theta)*std::cos(phi),
554  - std::cos(theta)*std::sin(phi),
555  std::sin(theta),
556 
557  std::cos(psi)*std::sin(phi) + std::sin(psi)*std::sin(theta)*std::cos(phi),
558  std::cos(psi)*std::cos(phi) - std::sin(psi)*std::sin(theta)*std::sin(phi),
559  -std::sin(psi)*std::cos(theta),
560 
561  std::sin(psi)*std::sin(phi) - std::cos(psi)*std::sin(theta)*std::cos(phi),
562  std::sin(psi)*std::cos(phi) + std::cos(psi)*std::sin(theta)*std::sin(phi),
563  std::cos(psi)*std::cos(theta)
564  );
565 
566 }
567 void convert( RotationZYX const & from, AxisAngle & to) {
568  // conversion to axis angle
569  // TODO better : temporary make conversion using Rotation3D
570  Rotation3D tmp;
571  convert(from,tmp);
572  convert(tmp,to);
573 }
574 void convert( RotationZYX const & from, EulerAngles & to) {
575  // conversion to Euler angle
576  // TODO better : temporary make conversion using Rotation3D
577  Rotation3D tmp;
578  convert(from,tmp);
579  convert(tmp,to);
580 }
581 void convert( RotationZYX const & from, Quaternion & to) {
582  double phi,theta,psi = 0;
583  from.GetComponents(phi,theta,psi);
584 
585  double sphi2 = std::sin(phi/2);
586  double cphi2 = std::cos(phi/2);
587  double stheta2 = std::sin(theta/2);
588  double ctheta2 = std::cos(theta/2);
589  double spsi2 = std::sin(psi/2);
590  double cpsi2 = std::cos(psi/2);
591  to.SetComponents( cphi2 * cpsi2 * ctheta2 - sphi2 * spsi2 * stheta2,
592  sphi2 * cpsi2 * stheta2 + cphi2 * spsi2 * ctheta2,
593  cphi2 * cpsi2 * stheta2 - sphi2 * spsi2 * ctheta2,
594  sphi2 * cpsi2 * ctheta2 + cphi2 * spsi2 * stheta2
595  );
596 }
597 
598 
599 // ----------------------------------------------------------------------
600 // conversions from RotationX
601 
602 void convert( RotationX const & from, Rotation3D & to)
603 {
604  // conversion from RotationX to Rotation3D
605 
606  const double c = from.CosAngle();
607  const double s = from.SinAngle();
608  to.SetComponents ( 1, 0, 0,
609  0, c, -s,
610  0, s, c );
611 }
612 
613 void convert( RotationX const & from, AxisAngle & to)
614 {
615  // conversion from RotationX to AxisAngle
616 
618  to.SetComponents ( axis, from.Angle() );
619 }
620 
621 void convert( RotationX const & from , EulerAngles & to )
622 {
623  // conversion from RotationX to EulerAngles
624  //TODO better: temporary make conversion using Rotation3D
625 
626  Rotation3D tmp;
627  convert(from,tmp);
628  convert(tmp,to);
629 }
630 
631 void convert( RotationX const & from, Quaternion & to)
632 {
633  // conversion from RotationX to Quaternion
634 
635  to.SetComponents (std::cos(from.Angle()/2), std::sin(from.Angle()/2), 0, 0);
636 }
637 
638 void convert( RotationX const & from , RotationZYX & to )
639 {
640  // conversion from RotationX to RotationZYX
641  to.SetComponents(0,0,from.Angle());
642 }
643 
644 
645 // ----------------------------------------------------------------------
646 // conversions from RotationY
647 
648 void convert( RotationY const & from, Rotation3D & to)
649 {
650  // conversion from RotationY to Rotation3D
651 
652  const double c = from.CosAngle();
653  const double s = from.SinAngle();
654  to.SetComponents ( c, 0, s,
655  0, 1, 0,
656  -s, 0, c );
657 }
658 
659 void convert( RotationY const & from, AxisAngle & to)
660 {
661  // conversion from RotationY to AxisAngle
662 
664  to.SetComponents ( axis, from.Angle() );
665 }
666 
667 void convert( RotationY const & from, EulerAngles & to )
668 {
669  // conversion from RotationY to EulerAngles
670  // TODO better: temporary make conversion using Rotation3D
671 
672  Rotation3D tmp;
673  convert(from,tmp);
674  convert(tmp,to);
675 }
676 
677 void convert( RotationY const & from , RotationZYX & to )
678 {
679  // conversion from RotationY to RotationZYX
680  to.SetComponents(0,from.Angle(),0);
681 }
682 
683 
684 void convert( RotationY const & from, Quaternion & to)
685 {
686  // conversion from RotationY to Quaternion
687 
688  to.SetComponents (std::cos(from.Angle()/2), 0, std::sin(from.Angle()/2), 0);
689 }
690 
691 
692 
693 // ----------------------------------------------------------------------
694 // conversions from RotationZ
695 
696 void convert( RotationZ const & from, Rotation3D & to)
697 {
698  // conversion from RotationZ to Rotation3D
699 
700  const double c = from.CosAngle();
701  const double s = from.SinAngle();
702  to.SetComponents ( c, -s, 0,
703  s, c, 0,
704  0, 0, 1 );
705 }
706 
707 void convert( RotationZ const & from, AxisAngle & to)
708 {
709  // conversion from RotationZ to AxisAngle
710 
712  to.SetComponents ( axis, from.Angle() );
713 }
714 
715 void convert( RotationZ const & from , EulerAngles & to )
716 {
717  // conversion from RotationZ to EulerAngles
718  // TODO better: temporary make conversion using Rotation3D
719 
720  Rotation3D tmp;
721  convert(from,tmp);
722  convert(tmp,to);
723 }
724 
725 void convert( RotationZ const & from , RotationZYX & to )
726 {
727  // conversion from RotationY to RotationZYX
728  to.SetComponents(from.Angle(),0,0);
729 }
730 
731 void convert( RotationZ const & from, Quaternion & to)
732 {
733  // conversion from RotationZ to Quaternion
734 
735  to.SetComponents (std::cos(from.Angle()/2), 0, 0, std::sin(from.Angle()/2));
736 }
737 
738 } //namespace gv_detail
739 } //namespace Math
740 } //namespace ROOT
DisplacementVector3D< CoordSystem, Tag > & SetCoordinates(const Scalar src[])
Set internal data based on a C-style array of 3 Scalar numbers.
Scalar Phi() const
Return Phi Euler angle.
Definition: EulerAngles.h:212
void SetComponents(IT begin, IT end)
Set the three Euler angles given a pair of pointers or iterators defining the beginning and end of an...
Definition: EulerAngles.h:152
const double pi
Float_t theta
Definition: shapesAnim.C:5
return c
Rotation class representing a 3D rotation about the Z axis by the angle of rotation.
Definition: RotationZ.h:43
Scalar Y() const
Cartesian Y, converting if necessary from internal coordinate system.
Rotation class with the (3D) rotation represented by angles describing first a rotation of an angle p...
Definition: RotationZYX.h:71
Scalar SinAngle() const
Sine or Cosine of the rotation angle.
Definition: RotationZ.h:108
AxisVector Axis() const
accesss to rotation axis
Definition: AxisAngle.h:183
Scalar CosAngle() const
Definition: RotationX.h:109
Rotation class with the (3D) rotation represented by a unit quaternion (u, i, j, k).
Definition: Quaternion.h:47
void GetComponents(IT begin, IT end) const
Get the axis and then the angle into data specified by an iterator begin and another to the end of th...
Definition: RotationZYX.h:152
Scalar J() const
Definition: Quaternion.h:173
double cos(double)
TFile * f
Scalar K() const
Definition: Quaternion.h:174
AxisAngle class describing rotation represented with direction axis (3D Vector) and an angle of rotat...
Definition: AxisAngle.h:41
double sqrt(double)
Scalar Angle() const
Angle of rotation.
Definition: RotationZ.h:103
double acos(double)
Rotation class representing a 3D rotation about the Y axis by the angle of rotation.
Definition: RotationY.h:43
void SetComponents(IT begin, IT end)
Set the axis and then the angle given a pair of pointers or iterators defining the beginning and end ...
Definition: AxisAngle.h:117
Scalar X() const
Cartesian X, converting if necessary from internal coordinate system.
#define M_PI_2
Definition: Math.h:42
Scalar SinAngle() const
Sine or Cosine of the rotation angle.
Definition: RotationY.h:108
Scalar CosAngle() const
Definition: RotationY.h:109
double sin(double)
Scalar Angle() const
access to rotation angle
Definition: AxisAngle.h:188
Scalar I() const
Definition: Quaternion.h:172
void Rectify()
Re-adjust components to eliminate small deviations from the axis being a unit vector and angles out o...
Definition: AxisAngle.cxx:47
VecExpr< UnaryOp< Fabs< T >, VecExpr< A, T, D >, T >, T, D > fabs(const VecExpr< A, T, D > &rhs)
ROOT::R::TRInterface & r
Definition: Object.C:4
#define M_PI
Definition: Rotated.cxx:105
Rotation class representing a 3D rotation about the X axis by the angle of rotation.
Definition: RotationX.h:43
Rotation class with the (3D) rotation represented by a 3x3 orthogonal matrix.
Definition: Rotation3D.h:65
TMarker * m
Definition: textangle.C:8
tuple w
Definition: qtexample.py:51
static void correctByPi(double &psi, double &phi)
Scalar Angle() const
Angle of rotation.
Definition: RotationX.h:103
void GetComponents(ForeignVector &v1, ForeignVector &v2, ForeignVector &v3) const
Get components into three vectors which will be the (orthonormal) columns of the rotation matrix...
Definition: Rotation3D.h:249
Scalar U() const
Access to the four quaternion components: U() is the coefficient of the identity Pauli matrix...
Definition: Quaternion.h:171
double asin(double)
REAL epsilon
Definition: triangle.c:617
Scalar Theta() const
Return Theta Euler angle.
Definition: EulerAngles.h:222
Float_t phi
Definition: shapesAnim.C:6
void convert(R1 const &, R2 const)
Definition: 3DConversions.h:41
Scalar Z() const
Cartesian Z, converting if necessary from internal coordinate system.
double atan2(double, double)
EulerAngles class describing rotation as three angles (Euler Angles).
Definition: EulerAngles.h:43
Scalar CosAngle() const
Definition: RotationZ.h:109
Scalar SinAngle() const
Sine or Cosine of the rotation angle.
Definition: RotationX.h:108
Bool_t axis
Definition: geodemo.C:37
void SetComponents(const ForeignVector &v1, const ForeignVector &v2, const ForeignVector &v3)
Set components from three orthonormal vectors (which must have methods x(), y() and z()) which will b...
Definition: Rotation3D.h:233
void Rectify()
Re-adjust components to eliminate small deviations from |Q| = 1 orthonormality.
Definition: Quaternion.cxx:35
float * q
Definition: THbookFile.cxx:87
void SetComponents(IT begin, IT end)
Set the three Euler angles given a pair of pointers or iterators defining the beginning and end of an...
Definition: RotationZYX.h:135
Plane3D::Scalar Scalar
Definition: Plane3D.cxx:29
void SetComponents(IT begin, IT end)
Set the four components given an iterator to the start of the desired data, and another to the end (4...
Definition: Quaternion.h:112
Scalar Angle() const
Angle of rotation.
Definition: RotationY.h:103
Scalar Psi() const
Return Psi Euler angle.
Definition: EulerAngles.h:232