source/geometry/affinetransformation.h

   1 #ifndef MSP_GEOMETRY_AFFINETRANSFORMATION_H_
   2 #define MSP_GEOMETRY_AFFINETRANSFORMATION_H_
   3
   4 #include <msp/linal/squarematrix.h>
   5 #include "angle.h"
   6
   7 namespace Msp {
   8 namespace Geometry {
   9
  10 template<typename T, unsigned D>
  11 class AffineTransformation;
  12
  13
  14 /**
  15 Helper class to provide specialized operations for AffineTransformation.
  16 */
  17 template<typename T, unsigned D>
  18 class AffineTransformationOps
  19 {
  20 protected:
  21         AffineTransformationOps() { }
  22 };
  23
  24 template<typename T>
  25 class AffineTransformationOps<T, 2>
  26 {
  27 protected:
  28         AffineTransformationOps() { }
  29
  30 public:
  31         static AffineTransformation<T, 2> rotation(const Angle<T> &);
  32 };
  33
  34 template<typename T>
  35 class AffineTransformationOps<T, 3>
  36 {
  37 protected:
  38         AffineTransformationOps() { }
  39
  40 public:
  41         static AffineTransformation<T, 3> rotation(const Angle<T> &, const LinAl::Vector<T, 3> &);
  42 };
  43
  44
  45 /**
  46 An affine transformation in D dimensions.  Affine transformations preserve
  47 straightness of lines and ratios of distances.  Angles and distances themselves
  48 may change.  Internally this is represented by a square matrix of size D+1.
  49 */
  50 template<typename T, unsigned D>
  51 class AffineTransformation: public AffineTransformationOps<T, D>
  52 {
  53         friend class AffineTransformationOps<T, D>;
  54
  55 private:
  56         LinAl::SquareMatrix<T, D+1> matrix;
  57
  58 public:
  59         AffineTransformation();
  60
  61         static AffineTransformation<T, D> translation(const LinAl::Vector<T, D> &);
  62         static AffineTransformation<T, D> scaling(const LinAl::Vector<T, D> &);
  63         static AffineTransformation<T, D> shear(const LinAl::Vector<T, D> &, const LinAl::Vector<T, D> &);
  64
  65         const LinAl::SquareMatrix<T, D+1> &get_matrix() const { return matrix; }
  66         operator const LinAl::SquareMatrix<T, D+1> &() const { return matrix; }
  67
  68         LinAl::Vector<T, D> transform(const LinAl::Vector<T, D> &) const;
  69         LinAl::Vector<T, D> transform_linear(const LinAl::Vector<T, D> &) const;
  70 };
  71
  72 template<typename T, unsigned D>
  73 inline AffineTransformation<T, D>::AffineTransformation()
  74 {
  75         this->matrix = LinAl::SquareMatrix<T, D+1>::identity();
  76 }
  77
  78
  79 template<typename T, unsigned D>
  80 AffineTransformation<T, D> AffineTransformation<T, D>::translation(const LinAl::Vector<T, D> &v)
  81 {
  82         AffineTransformation<T, D> r;
  83         for(unsigned i=0; i<D; ++i)
  84                 r.matrix(i, D) = v[i];
  85         return r;
  86 }
  87
  88 template<typename T, unsigned D>
  89 AffineTransformation<T, D> AffineTransformation<T, D>::scaling(const LinAl::Vector<T, D> &factors)
  90 {
  91         AffineTransformation<T, D> r;
  92         for(unsigned i=0; i<D; ++i)
  93                 r.matrix(i, i) = factors[i];
  94         return r;
  95 }
  96
  97 template<typename T, unsigned D>
  98 AffineTransformation<T, D> AffineTransformation<T, D>::shear(const LinAl::Vector<T, D> &normal, const LinAl::Vector<T, D> &shift)
  99 {
 100         AffineTransformation<T, D> r;
 101         for(unsigned i=0; i<D; ++i)
 102                 for(unsigned j=0; j<D; ++j)
 103                         r.matrix(i, j) += normal[j]*shift[i];
 104         return r;
 105 }
 106
 107 template<typename T>
 108 AffineTransformation<T, 2> AffineTransformationOps<T, 2>::rotation(const Angle<T> &angle)
 109 {
 110         AffineTransformation<T, 2> r;
 111         T c = cos(angle);
 112         T s = sin(angle);
 113         r.matrix(0, 0) = c;
 114         r.matrix(0, 1) = -s;
 115         r.matrix(1, 0) = s;
 116         r.matrix(1, 1) = c;
 117         return r;
 118 }
 119
 120 template<typename T>
 121 AffineTransformation<T, 3> AffineTransformationOps<T, 3>::rotation(const Angle<T> &angle, const LinAl::Vector<T, 3> &axis)
 122 {
 123         AffineTransformation<T, 3> r;
 124         LinAl::Vector<T, 3> axn = normalize(axis);
 125         T c = cos(angle);
 126         T s = sin(angle);
 127         // http://en.wikipedia.org/wiki/Rotation_matrix#Rotation_matrix_from_axis_and_angle
 128         r.matrix(0, 0) = c+axn.x*axn.x*(1-c);
 129         r.matrix(0, 1) = axn.x*axn.y*(1-c)-axn.z*s;
 130         r.matrix(0, 2) = axn.x*axn.z*(1-c)+axn.y*s;
 131         r.matrix(1, 0) = axn.y*axn.x*(1-c)+axn.z*s;
 132         r.matrix(1, 1) = c+axn.y*axn.y*(1-c);
 133         r.matrix(1, 2) = axn.y*axn.z*(1-c)-axn.x*s;
 134         r.matrix(2, 0) = axn.z*axn.x*(1-c)-axn.y*s;
 135         r.matrix(2, 1) = axn.z*axn.y*(1-c)+axn.x*s;
 136         r.matrix(2, 2) = c+axn.z*axn.z*(1-c);
 137         return r;
 138 }
 139
 140
 141 template<typename T, unsigned N>
 142 inline LinAl::Vector<T, N+1> augment_vector(const LinAl::Vector<T, N> &v, T s)
 143 {
 144         LinAl::Vector<T, N+1> r;
 145         for(unsigned i=0; i<N; ++i)
 146                 r[i] = v[i];
 147         r[N] = s;
 148         return r;
 149 }
 150
 151 template<typename T, unsigned N>
 152 inline LinAl::Vector<T, N-1> reduce_vector(const LinAl::Vector<T, N> &v)
 153 {
 154         LinAl::Vector<T, N-1> r;
 155         for(unsigned i=0; i<N-1; ++i)
 156                 r[i] = v[i];
 157         return r;
 158 }
 159
 160 template<typename T, unsigned N>
 161 inline LinAl::Vector<T, N-1> divide_vector(const LinAl::Vector<T, N> &v)
 162 {
 163         LinAl::Vector<T, N-1> r;
 164         for(unsigned i=0; i<N-1; ++i)
 165                 r[i] = v[i]/v[N-1];
 166         return r;
 167 }
 168
 169
 170 template<typename T, unsigned D>
 171 inline LinAl::Vector<T, D> AffineTransformation<T, D>::transform(const LinAl::Vector<T, D> &v) const
 172 {
 173         return reduce_vector(matrix*augment_vector(v, T(1)));
 174 }
 175
 176 template<typename T, unsigned D>
 177 inline LinAl::Vector<T, D> AffineTransformation<T, D>::transform_linear(const LinAl::Vector<T, D> &v) const
 178 {
 179         return reduce_vector(matrix*augment_vector(v, T(0)));
 180 }
 181
 182 } // namespace Geometry
 183 } // namespace Msp
 184
 185 #endif