rotmatrix.h

Go to the documentation of this file.
00001 /*
00002  *  This file is part of Healpix_cxx.
00003  *
00004  *  Healpix_cxx is free software; you can redistribute it and/or modify
00005  *  it under the terms of the GNU General Public License as published by
00006  *  the Free Software Foundation; either version 2 of the License, or
00007  *  (at your option) any later version.
00008  *
00009  *  Healpix_cxx is distributed in the hope that it will be useful,
00010  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
00011  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00012  *  GNU General Public License for more details.
00013  *
00014  *  You should have received a copy of the GNU General Public License
00015  *  along with Healpix_cxx; if not, write to the Free Software
00016  *  Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301  USA
00017  *
00018  *  For more information about HEALPix, see http://healpix.jpl.nasa.gov
00019  */
00020 
00021 /*
00022  *  Healpix_cxx is being developed at the Max-Planck-Institut fuer Astrophysik
00023  *  and financially supported by the Deutsches Zentrum fuer Luft- und Raumfahrt
00024  *  (DLR).
00025  */
00026 
00027 /*! \file rotmatrix.h
00028  *  Class for rotation transforms in 3D space
00029  *
00030  *  Copyright (C) 2003 Max-Planck-Society
00031  *  \author Martin Reinecke
00032  */
00033 
00034 #ifndef PLANCK_ROTMATRIX_H
00035 #define PLANCK_ROTMATRIX_H
00036 
00037 #include <iostream>
00038 #include "cxxutils.h"
00039 #include "vec3.h"
00040 
00041 /*! \defgroup rotmatrixgroup Rotation matrices */
00042 /*! \{ */
00043 
00044 /*! Class for rotation transforms in 3D space */
00045 class rotmatrix
00046   {
00047   public:
00048     double entry[3][3];
00049 
00050     rotmatrix () {}
00051 
00052     /*! Constructs a rotation matrix from its nine entries */
00053     rotmatrix (double a00, double a01, double a02,
00054                double a10, double a11, double a12,
00055                double a20, double a21, double a22)
00056       {
00057       entry[0][0]=a00; entry[0][1]=a01; entry[0][2]=a02;
00058       entry[1][0]=a10; entry[1][1]=a11; entry[1][2]=a12;
00059       entry[2][0]=a20; entry[2][1]=a21; entry[2][2]=a22;
00060       }
00061 
00062     /*! Constructs a rotation matrix so that \a a is the first column,
00063         \a b is the second column and \a c is the third column.
00064         \note The vectors \a a, \a b and \a c must form an orthonormal system!
00065      */
00066     rotmatrix (const vec3 &a, const vec3 &b, const vec3 &c);
00067 
00068     /*! Sets the matrix to the identity matrix. */
00069     void SetToIdentity ();
00070     /*! Sets all matrix elements to zero. */
00071     void SetToZero ();
00072     /*! Transposes the matrix. */
00073     void Transpose ();
00074 
00075     /*! Extracts a unit-length rotation axis \a axis and a rotation angle
00076         \a angle from the matrix. */
00077     void toAxisAngle (vec3 &axis, double &angle) const;
00078 
00079     /*! Constructs a matrix which causes a rotation by \a angle around
00080         \a axis. \a axis must have unit length. */
00081     void Make_Axis_Rotation_Transform (const vec3 &axis, double angle)
00082       {
00083       double sa=sin(angle), ca=cos(angle);
00084       double ica=1-ca;
00085       entry[0][0] = axis.x*axis.x*ica + ca;
00086       entry[1][1] = axis.y*axis.y*ica + ca;
00087       entry[2][2] = axis.z*axis.z*ica + ca;
00088       double t1 = axis.x*axis.y*ica, t2 = axis.z*sa;
00089       entry[1][0] = t1 + t2;
00090       entry[0][1] = t1 - t2;
00091       t1 = axis.x*axis.z*ica; t2 = axis.y*sa;
00092       entry[2][0] = t1 - t2;
00093       entry[0][2] = t1 + t2;
00094       t1 = axis.y*axis.z*ica; t2 = axis.x*sa;
00095       entry[1][2] = t1 - t2;
00096       entry[2][1] = t1 + t2;
00097       }
00098 
00099     /*! Creates a rotation matrix \a A, which performs the following operations
00100         on a vector \a v, when \a Av is calculated:
00101         -# rotate \a v around the z-axis by \a gamma,
00102         -# rotate \a v' around the y-axis by \a beta,
00103         -# rotate \a v'' around the z-axis by \a alpha.
00104 
00105         \note \a alpha, \a beta and \a gamma are given in radians,
00106               the rotations are right handed.
00107 
00108         \note This transformation rotates the \e vectors, not the coordinate
00109               axes! */
00110     void Make_CPAC_Euler_Matrix (double alpha, double beta, double gamma);
00111 
00112     /*! Extracts the Euler angles \a alpha, \a beta and \a gamma from the
00113         matrix. For their definition see Make_CPAC_Euler_Matrix().
00114 
00115         \note In case of ambiguity \a alpha will be 0. */
00116     void Extract_CPAC_Euler_Angles
00117       (double &alpha, double &beta, double &gamma) const;
00118 
00119     /*! Returns the vector \a vec, transformed by the matrix. */
00120     vec3 Transform (const vec3 &vec) const
00121       {
00122       return vec3
00123         (vec.x*entry[0][0] + vec.y*entry[0][1] + vec.z*entry[0][2],
00124          vec.x*entry[1][0] + vec.y*entry[1][1] + vec.z*entry[1][2],
00125          vec.x*entry[2][0] + vec.y*entry[2][1] + vec.z*entry[2][2]);
00126       }
00127     /*! Returns the vector \a vec, transformed by the matrix, in \a vec2. */
00128     void Transform (const vec3 &vec, vec3 &vec2) const
00129       {
00130       vec2.x = vec.x*entry[0][0] + vec.y*entry[0][1] + vec.z*entry[0][2];
00131       vec2.y = vec.x*entry[1][0] + vec.y*entry[1][1] + vec.z*entry[1][2];
00132       vec2.z = vec.x*entry[2][0] + vec.y*entry[2][1] + vec.z*entry[2][2];
00133       }
00134   };
00135 
00136 /*! Returns \a a * \a b.
00137     \relates rotmatrix */
00138 rotmatrix operator* (const rotmatrix &a, const rotmatrix &b);
00139 /*! Returns \a a * \a b in \a res.
00140     \relates rotmatrix */
00141 inline void matmult (const rotmatrix &a, const rotmatrix &b, rotmatrix &res)
00142   {
00143   for (int i=0; i<3; ++i)
00144     for (int j=0; j<3; ++j)
00145       res.entry[i][j] = a.entry[i][0] * b.entry[0][j]
00146                       + a.entry[i][1] * b.entry[1][j]
00147                       + a.entry[i][2] * b.entry[2][j];
00148   }
00149 
00150 /*! Returns \a a^T * \a b in \a res.
00151     \relates rotmatrix */
00152 void TransposeTimes (const rotmatrix &a, const rotmatrix &b, rotmatrix &res);
00153 
00154 /*! Writes \a mat to \a os.
00155     \relates rotmatrix */
00156 std::ostream &operator<< (std::ostream &os, const rotmatrix &mat);
00157 
00158 /*! \} */
00159 
00160 #endif