--- /dev/null
+#ifndef MSP_GEOMETRY_AFFINETRANSFORM_H_
+#define MSP_GEOMETRY_AFFINETRANSFORM_H_
+
+#include <msp/linal/matrix.h>
+#include "angle.h"
+#include "boundingbox.h"
+#include "ray.h"
+
+namespace Msp {
+namespace Geometry {
+
+template<typename T, unsigned D>
+class AffineTransform;
+
+
+/**
+Helper class to provide specialized operations for AffineTransform.
+*/
+template<typename T, unsigned D>
+class AffineTransformOps
+{
+protected:
+ AffineTransformOps() { }
+};
+
+template<typename T>
+class AffineTransformOps<T, 2>
+{
+protected:
+ AffineTransformOps() { }
+
+public:
+ static AffineTransform<T, 2> rotation(const Angle<T> &);
+};
+
+template<typename T>
+class AffineTransformOps<T, 3>
+{
+protected:
+ AffineTransformOps() { }
+
+public:
+ static AffineTransform<T, 3> rotation(const Angle<T> &, const LinAl::Vector<T, 3> &);
+};
+
+
+/**
+An affine transformation in D dimensions. Affine transformations preserve
+straightness of lines and ratios of distances. Angles and distances themselves
+may change. Internally this is represented by a square matrix of size D+1.
+*/
+template<typename T, unsigned D>
+class AffineTransform: public AffineTransformOps<T, D>
+{
+ friend class AffineTransformOps<T, D>;
+
+private:
+ LinAl::Matrix<T, D+1, D+1> matrix;
+
+public:
+ AffineTransform();
+
+ static AffineTransform<T, D> translation(const LinAl::Vector<T, D> &);
+ static AffineTransform<T, D> scaling(const LinAl::Vector<T, D> &);
+ static AffineTransform<T, D> shear(const LinAl::Vector<T, D> &, const LinAl::Vector<T, D> &);
+
+ AffineTransform &operator*=(const AffineTransform &);
+ AffineTransform &invert();
+
+ const LinAl::Matrix<T, D+1, D+1> &get_matrix() const { return matrix; }
+ operator const LinAl::Matrix<T, D+1, D+1> &() const { return matrix; }
+
+ LinAl::Vector<T, D> transform(const LinAl::Vector<T, D> &) const;
+ LinAl::Vector<T, D> transform_linear(const LinAl::Vector<T, D> &) const;
+ Ray<T, D> transform(const Ray<T, D> &) const;
+ BoundingBox<T, D> transform(const BoundingBox<T, D> &) const;
+};
+
+template<typename T, unsigned D>
+inline AffineTransform<T, D>::AffineTransform()
+{
+ this->matrix = LinAl::Matrix<T, D+1, D+1>::identity();
+}
+
+
+template<typename T, unsigned D>
+AffineTransform<T, D> AffineTransform<T, D>::translation(const LinAl::Vector<T, D> &v)
+{
+ AffineTransform<T, D> r;
+ for(unsigned i=0; i<D; ++i)
+ r.matrix(i, D) = v[i];
+ return r;
+}
+
+template<typename T, unsigned D>
+AffineTransform<T, D> AffineTransform<T, D>::scaling(const LinAl::Vector<T, D> &factors)
+{
+ AffineTransform<T, D> r;
+ for(unsigned i=0; i<D; ++i)
+ r.matrix(i, i) = factors[i];
+ return r;
+}
+
+template<typename T, unsigned D>
+AffineTransform<T, D> AffineTransform<T, D>::shear(const LinAl::Vector<T, D> &normal, const LinAl::Vector<T, D> &shift)
+{
+ AffineTransform<T, D> r;
+ for(unsigned i=0; i<D; ++i)
+ for(unsigned j=0; j<D; ++j)
+ r.matrix(i, j) += normal[j]*shift[i];
+ return r;
+}
+
+template<typename T>
+AffineTransform<T, 2> AffineTransformOps<T, 2>::rotation(const Angle<T> &angle)
+{
+ AffineTransform<T, 2> r;
+ T c = cos(angle);
+ T s = sin(angle);
+ r.matrix(0, 0) = c;
+ r.matrix(0, 1) = -s;
+ r.matrix(1, 0) = s;
+ r.matrix(1, 1) = c;
+ return r;
+}
+
+template<typename T>
+AffineTransform<T, 3> AffineTransformOps<T, 3>::rotation(const Angle<T> &angle, const LinAl::Vector<T, 3> &axis)
+{
+ AffineTransform<T, 3> r;
+ LinAl::Vector<T, 3> axn = normalize(axis);
+ T c = cos(angle);
+ T s = sin(angle);
+ // http://en.wikipedia.org/wiki/Rotation_matrix#Rotation_matrix_from_axis_and_angle
+ r.matrix(0, 0) = c+axn.x*axn.x*(1-c);
+ r.matrix(0, 1) = axn.x*axn.y*(1-c)-axn.z*s;
+ r.matrix(0, 2) = axn.x*axn.z*(1-c)+axn.y*s;
+ r.matrix(1, 0) = axn.y*axn.x*(1-c)+axn.z*s;
+ r.matrix(1, 1) = c+axn.y*axn.y*(1-c);
+ r.matrix(1, 2) = axn.y*axn.z*(1-c)-axn.x*s;
+ r.matrix(2, 0) = axn.z*axn.x*(1-c)-axn.y*s;
+ r.matrix(2, 1) = axn.z*axn.y*(1-c)+axn.x*s;
+ r.matrix(2, 2) = c+axn.z*axn.z*(1-c);
+ return r;
+}
+
+template<typename T, unsigned D>
+inline AffineTransform<T, D> &AffineTransform<T, D>::operator*=(const AffineTransform<T, D> &other)
+{
+ matrix *= other.get_matrix();
+ return *this;
+}
+
+template<typename T, unsigned D>
+inline AffineTransform<T, D> operator*(const AffineTransform<T, D> &at1, const AffineTransform<T, D> &at2)
+{
+ AffineTransform<T, D> r = at1;
+ return r *= at2;
+}
+
+template<typename T, unsigned D>
+inline AffineTransform<T, D> &AffineTransform<T, D>::invert()
+{
+ matrix.invert();
+ return *this;
+}
+
+template<typename T, unsigned D>
+inline AffineTransform<T, D> invert(const AffineTransform<T, D> &at)
+{
+ AffineTransform<T, D> r = at;
+ return r.invert();
+}
+
+template<typename T, unsigned D>
+inline LinAl::Vector<T, D> AffineTransform<T, D>::transform(const LinAl::Vector<T, D> &v) const
+{
+ return (matrix*compose(v, T(1))).template slice<D>(0);
+}
+
+template<typename T, unsigned D>
+inline LinAl::Vector<T, D> AffineTransform<T, D>::transform_linear(const LinAl::Vector<T, D> &v) const
+{
+ return (matrix*compose(v, T(0))).template slice<D>(0);
+}
+
+template<typename T, unsigned D>
+inline Ray<T, D> AffineTransform<T, D>::transform(const Ray<T, D> &ray) const
+{
+ LinAl::Vector<T, D> dir = transform_linear(ray.get_direction());
+ return Ray<T, D>(transform(ray.get_start()), dir, ray.get_limit()*dir.norm());
+}
+
+template<typename T, unsigned D>
+inline BoundingBox<T, D> AffineTransform<T, D>::transform(const BoundingBox<T, D> &bbox) const
+{
+ LinAl::Vector<T, D> min_pt;
+ LinAl::Vector<T, D> max_pt;
+ for(unsigned i=0; i<(1<<D); ++i)
+ {
+ LinAl::Vector<T, D> point;
+ for(unsigned j=0; j<D; ++j)
+ point[j] = ((i>>j)&1 ? bbox.get_maximum_coordinate(j) : bbox.get_minimum_coordinate(j));
+
+ point = transform(point);
+
+ if(i==0)
+ {
+ min_pt = point;
+ max_pt = point;
+ }
+ else
+ {
+ for(unsigned j=0; j<D; ++j)
+ {
+ min_pt[j] = std::min(min_pt[j], point[j]);
+ max_pt[j] = std::max(max_pt[j], point[j]);
+ }
+ }
+ }
+
+ return BoundingBox<T, D>(min_pt, max_pt);
+}
+
+} // namespace Geometry
+} // namespace Msp
+
+#endif
projects / libs / math.git / blobdiff