8 #ifndef IMPALGEBRA_ROTATION_3D_H
9 #define IMPALGEBRA_ROTATION_3D_H
11 #include <IMP/algebra/algebra_config.h>
22 IMPALGEBRA_BEGIN_NAMESPACE
24 #if !defined(IMP_DOXYGEN) && !defined(SWIG)
26 Rotation3D
compose(
const Rotation3D &a,
const Rotation3D &b) ;
48 mutable bool has_cache_;
52 void fill_cache()
const {
54 "Attempting to apply uninitialized rotation");
56 double v0s=get_squared(v_[0]);
57 double v1s=get_squared(v_[1]);
58 double v2s=get_squared(v_[2]);
59 double v3s=get_squared(v_[3]);
60 double v12= v_[1]*v_[2];
61 double v01=v_[0]*v_[1];
62 double v02=v_[0]*v_[2];
63 double v23=v_[2]*v_[3];
64 double v03=v_[0]*v_[3];
65 double v13=v_[1]*v_[3];
66 matrix_[0]=
Vector3D(v0s+v1s-v2s-v3s,
79 v_(v.get_unit_vector()),
87 Rotation3D(
double a,
double b,
double c,
double d): v_(a,b,c,d),
90 "Attempting to construct a rotation from a "
91 <<
" non-quaternion value. The coefficient vector"
92 <<
" must have a length of 1. Got: "
93 << a <<
" " << b <<
" " << c <<
" " << d
94 <<
" gives " << v_.get_squared_magnitude());
107 "Attempting to access uninitialized rotation");
108 return Vector3D((v_[0]*v_[0]+v_[1]*v_[1]-v_[2]*v_[2]-v_[3]*v_[3])*o[0]
109 + 2*(v_[1]*v_[2]-v_[0]*v_[3])*o[1]
110 + 2*(v_[1]*v_[3]+v_[0]*v_[2])*o[2],
111 2*(v_[1]*v_[2]+v_[0]*v_[3])*o[0]
112 + (v_[0]*v_[0]-v_[1]*v_[1]+v_[2]*v_[2]-v_[3]*v_[3])*o[1]
113 + 2*(v_[2]*v_[3]-v_[0]*v_[1])*o[2],
114 2*(v_[1]*v_[3]-v_[0]*v_[2])*o[0]
115 + 2*(v_[2]*v_[3]+v_[0]*v_[1])*o[1]
116 + (v_[0]*v_[0]-v_[1]*v_[1]-v_[2]*v_[2]+v_[3]*v_[3])*o[2]);
120 double get_rotated_one_coordinate_no_cache(
const Vector3D &o,
121 unsigned int coord)
const {
123 "Attempting to apply uninitialized rotation");
126 return (v_[0]*v_[0]+v_[1]*v_[1]-v_[2]*v_[2]-v_[3]*v_[3])*o[0]
127 + 2*(v_[1]*v_[2]-v_[0]*v_[3])*o[1]
128 + 2*(v_[1]*v_[3]+v_[0]*v_[2])*o[2];
131 return 2*(v_[1]*v_[2]+v_[0]*v_[3])*o[0]
132 + (v_[0]*v_[0]-v_[1]*v_[1]+v_[2]*v_[2]-v_[3]*v_[3])*o[1]
133 + 2*(v_[2]*v_[3]-v_[0]*v_[1])*o[2];
137 return 2*(v_[1]*v_[3]-v_[0]*v_[2])*o[0]
138 + 2*(v_[2]*v_[3]+v_[0]*v_[1])*o[1]
139 + (v_[0]*v_[0]-v_[1]*v_[1]-v_[2]*v_[2]+v_[3]*v_[3])*o[2];
142 IMP_THROW(
"Out of range coordinate " << coord,
150 "Attempting to apply uninitialized rotation");
151 if (!has_cache_) fill_cache();
158 double get_rotated_one_coordinate(
const Vector3D &o,
159 unsigned int coord)
const {
161 "Attempting to apply uninitialized rotation");
162 if (!has_cache_) fill_cache();
163 return o*matrix_[coord];
168 return get_rotated(v);
170 Vector3D get_rotation_matrix_row(
int i)
const {
172 if (!has_cache_) fill_cache();
181 "Attempting to invert uninitialized rotation");
182 Rotation3D ret(v_[0], -v_[1], -v_[2], -v_[3]);
193 "Attempting to access uninitialized rotation");
199 "Attempting to compose uninitialized rotation");
206 "Attempting to compose uninitialized rotation");
246 double ans= std::min(dot, odot);
247 const double s2=std::sqrt(2.0);
249 return std::max(std::min(ret, 1.0), 0.0);
277 IMPALGEBRAEXPORT Rotation3D
279 double m10,
double m11,
double m12,
280 double m20,
double m21,
double m22);
328 return Rotation3D(uv[0], uv[1], uv[2], uv[3]);
334 inline Rotation3D
compose(
const Rotation3D &a,
const Rotation3D &b) {
335 return Rotation3D(a.v_[0]*b.v_[0] - a.v_[1]*b.v_[1]
336 - a.v_[2]*b.v_[2] - a.v_[3]*b.v_[3],
337 a.v_[0]*b.v_[1] + a.v_[1]*b.v_[0]
338 + a.v_[2]*b.v_[3] - a.v_[3]*b.v_[2],
339 a.v_[0]*b.v_[2] - a.v_[1]*b.v_[3]
340 + a.v_[2]*b.v_[0] + a.v_[3]*b.v_[1],
341 a.v_[0]*b.v_[3] + a.v_[1]*b.v_[2]
342 - a.v_[2]*b.v_[1] + a.v_[3]*b.v_[0]);
410 FixedZYZ(
double rot,
double tilt,
double psi)
411 {v_[0]=rot; v_[1]= tilt; v_[2]=psi;}
412 double get_rot()
const {
415 double get_tilt()
const {
418 double get_psi()
const {
422 {out << v_[0] <<
" " << v_[1]
436 FixedZXZ(
double psi,
double theta,
double phi)
437 {v_[0]=psi; v_[1]= theta; v_[2]=phi;}
438 double get_psi()
const {
441 double get_theta()
const {
444 double get_phi()
const {
448 {out << v_[0] <<
" " << v_[1]
460 FixedXYZ(
double x,
double y,
double z)
461 {v_[0]=x; v_[1]= y; v_[2]=z;}
462 double get_x()
const {
465 double get_y()
const {
468 double get_z()
const {
472 out << v_[0] <<
" " << v_[1] <<
" " << v_[2];
513 if (bq*aq < 0) bq=-bq;
514 return f*aq+(1-f)*bq;
539 typedef std::pair<Vector3D,double> AxisAnglePair;
544 IMPALGEBRA_END_NAMESPACE