transform.h

Go to the documentation of this file.

/*
    This file is part of Mitsuba, a physically based rendering system.

    Copyright (c) 2007-2014 by Wenzel Jakob and others.

    Mitsuba is free software; you can redistribute it and/or modify
    it under the terms of the GNU General Public License Version 3
    as published by the Free Software Foundation.

    Mitsuba is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with this program. If not, see <http://www.gnu.org/licenses/>.
*/

#pragma once
#if !defined(__MITSUBA_CORE_TRANSFORM_H_)
#define __MITSUBA_CORE_TRANSFORM_H_

#include <mitsuba/core/matrix.h>
#include <mitsuba/core/ray.h>

MTS_NAMESPACE_BEGIN

/**
 * \brief Encapsulates a 4x4 linear transformation and its inverse
 * \ingroup libcore
 * \ingroup libpython
 */
struct MTS_EXPORT_CORE Transform {
public:
        /// Create an identity transformation
        Transform() {
                m_transform.setIdentity();
                m_invTransform.setIdentity();
        }

        /// Unserialize a transformation from a stream
        inline Transform(Stream *stream) {
                m_transform = Matrix4x4(stream);
                m_invTransform = Matrix4x4(stream);
        }

        /** \brief Create a transform from the given matrix
         * and calculate the inverse
         */
        Transform(const Matrix4x4 &trafo)
                : m_transform(trafo) {
                bool success = m_transform.invert(m_invTransform);
                if (!success)
                        SLog(EError, "Unable to invert singular matrix %s", trafo.toString().c_str());
        }

        /// Create a transform from the given matrices
        Transform(const Matrix4x4 &trafo, const Matrix4x4 &invTrafo)
                : m_transform(trafo), m_invTransform(invTrafo) {
        }

        /// Return the inverse transform
        Transform inverse() const {
                return Transform(m_invTransform, m_transform);
        }

        /// Matrix-matrix multiplication
        Transform operator*(const Transform &t) const;

        /// Return the determinant of the upper left 3x3 submatrix
        inline Float det3x3() const {
                return m_transform.det3x3();
        }

        /// Test for a scale component
        inline bool hasScale() const {
                for (int i=0; i<3; ++i) {
                        for (int j=i; j<3; ++j) {
                                Float sum = 0;
                                for (int k=0; k<3; ++k)
                                        sum += m_transform.m[i][k] * m_transform.m[j][k];

                                if (i == j && std::abs(sum-1) > 1e-3f)
                                        return true;
                                else if (i != j && std::abs(sum) > 1e-3f)
                                        return true;
                        }
                }
                return false;
        }

        /// Test if this is the identity transformation
        inline bool isIdentity() const {
                return m_transform.isIdentity();
        }

        /// Test if this is the zero transformation
        inline bool isZero() const {
                return m_transform.isZero();
        }

        /**
         * \brief Matrix-vector multiplication for points in 3d space
         *
         * \remark In the Python bindings, this is function implemented as
         * the multiplication operator (\c __mul__).
         */
        inline Point operator()(const Point &p) const {
                Float x = m_transform.m[0][0] * p.x + m_transform.m[0][1] * p.y
                        + m_transform.m[0][2] * p.z + m_transform.m[0][3];
                Float y = m_transform.m[1][0] * p.x + m_transform.m[1][1] * p.y
                        + m_transform.m[1][2] * p.z + m_transform.m[1][3];
                Float z = m_transform.m[2][0] * p.x + m_transform.m[2][1] * p.y
                        + m_transform.m[2][2] * p.z + m_transform.m[2][3];
                Float w = m_transform.m[3][0] * p.x + m_transform.m[3][1] * p.y
                        + m_transform.m[3][2] * p.z + m_transform.m[3][3];
#ifdef MTS_DEBUG
                if (w == 0)
                        SLog(EWarn, "w==0 in Transform::operator(Point &)");
#endif
                if (w == 1.0f)
                        return Point(x, y, z);
                else
                        return Point(x, y, z) / w;
        }

        /// Transform a point by an affine / non-projective matrix
        inline Point transformAffine(const Point &p) const {
                Float x = m_transform.m[0][0] * p.x + m_transform.m[0][1] * p.y
                        + m_transform.m[0][2] * p.z + m_transform.m[0][3];
                Float y = m_transform.m[1][0] * p.x + m_transform.m[1][1] * p.y
                        + m_transform.m[1][2] * p.z + m_transform.m[1][3];
                Float z = m_transform.m[2][0] * p.x + m_transform.m[2][1] * p.y
                        + m_transform.m[2][2] * p.z + m_transform.m[2][3];
                return Point(x,y,z);
        }

        /// Transform a point by a affine / non-projective matrix (no temporaries)
        inline void transformAffine(const Point &p, Point &dest) const {
                dest.x = m_transform.m[0][0] * p.x + m_transform.m[0][1] * p.y
                       + m_transform.m[0][2] * p.z + m_transform.m[0][3];
                dest.y = m_transform.m[1][0] * p.x + m_transform.m[1][1] * p.y
                       + m_transform.m[1][2] * p.z + m_transform.m[1][3];
                dest.z = m_transform.m[2][0] * p.x + m_transform.m[2][1] * p.y
                       + m_transform.m[2][2] * p.z + m_transform.m[2][3];
        }

        /**
         * \brief Matrix-vector multiplication for points in 3d space (no temporaries)
         * \remark This function is not available in the Python bindings
         */
    inline void operator()(const Point &p, Point &dest) const {
                dest.x = m_transform.m[0][0] * p.x + m_transform.m[0][1] * p.y
                       + m_transform.m[0][2] * p.z + m_transform.m[0][3];
                dest.y = m_transform.m[1][0] * p.x + m_transform.m[1][1] * p.y
                       + m_transform.m[1][2] * p.z + m_transform.m[1][3];
                dest.z = m_transform.m[2][0] * p.x + m_transform.m[2][1] * p.y
                       + m_transform.m[2][2] * p.z + m_transform.m[2][3];
                Float w = m_transform.m[3][0] * p.x + m_transform.m[3][1] * p.y
                        + m_transform.m[3][2] * p.z + m_transform.m[3][3];

#ifdef MTS_DEBUG
                if (w == 0)
                        SLog(EWarn, "w==0 in Transform::operator(Point &, Point &)");
#endif
                if (w != 1.0f)
                        dest /= w;
        }

        /**
         * \brief Matrix-vector multiplication for vectors in 3d space
         * \remark In the Python bindings, this is function implemented as
         * the multiplication operator (\c __mul__).
         */
    inline Vector operator()(const Vector &v) const {
                Float x = m_transform.m[0][0] * v.x + m_transform.m[0][1] * v.y
                        + m_transform.m[0][2] * v.z;
                Float y = m_transform.m[1][0] * v.x + m_transform.m[1][1] * v.y
                        + m_transform.m[1][2] * v.z;
                Float z = m_transform.m[2][0] * v.x + m_transform.m[2][1] * v.y
                        + m_transform.m[2][2] * v.z;
                return Vector(x, y, z);
        }

        /**
         * \brief Matrix-vector multiplication for vectors in 3d space (no temporaries)
         * \remark This function is not available in the Python bindings
         */
    inline void operator()(const Vector &v, Vector &dest) const {
                dest.x = m_transform.m[0][0] * v.x + m_transform.m[0][1] * v.y
                       + m_transform.m[0][2] * v.z;
                dest.y = m_transform.m[1][0] * v.x + m_transform.m[1][1] * v.y
                       + m_transform.m[1][2] * v.z;
                dest.z = m_transform.m[2][0] * v.x + m_transform.m[2][1] * v.y
                       + m_transform.m[2][2] * v.z;
        }

        /**
         * \brief Matrix-normal multiplication
         * \remark In the Python bindings, this is function implemented as
         * the multiplication operator (\c __mul__).
         */
    inline Normal operator()(const Normal &v) const {
                Float x = m_invTransform.m[0][0] * v.x + m_invTransform.m[1][0] * v.y
                        + m_invTransform.m[2][0] * v.z;
                Float y = m_invTransform.m[0][1] * v.x + m_invTransform.m[1][1] * v.y
                        + m_invTransform.m[2][1] * v.z;
                Float z = m_invTransform.m[0][2] * v.x + m_invTransform.m[1][2] * v.y
                        + m_invTransform.m[2][2] * v.z;
                return Normal(x, y, z);
        }

        /**
         * \brief Matrix-normal multiplication (no temporaries)
         * \remark This function is not available in the Python bindings
         */
    inline void operator()(const Normal &v, Normal &dest) const {
                dest.x = m_invTransform.m[0][0] * v.x + m_invTransform.m[1][0] * v.y
                       + m_invTransform.m[2][0] * v.z;
                dest.y = m_invTransform.m[0][1] * v.x + m_invTransform.m[1][1] * v.y
                       + m_invTransform.m[2][1] * v.z;
                dest.z = m_invTransform.m[0][2] * v.x + m_invTransform.m[1][2] * v.y
                       + m_invTransform.m[2][2] * v.z;
        }

        /**
         * \brief 4D matrix-vector multiplication
         * \remark In the Python bindings, this is function implemented as
         * the multiplication operator (\c __mul__).
         */
        inline Vector4 operator()(const Vector4 &v) const {
                Float x = m_transform.m[0][0] * v.x + m_transform.m[0][1] * v.y
                        + m_transform.m[0][2] * v.z + m_transform.m[0][3] * v.w;
                Float y = m_transform.m[1][0] * v.x + m_transform.m[1][1] * v.y
                        + m_transform.m[1][2] * v.z + m_transform.m[1][3] * v.w;
                Float z = m_transform.m[2][0] * v.x + m_transform.m[2][1] * v.y
                        + m_transform.m[2][2] * v.z + m_transform.m[2][3] * v.w;
                Float w = m_transform.m[3][0] * v.x + m_transform.m[3][1] * v.y
                        + m_transform.m[3][2] * v.z + m_transform.m[3][3] * v.w;
                return Vector4(x,y,z,w);
        }

        /**
         * \brief 4D matrix-vector multiplication (no temporaries)
         * \remark This function is not available in the Python bindings
         */
        inline void operator()(const Vector4 &v, Vector4 &dest) const {
                dest.x = m_transform.m[0][0] * v.x + m_transform.m[0][1] * v.y
                       + m_transform.m[0][2] * v.z + m_transform.m[0][3] * v.w;
                dest.y = m_transform.m[1][0] * v.x + m_transform.m[1][1] * v.y
                       + m_transform.m[1][2] * v.z + m_transform.m[1][3] * v.w;
                dest.z = m_transform.m[2][0] * v.x + m_transform.m[2][1] * v.y
                       + m_transform.m[2][2] * v.z + m_transform.m[2][3] * v.w;
                dest.w = m_transform.m[3][0] * v.x + m_transform.m[3][1] * v.y
                       + m_transform.m[3][2] * v.z + m_transform.m[3][3] * v.w;
        }

        /**
         * \brief Transform a ray. Assumes that there is no scaling (no temporaries)
         * \remark This function is not available in the Python bindings
         */
        inline void operator()(const Ray &a, Ray &b) const {
                b.mint = a.mint;
                b.maxt = a.maxt;
                operator()(a.o, b.o);
                operator()(a.d, b.d);
#ifdef MTS_DEBUG_FP
                bool state = disableFPExceptions();
#endif
                /* Re-compute the reciprocal */
                b.dRcp.x = 1.0f / b.d.x;
                b.dRcp.y = 1.0f / b.d.y;
                b.dRcp.z = 1.0f / b.d.z;
#ifdef MTS_DEBUG_FP
                restoreFPExceptions(state);
#endif
                b.time = a.time;
        }

        /**
         * \brief Transform a ray
         * \remark In the Python bindings, this is function implemented as
         * the multiplication operator (\c __mul__).
         */
        inline Ray operator()(const Ray &ray) const {
                Ray result;
                operator()(ray, result);
                return result;
        }

        /// Transform a ray by an affine / non-projective matrix (no temporaries)
        inline void transformAffine(const Ray &a, Ray &b) const {
                b.mint = a.mint;
                b.maxt = a.maxt;
                transformAffine(a.o, b.o);
                operator()(a.d, b.d);
#ifdef MTS_DEBUG_FP
                bool state = disableFPExceptions();
#endif
                /* Re-compute the reciprocal */
                b.dRcp.x = 1.0f / b.d.x;
                b.dRcp.y = 1.0f / b.d.y;
                b.dRcp.z = 1.0f / b.d.z;
#ifdef MTS_DEBUG_FP
                restoreFPExceptions(state);
#endif
        }

        /// Transform a ray by an affine / non-projective matrix (no temporaries)
        inline Ray transformAffine(const Ray &ray) const {
                Ray result;
                transformAffine(ray, result);
                return result;
        }

        /// Return the underlying matrix
        inline const Matrix4x4 &getMatrix() const { return m_transform; }

        /// Return the underlying inverse matrix (const version)
        inline const Matrix4x4 &getInverseMatrix() const { return m_invTransform; }

        /// Create a translation transformation
        static Transform translate(const Vector &v);

        /// Create a rotation transformation around an arbitrary axis. The angle is specified in degrees
        static Transform rotate(const Vector &axis, Float angle);

        /// Create a scale transformation
        static Transform scale(const Vector &v);

        /** \brief Create a perspective transformation.
         *   (Maps [near, far] to [0, 1])
         * \param fov Field of view in degrees
         * \param clipNear Near clipping plane
         * \param clipFar Far clipping plane
         */
        static Transform perspective(Float fov, Float clipNear, Float clipFar);

        /** \brief Create a perspective transformation for OpenGL.
         *   (Maps [-near, -far] to [-1, 1])
         * \param fov Field of view in degrees
         * \param clipNear Near clipping plane distance
         * \param clipFar Far clipping plane distance
         */
        static Transform glPerspective(Float fov, Float clipNear, Float clipFar);

        /** \brief Create a perspective transformation for OpenGL.
         * \param left Left clipping plane coordinate
         * \param right Right clipping plane coordinate
         * \param top Top clipping plane coordinate
         * \param bottom Bottom clipping plane coordinate
         * \param nearVal Near clipping plane distance
         * \param farVal Far clipping plane distance
         */
        static Transform glFrustum(Float left, Float right, Float bottom, Float top, Float nearVal, Float farVal);

        /** \brief Create an orthographic transformation, which maps Z to [0,1]
         * and leaves the X and Y coordinates untouched.
         * \param clipNear Near clipping plane
         * \param clipFar Far clipping plane
         */
        static Transform orthographic(Float clipNear, Float clipFar);

        /** \brief Create an orthographic transformation for OpenGL
         * \param clipNear Near clipping plane
         * \param clipFar Far clipping plane
         */
        static Transform glOrthographic(Float clipNear, Float clipFar);

        /** \brief Create an orthographic transformation for OpenGL
         *
         * Slightly extended variant which also handles non-unity clipping
         * planes on the X and Y axes and matches the 'glOrtho' spec.
         */
        static Transform glOrthographic(Float clipLeft, Float clipRight,
                Float clipBottom, Float clipTop, Float clipNear, Float clipFar);

        /** \brief Create a look-at camera transformation
         * \param p Camera position
         * \param t Target vector
         * \param u Up vector
         */
        static Transform lookAt(const Point &p, const Point &t, const Vector &u);

        /** \brief Create an orthogonal transformation that takes
         * the standard to the supplied frame
         */
        static Transform fromFrame(const Frame &frame);

        /// Serialize a transformation to a stream
        inline void serialize(Stream *stream) const {
                m_transform.serialize(stream);
                m_invTransform.serialize(stream);
        }

        /// Equality comparison operator
        inline bool operator==(const Transform &trafo) const {
                return m_transform == trafo.m_transform;
        }

        /// Inequality comparison operator
        inline bool operator!=(const Transform &trafo) const {
                return m_transform != trafo.m_transform;
        }

        /// Return a string representation
        std::string toString() const;
private:
        Matrix4x4 m_transform;
        Matrix4x4 m_invTransform;
};

MTS_NAMESPACE_END

#endif /* __MITSUBA_CORE_TRANSFORM_H_ */