1 | // -*- C++ -*-
|
---|
2 | // ---------------------------------------------------------------------------
|
---|
3 | //
|
---|
4 | // This file is a part of the CLHEP - a Class Library for High Energy Physics.
|
---|
5 | //
|
---|
6 | // This is the implementation of those methods of the HepRotation class which
|
---|
7 | // were introduced when ZOOM PhysicsVectors was merged in, and which involve
|
---|
8 | // the angle/axis representation of a Rotation.
|
---|
9 | //
|
---|
10 |
|
---|
11 | #ifdef GNUPRAGMA
|
---|
12 | #pragma implementation
|
---|
13 | #endif
|
---|
14 |
|
---|
15 | #include "CLHEP/Vector/defs.h"
|
---|
16 | #include "CLHEP/Vector/Rotation.h"
|
---|
17 | #include "CLHEP/Units/PhysicalConstants.h"
|
---|
18 |
|
---|
19 | #include <iostream>
|
---|
20 | #include <cmath>
|
---|
21 |
|
---|
22 | namespace CLHEP {
|
---|
23 |
|
---|
24 | // ---------- Constructors and Assignment:
|
---|
25 |
|
---|
26 | // axis and angle
|
---|
27 |
|
---|
28 | HepRotation & HepRotation::set( const Hep3Vector & axis, double delta ) {
|
---|
29 |
|
---|
30 | register double sinDelta = sin(delta), cosDelta = cos(delta);
|
---|
31 | register double oneMinusCosDelta = 1.0 - cosDelta;
|
---|
32 |
|
---|
33 | Hep3Vector u = axis.unit();
|
---|
34 |
|
---|
35 | register double uX = u.getX();
|
---|
36 | register double uY = u.getY();
|
---|
37 | register double uZ = u.getZ();
|
---|
38 |
|
---|
39 | rxx = oneMinusCosDelta * uX * uX + cosDelta;
|
---|
40 | rxy = oneMinusCosDelta * uX * uY - sinDelta * uZ;
|
---|
41 | rxz = oneMinusCosDelta * uX * uZ + sinDelta * uY;
|
---|
42 |
|
---|
43 | ryx = oneMinusCosDelta * uY * uX + sinDelta * uZ;
|
---|
44 | ryy = oneMinusCosDelta * uY * uY + cosDelta;
|
---|
45 | ryz = oneMinusCosDelta * uY * uZ - sinDelta * uX;
|
---|
46 |
|
---|
47 | rzx = oneMinusCosDelta * uZ * uX - sinDelta * uY;
|
---|
48 | rzy = oneMinusCosDelta * uZ * uY + sinDelta * uX;
|
---|
49 | rzz = oneMinusCosDelta * uZ * uZ + cosDelta;
|
---|
50 |
|
---|
51 | return *this;
|
---|
52 |
|
---|
53 | } // HepRotation::set(axis, delta)
|
---|
54 |
|
---|
55 | HepRotation::HepRotation ( const Hep3Vector & axis, double delta )
|
---|
56 | {
|
---|
57 | set( axis, delta );
|
---|
58 | }
|
---|
59 | HepRotation & HepRotation::set( const HepAxisAngle & ax ) {
|
---|
60 | return set ( ax.axis(), ax.delta() );
|
---|
61 | }
|
---|
62 | HepRotation::HepRotation ( const HepAxisAngle & ax )
|
---|
63 | {
|
---|
64 | set ( ax.axis(), ax.delta() );
|
---|
65 | }
|
---|
66 |
|
---|
67 | |
---|
68 |
|
---|
69 |
|
---|
70 | double HepRotation::delta() const {
|
---|
71 |
|
---|
72 | double cosdelta = (rxx + ryy + rzz - 1.0) / 2.0;
|
---|
73 | if (cosdelta > 1.0) {
|
---|
74 | return 0;
|
---|
75 | } else if (cosdelta < -1.0) {
|
---|
76 | return CLHEP::pi;
|
---|
77 | } else {
|
---|
78 | return acos( cosdelta ); // Already safe due to the cosdelta > 1 check
|
---|
79 | }
|
---|
80 |
|
---|
81 | } // delta()
|
---|
82 |
|
---|
83 | Hep3Vector HepRotation::axis () const {
|
---|
84 |
|
---|
85 | // Determine 2*sin(delta) times the u components (I call this uX, uY, Uz)
|
---|
86 | // Normalization is not needed; it will be done when returning the 3-Vector
|
---|
87 |
|
---|
88 | double Uz = ryx - rxy;
|
---|
89 | double Uy = rxz - rzx;
|
---|
90 | double Ux = rzy - ryz;
|
---|
91 |
|
---|
92 | if ( (Uz==0) && (Uy==0) && (Ux==0) ) {
|
---|
93 | if ( rzz>0 ) {
|
---|
94 | return Hep3Vector(0,0,1);
|
---|
95 | } else if ( ryy>0 ) {
|
---|
96 | return Hep3Vector(0,1,0);
|
---|
97 | } else {
|
---|
98 | return Hep3Vector(1,0,0);
|
---|
99 | }
|
---|
100 | } else {
|
---|
101 | return Hep3Vector( Ux, Uy, Uz ).unit();
|
---|
102 | }
|
---|
103 |
|
---|
104 | } // axis()
|
---|
105 |
|
---|
106 | HepAxisAngle HepRotation::axisAngle() const {
|
---|
107 |
|
---|
108 | return HepAxisAngle (axis(), delta());
|
---|
109 |
|
---|
110 | } // axisAngle()
|
---|
111 |
|
---|
112 |
|
---|
113 | void HepRotation::setAxis (const Hep3Vector & axis) {
|
---|
114 | set ( axis, delta() );
|
---|
115 | }
|
---|
116 |
|
---|
117 | void HepRotation::setDelta (double delta) {
|
---|
118 | set ( axis(), delta );
|
---|
119 | }
|
---|
120 |
|
---|
121 | } // namespace CLHEP
|
---|