quat.hpp
上传用户:yhdzpy8989
上传日期:2007-06-13
资源大小:13604k
文件大小:10k
源码类别:

生物技术

开发平台:

C/C++

  1. /*
  2.  * ===========================================================================
  3.  * PRODUCTION $Log: quat.hpp,v $
  4.  * PRODUCTION Revision 1000.2  2004/06/01 19:49:02  gouriano
  5.  * PRODUCTION PRODUCTION: UPGRADED [GCC34_MSVC7] Dev-tree R1.7
  6.  * PRODUCTION
  7.  * ===========================================================================
  8.  */
  9. #ifndef GUI_MATH___QUAT___HPP
  10. #define GUI_MATH___QUAT___HPP
  11. /*  $Id: quat.hpp,v 1000.2 2004/06/01 19:49:02 gouriano Exp $
  12.  * ===========================================================================
  13.  *
  14.  *                            PUBLIC DOMAIN NOTICE
  15.  *               National Center for Biotechnology Information
  16.  *
  17.  *  This software / database is a "United States Government Work" under the
  18.  *  terms of the United States Copyright Act.  It was written as part of
  19.  *  the author's official duties as a United States Government employee and
  20.  *  thus cannot be copyrighted.  This software / database is freely available
  21.  *  to the public for use. The National Library of Medicine and the U.S.
  22.  *  Government have not placed any restriction on its use or reproduction.
  23.  *
  24.  *  Although all reasonable efforts have been taken to ensure the accuracy
  25.  *  and reliability of the software and data, the NLM and the U.S.
  26.  *  Government do not and cannot warrant the performance or results that
  27.  *  may be obtained by using this software or data. The NLM and the U.S.
  28.  *  Government disclaim all warranties, express or implied, including
  29.  *  warranties of performance, merchantability or fitness for any particular
  30.  *  purpose.
  31.  *
  32.  *  Please cite the author in any work or product based on this material.
  33.  *
  34.  * ===========================================================================
  35.  *
  36.  * Authors:  Mike DiCuccio
  37.  *
  38.  * File Description:
  39.  *
  40.  */
  41. #include <corelib/ncbistd.hpp>
  42. #include <math.h>
  43. #include <gui/math/vect3.hpp>
  44. #include <gui/math/matrix4.hpp>
  45. /** @addtogroup GUI_MATH
  46.  *
  47.  * @{
  48.  */
  49. BEGIN_NCBI_SCOPE
  50. template <class T>
  51. class CQuat
  52. {
  53. public:
  54.     // default ctor
  55.     CQuat(void);
  56.     // conversion ctor
  57.     CQuat(T x, T y, T z, T w);
  58.     // conversion ctor: rotation of degrees about axis
  59.     CQuat(const CVect3<T>& axis, T degrees);
  60.     // operator: index
  61.     T&          operator[] (int i)          { return m_Xyzw[i]; }
  62.     const T&    operator[] (int i) const    { return m_Xyzw[i]; }
  63.     // operator: *=
  64.     CQuat<T>&    operator *= (const CQuat<T>&);
  65.     // operator: +=
  66.     CQuat<T>&    operator += (const CQuat<T>&);
  67.     // operator: -=
  68.     CQuat<T>&    operator -= (const CQuat<T>&);
  69.     // access data through named functions
  70.     T&          X(void)         { return m_Xyzw[0]; }
  71.     const T&    X(void) const   { return m_Xyzw[0]; }
  72.     T&          Y(void)         { return m_Xyzw[1]; }
  73.     const T&    Y(void) const   { return m_Xyzw[1]; }
  74.     T&          Z(void)         { return m_Xyzw[2]; }
  75.     const T&    Z(void) const   { return m_Xyzw[2]; }
  76.     T&          W(void)         { return m_Xyzw[3]; }
  77.     const T&    W(void) const   { return m_Xyzw[3]; }
  78.     //
  79.     // other named functions
  80.     // create a rotation matric from a quaternion
  81.     void        ToMatrix(CMatrix4<T>&) const;
  82.     CMatrix4<T> ToMatrix(void) const;
  83.     // get the conjugate of a quaternion
  84.     CQuat<T>    Conjugate(void) const;
  85.     // return the magnitude of a quaternion
  86.     T           Length(void) const;
  87.     // return the magnitude squared of a quaternion
  88.     T           Length2(void) const;
  89.     // inverse of a quaterion(1 / q) = q_conj / (q * q_conj)
  90.     CQuat<T>    Inverse(void) const;
  91.     // rotate a vector using the current quaternion
  92.     void        Rotate(CVect3<float>& point) const;
  93. protected:  // data
  94.     // the data is stored as a 4 - element array
  95.     T           m_Xyzw[4];
  96. protected:  // internals
  97. };
  98. //
  99. // global operators
  100. //
  101. //
  102. // operator: compare - equals
  103. template <class T>
  104. inline bool
  105. operator== (const CQuat<T>& first, const CQuat<T>& second)
  106. {
  107.     return (first[0] == second[0]  &&  
  108.             first[1] == second[1]  &&  
  109.             first[2] == second[2]  &&  
  110.             first[3] == second[3]);
  111. }
  112. //
  113. // operator: compare - not - equals
  114. template <class T>
  115. inline bool
  116. operator!= (const CQuat<T>& first, const CQuat<T>& second)
  117. {
  118.     return !(first == second);
  119. }
  120. //
  121. // operator: multiply
  122. template <class T>
  123. inline CQuat<T>
  124. operator* (const CQuat<T>& first, const CQuat<T>& second)
  125. {
  126.     return (CQuat<T> (first) *= second);
  127. }
  128. //
  129. // operator: multiply
  130. template <class T>
  131. inline CQuat<T>
  132. operator* (const CQuat<T>& first, T second)
  133. {
  134.     return (CQuat<T> (first.X() * second,
  135.                         first.Y() * second,
  136.                         first.Z() * second,
  137.                         first.W() * second));
  138. }
  139. //
  140. // operator: multiply
  141. template <class T>
  142. inline CQuat<T>
  143. operator* (T first, const CQuat<T>& second)
  144. {
  145.     return (CQuat<T> (second.X() * first,
  146.                         second.Y() * first,
  147.                         second.Z() * first,
  148.                         second.W() * first));
  149. }
  150. //
  151. // operator: divide
  152. template <class T>
  153. inline CQuat<T>
  154. operator/ (const CQuat<T>& first, T second)
  155. {
  156.     return (CQuat<T> (first.X() / second,
  157.                         first.Y() / second,
  158.                         first.Z() / second,
  159.                         first.W() / second));
  160. }
  161. //
  162. // operator: add
  163. template <class T>
  164. inline CQuat<T>
  165. operator+ (const CQuat<T>& first, const CQuat<T>& second)
  166. {
  167.     return (CQuat<T> (first) += second);
  168. }
  169. //
  170. // operator: subtract
  171. template <class T>
  172. inline CQuat<T>
  173. operator- (const CQuat<T>& first, const CQuat<T>& second)
  174. {
  175.     return (CQuat<T> (first) -= second);
  176. }
  177. //
  178. // default ctor
  179. template <class T>
  180. inline
  181. CQuat<T>::CQuat()
  182. {
  183.     m_Xyzw[0] = T(0);
  184.     m_Xyzw[1] = T(0);
  185.     m_Xyzw[2] = T(0);
  186.     m_Xyzw[3] = T(0);
  187. }
  188. //
  189. // conversion ctor
  190. template <class T>
  191. inline
  192. CQuat<T>::CQuat(T x, T y, T z, T w)
  193. {
  194.     m_Xyzw[0] = x;
  195.     m_Xyzw[1] = y;
  196.     m_Xyzw[2] = z;
  197.     m_Xyzw[3] = w;
  198. }
  199. //
  200. // conversion ctor: rotation quat
  201. // we want to create a rotation quat about an axis
  202. template <class T>
  203. inline
  204. CQuat<T>::CQuat(const CVect3<T>& axis, T theta)
  205. {
  206.     theta *= math::pi / 180.0f;
  207.     double cos_theta = cos(theta / 2.0);
  208.     double sin_theta = sin(theta / 2.0);
  209.     m_Xyzw[0] = axis[0]*sin_theta;
  210.     m_Xyzw[1] = axis[1]*sin_theta;
  211.     m_Xyzw[2] = axis[2]*sin_theta;
  212.     m_Xyzw[3] = cos_theta;
  213. }
  214. //
  215. // operator: *=
  216. template <class T>
  217. inline CQuat<T>&
  218. CQuat<T>::operator*= (const CQuat<T>& q)
  219. {
  220.     T new_x = W() * q.X() +  X() * q.W() + ( Y() * q.Z() - Z() * q.Y() );
  221.     T new_y = W() * q.Y() +  Y() * q.W() + ( Z() * q.X() - X() * q.Z() );
  222.     T new_z = W() * q.Z() +  Z() * q.W() + ( X() * q.Y() - Y() * q.X() );
  223.     T new_w = W() * q.W() - (X() * q.X() + Y() * q.Y() + Z() * q.Z());
  224.     m_Xyzw[0] = new_x;
  225.     m_Xyzw[1] = new_y;
  226.     m_Xyzw[2] = new_z;
  227.     m_Xyzw[3] = new_w;
  228.     return *this;
  229. }
  230. //
  231. // operator: -=
  232. template <class T>
  233. inline CQuat<T>&
  234. CQuat<T>::operator-= (const CQuat<T>& q)
  235. {
  236.     X() -= q.X();
  237.     Y() -= q.Y();
  238.     Z() -= q.Z();
  239.     W() -= q.W();
  240.     return *this;
  241. }
  242. //
  243. // operator: +=
  244. template <class T>
  245. inline CQuat<T>&
  246. CQuat<T>::operator+= (const CQuat<T>& q)
  247. {
  248.     X() += q.X();
  249.     Y() += q.Y();
  250.     Z() += q.Z();
  251.     W() += q.W();
  252.     return *this;
  253. }
  254. //
  255. // create a rotation matrix from a quaternion
  256. template <class T>
  257. inline CMatrix4<T>
  258. CQuat<T>::ToMatrix() const
  259. {
  260.     CMatrix4<T> m;
  261.     ToMatrix(m);
  262.     return m;
  263. }
  264. template <class T>
  265. inline void
  266. CQuat<T>::ToMatrix(CMatrix4<T>& m) const
  267. {
  268.     T Nq = m_Xyzw[0]*m_Xyzw[0] +
  269.         m_Xyzw[1]*m_Xyzw[1] +
  270.         m_Xyzw[2]*m_Xyzw[2] +
  271.         m_Xyzw[3]*m_Xyzw[3];
  272.     T s = (Nq > 0.0) ? (2.0 / Nq) : 0.0;
  273.     T xs = m_Xyzw[0]*s;
  274.     T ys = m_Xyzw[1]*s;
  275.     T zs = m_Xyzw[2]*s;
  276.     T wx = m_Xyzw[3]*xs;
  277.     T wy = m_Xyzw[3]*ys;
  278.     T wz = m_Xyzw[3]*zs;
  279.     T xx = m_Xyzw[0]*xs;
  280.     T xy = m_Xyzw[0]*ys;
  281.     T xz = m_Xyzw[0]*zs;
  282.     T yy = m_Xyzw[1]*ys;
  283.     T yz = m_Xyzw[1]*zs;
  284.     T zz = m_Xyzw[2]*zs;
  285.     m(0,0) = 1.0 - (yy + zz);
  286.     m(1,0) = xy + wz;
  287.     m(2,0) = xz - wy;
  288.     m(0,1) = xy - wz;
  289.     m(1,1) = 1.0 - (xx + zz);
  290.     m(2,1) = yz + wx;
  291.     m(0,2) = xz + wy;
  292.     m(1,2) = yz - wx;
  293.     m(2,2) = 1.0 - (xx + yy);
  294.     m(0,3) = 0.0;
  295.     m(1,3) = 0.0;
  296.     m(2,3) = 0.0;
  297.     m(3,0) = 0.0;
  298.     m(3,1) = 0.0;
  299.     m(3,2) = 0.0;
  300.     m(3,3) = 1.0;
  301. }
  302. //
  303. // return the conjugate of a quaternion
  304. template <class T>
  305. inline CQuat<T>
  306. CQuat<T>::Conjugate() const
  307. {
  308.     CQuat<T> q(*this);
  309.     q.X() = -q.X();
  310.     q.Y() = -q.Y();
  311.     q.Z() = -q.Z();
  312.     return q;
  313. }
  314. //
  315. // get the length(magnitude) of a quaternion
  316. template <class T>
  317. inline T
  318. CQuat<T>::Length() const
  319. {
  320.     return sqrt(Length2());
  321. }
  322. //
  323. // get the length(magnitude) squared of a quaternion
  324. template <class T>
  325. inline T
  326. CQuat<T>::Length2() const
  327. {
  328.     return (m_Xyzw[0]*m_Xyzw[0] +
  329.             m_Xyzw[1]*m_Xyzw[1] +
  330.             m_Xyzw[2]*m_Xyzw[2] +
  331.             m_Xyzw[3]*m_Xyzw[3]);
  332. }
  333. //
  334. // get the inverse of a quaternion
  335. template <class T>
  336. inline CQuat<T>
  337. CQuat<T>::Inverse() const
  338. {
  339.     return (Conjugate() / Length2());
  340. }
  341. //
  342. // rotate a point using a quaternion
  343. template <class T>
  344. inline void
  345. CQuat<T>::Rotate(CVect3<float>& point) const
  346. {
  347.     // formula is qPq', where
  348.     // P = (0,p)
  349.     CQuat<T> p(point.X(), point.Y(), point.Z(), 0);
  350.     // there's a little extra math done in here
  351.     // W() is 0, so we could eliminate this
  352.     p = *this * p;
  353.     p *= Inverse();
  354.     point.X() = p.X();
  355.     point.Y() = p.Y();
  356.     point.Z() = p.Z();
  357.     // w remains 0 throughout
  358. }
  359. //
  360. // streams - based output
  361. template <class T>
  362. ostream&
  363. operator<< (ostream& o, const CQuat<T>& quat)
  364. {
  365.     o << "[(" << quat.X() << ", " << quat.Y() << ", " << quat.Z()
  366.         << "), " << quat.W() << "]";
  367.     return o;
  368. }
  369. END_NCBI_SCOPE
  370. /* @} */
  371. /*
  372.  * ===========================================================================
  373.  * $Log: quat.hpp,v $
  374.  * Revision 1000.2  2004/06/01 19:49:02  gouriano
  375.  * PRODUCTION: UPGRADED [GCC34_MSVC7] Dev-tree R1.7
  376.  *
  377.  * Revision 1.7  2004/05/11 18:53:50  dicuccio
  378.  * Added doxygen modules info
  379.  *
  380.  * Revision 1.6  2004/05/03 13:10:01  dicuccio
  381.  * Removed gui/gui.hpp --> not necessary for math projects
  382.  *
  383.  * Revision 1.5  2004/05/03 12:43:35  dicuccio
  384.  * Added #include for gui/gui.hpp
  385.  *
  386.  * Revision 1.4  2004/04/26 17:29:06  ucko
  387.  * Length(): Use plain sqrt rather than math::Sqrt, which seems not to exist.
  388.  *
  389.  * Revision 1.3  2004/03/11 17:17:07  dicuccio
  390.  * FIxed include guards
  391.  *
  392.  * Revision 1.2  2003/12/22 19:14:21  dicuccio
  393.  * Fixed lots of bugs in referencing non-existent functions
  394.  *
  395.  * Revision 1.1  2003/06/09 19:30:50  dicuccio
  396.  * Initial revision
  397.  *
  398.  * ===========================================================================
  399.  */
  400. #endif  // GUI_MATH___QUAT___HPP