source: GTP/trunk/App/Demos/Vis/CHC_revisited/Vector3.h @ 2752

#ifndef _Vector3_h__
#define _Vector3_h__

#include <iostream>
#include "common.h"
#include <math.h>


// Forward-declare some other classes.
class Matrix4x4;


namespace CHCDemo
{

/** Vector class.
*/
class Vector3
{
public:

        ///////////
        //-- members;

        float x; 
        float y;
        float z;


        /////////////////////////

        void  SetX(float q) { x = q; }
        void  SetY(float q) { y = q; }
        void  SetZ(float q) { z = q; }

        float GetX() const { return x; }
        float GetY() const { return y; }
        float GetZ() const { return z; }


        //////////////
        //-- constructors

        Vector3() { }

        Vector3(float x, float y, float z):
        x(x), y(y), z(z)
        {
        }
        Vector3(float v):
        x(v), y(v), z(v)
        {
        }
        /** Copy constructor.
        */
        Vector3(const Vector3 &v): x(v.x), y(v.y), z(v.z)
        {}

        // Functions to get at the std::vector components
        float& operator[] (const int inx) 
        {
                return (&x)[inx];
        }

        operator const float*() const 
        { 
                return (const float*) this; 
        }

        const float &operator[] (const int inx) const 
        {
                return *(&x + inx);
        }

        void ExtractVerts(float *px, float *py, int which) const;
  
        void SetValue(float a, float b, float c)
        {  
                x = a; y = b; z = c; 
        }

        void SetValue(float a) 
        { 
                x = y = z = a; 
        }
 
        /** Returns the axis, where the std::vector has the largest value.
        */
        int DrivingAxis(void) const;
        /** returns the axis, where the std::vector has the smallest value
        */
        int TinyAxis(void) const;
        /** Returns largest component in this vector
        */
        inline float MaxComponent(void) const 
        {
                return (x > y) ? ( (x > z) ? x : z) : ( (y > z) ? y : z);
        }
        /** Returns copy of this vector where all components are positiv.
        */
        inline Vector3 Abs(void) const 
        {
                return Vector3(fabs(x), fabs(y), fabs(z));
        }
        /** normalizes the std::vector of unit size corresponding to given std::vector
        */
        inline void Normalize();
        /** Returns false if this std::vector has a nan component.
        */
        bool CheckValidity() const;
        /**
        Find a right handed coordinate system with (*this) being
        the z-axis. For a right-handed system, U x V = (*this) holds.
        This implementation is here to avoid inconsistence and confusion
        when construction coordinate systems using ArbitraryNormal():
        In fact:
        V = ArbitraryNormal(N);
        U = CrossProd(V,N);
        constructs a right-handed coordinate system as well, BUT:
        1) bugs can be introduced if one mistakenly constructs a 
        left handed sytems e.g. by doing
        U = ArbitraryNormal(N);
        V = CrossProd(U,N);
        2) this implementation gives non-negative base vectors
        for (*this)==(0,0,1) |  (0,1,0) | (1,0,0), which is
        good for debugging and is not the case with the implementation
        using ArbitraryNormal()

        ===> Using ArbitraryNormal() for constructing coord systems 
        is obsoleted by this method (<JK> 12/20/03).   
        */
        void RightHandedBase(Vector3& U, Vector3& V) const;
        
        
        // Unary operators

        Vector3 operator+ () const;
        Vector3 operator- () const;

        // Assignment operators
        
        Vector3& operator+= (const Vector3 &A);
        Vector3& operator-= (const Vector3 &A);
        Vector3& operator*= (const Vector3 &A);
        Vector3& operator*= (float A);
        Vector3& operator/= (float A);


        //////////
        //-- friends
        
        /**
        ===> Using ArbitraryNormal() for constructing coord systems 
        ===> is obsoleted by RightHandedBase() method (<JK> 12/20/03).   

        Return an arbitrary normal to `v'.
        In fact it tries v x (0,0,1) an if the result is too small,
        it definitely does v x (0,1,0). It will always work for
        non-degenareted std::vector and is much faster than to use
        TangentVectors.

        @param v(in) The std::vector we want to find normal for.
        @return The normal std::vector to v.
        */
        friend inline Vector3 ArbitraryNormal(const Vector3 &v); 

        /// Transforms a std::vector to the global coordinate frame.
        /**
        Given a local coordinate frame (U,V,N) (i.e. U,V,N are 
        the x,y,z axes of the local coordinate system) and
        a std::vector 'loc' in the local coordiante system, this
        function returns a the coordinates of the same std::vector
        in global frame (i.e. frame (1,0,0), (0,1,0), (0,0,1).
        */
        friend inline Vector3 ToGlobalFrame(const Vector3& loc,
                                                const Vector3& U,
                                                                                const Vector3& V,
                                                                                const Vector3& N);

        /// Transforms a std::vector to a local coordinate frame.
        /**
        Given a local coordinate frame (U,V,N) (i.e. U,V,N are 
        the x,y,z axes of the local coordinate system) and
        a std::vector 'loc' in the global coordiante system, this
        function returns a the coordinates of the same std::vector
        in the local frame.
        */
        friend inline Vector3 ToLocalFrame(const Vector3& loc,
                                           const Vector3& U,
                                                                           const Vector3& V,
                                                                           const Vector3& N);

        /// the magnitude=size of the std::vector
        friend inline float Magnitude(const Vector3 &v);
        /// the squared magnitude of the std::vector .. for efficiency in some cases
        friend inline float SqrMagnitude(const Vector3 &v);
        /// Magnitude(v1-v2)
        friend inline float Distance(const Vector3 &v1, const Vector3 &v2);
        /// SqrMagnitude(v1-v2)
        friend inline float SqrDistance(const Vector3 &v1, const Vector3 &v2);
        /// creates the std::vector of unit size corresponding to given std::vector
        friend inline Vector3 Normalize(const Vector3 &A);
        /// // Rotate a direction vector
        friend Vector3 PlaneRotate(const Matrix4x4 &, const Vector3 &);

        // construct view vectors .. DirAt is the main viewing direction
        // Viewer is the coordinates of viewer location, UpL is the std::vector.
        friend void ViewVectors(const Vector3 &DirAt, 
                                                        const Vector3 &Viewer,
                                                        const Vector3 &UpL, 
                                                        Vector3 &ViewV,
                                                        Vector3 &ViewU, 
                                                        Vector3 &ViewN);

        // Given the intersection point `P', you have available normal `N'
        // of unit length. Let us suppose the incoming ray has direction `D'.
        // Then we can construct such two vectors `U' and `V' that
        // `U',`N', and `D' are coplanar, and `V' is perpendicular
        // to the vectors `N','D', and `V'. Then 'N', 'U', and 'V' create
        // the orthonormal base in space R3.
        friend void TangentVectors(Vector3 &U, Vector3 &V, // output
                                                           const Vector3 &normal, // input
                                                           const Vector3 &dirIncoming);


        // Binary operators
        
        friend inline Vector3 operator+ (const Vector3 &A, const Vector3 &B);
        friend inline Vector3 operator- (const Vector3 &A, const Vector3 &B);
        friend inline Vector3 operator* (const Vector3 &A, const Vector3 &B);
        friend inline Vector3 operator* (const Vector3 &A, float B);
        friend inline Vector3 operator* (float A, const Vector3 &B);
        friend Vector3 operator* (const Matrix4x4 &, const Vector3 &);
        friend inline Vector3 operator/ (const Vector3 &A, const Vector3 &B);

        friend inline int operator< (const Vector3 &A, const Vector3 &B);
        friend inline int operator<= (const Vector3 &A, const Vector3 &B);

        friend inline Vector3 operator/ (const Vector3 &A, float B);
        friend inline int operator== (const Vector3 &A, const Vector3 &B);
        friend inline float DotProd(const Vector3 &A, const Vector3 &B);
        friend inline Vector3 CrossProd (const Vector3 &A, const Vector3 &B);

        friend std::ostream& operator<< (std::ostream &s, const Vector3 &A);
        friend std::istream& operator>> (std::istream &s, Vector3 &A);

        friend void Minimize(Vector3 &min, const Vector3 &candidate);
        friend void Maximize(Vector3 &max, const Vector3 &candidate);

        friend inline int EpsilonEqualV3(const Vector3 &v1, const Vector3 &v2, float thr);
        friend inline int EpsilonEqualV3(const Vector3 &v1, const Vector3 &v2);

        friend Vector3 CosineRandomVector(const Vector3 &normal);
        friend Vector3 UniformRandomVector(const Vector3 &normal);
        friend Vector3 UniformRandomVector();

        friend Vector3 UniformRandomVector(float r1, float r2);

        friend Vector3 CosineRandomVector(float r1, float r2, const Vector3 &normal);
};



// forward declaration
extern Vector3 UniformRandomVector(const Vector3 &normal);
extern Vector3 UniformRandomVector();
extern Vector3 UniformRandomVector(const float r1, const float r2);

inline Vector3 ArbitraryNormal(const Vector3 &N)
{
        float dist2 = N.x * N.x + N.y * N.y;

        if (dist2 > 0.0001f) 
        {
                float inv_size = 1.0f / sqrtf(dist2);
                return Vector3(N.y * inv_size, -N.x * inv_size, 0); // N x (0,0,1)
        }

        float inv_size = 1.0f / sqrtf(N.z * N.z + N.x * N.x);
        return Vector3(-N.z * inv_size, 0, N.x * inv_size); // N x (0,1,0)
}


inline void Vector3::RightHandedBase(Vector3& U, Vector3& V) const
{
        // HACK
        V = ArbitraryNormal(*this);
        U = CrossProd(V, *this);
}


inline Vector3 ToGlobalFrame(const Vector3 &loc,
                                                         const Vector3 &U,
                                                         const Vector3 &V,
                                                         const Vector3 &N)
{
        return loc.x * U + loc.y * V + loc.z * N;
}


inline Vector3 ToLocalFrame(const Vector3 &loc,
                                                        const Vector3 &U,
                                                        const Vector3 &V,
                                                        const Vector3 &N)
{
        return Vector3(loc.x * U.x + loc.y * U.y + loc.z * U.z,
                                   loc.x * V.x + loc.y * V.y + loc.z * V.z,
                                   loc.x * N.x + loc.y * N.y + loc.z * N.z);
}


inline float Magnitude(const Vector3 &v) 
{
        return sqrtf(v.x * v.x + v.y * v.y + v.z * v.z);
}


inline float SqrMagnitude(const Vector3 &v) 
{
  return v.x * v.x + v.y * v.y + v.z * v.z;
}


inline float Distance(const Vector3 &v1, const Vector3 &v2)
{
  return sqrtf(sqrt(v1.x - v2.x) + sqrt(v1.y - v2.y) + sqrt(v1.z - v2.z));
}


inline float SqrDistance(const Vector3 &v1, const Vector3 &v2)
{
  return sqrt(v1.x - v2.x) + sqrt(v1.y - v2.y) + sqrt(v1.z - v2.z);
}


inline Vector3 Normalize(const Vector3 &A)
{
  return A * (1.0f / Magnitude(A));
}


inline float DotProd(const Vector3 &A, const Vector3 &B)
{
  return A.x * B.x + A.y * B.y + A.z * B.z;
}


// angle between two vectors with respect to a surface normal in the
// range [0 .. 2 * pi]
inline float Angle(const Vector3 &A, const Vector3 &B, const Vector3 &norm)
{
        Vector3 cross = CrossProd(A, B);

        float signedAngle;

        if (DotProd(cross, norm) > 0)
                signedAngle = atan2(-Magnitude(CrossProd(A, B)), DotProd(A, B));
        else
                signedAngle = atan2(Magnitude(CrossProd(A, B)), DotProd(A, B));

        if (signedAngle < 0)
                return 2 * M_PI + signedAngle;

        return signedAngle;
}


inline Vector3 Vector3::operator+() const
{
        return *this;
}


inline Vector3 Vector3::operator-() const
{
  return Vector3(-x, -y, -z);
}


inline Vector3 &Vector3::operator+=(const Vector3 &A)
{
        x += A.x;  y += A.y;  z += A.z;
        return *this;
}


inline Vector3& Vector3::operator-=(const Vector3 &A)
{
        x -= A.x;  y -= A.y;  z -= A.z;
        
        return *this;
}


inline Vector3& Vector3::operator*= (float A)
{
        x *= A;  y *= A;  z *= A;
        return *this;
}


inline Vector3& Vector3::operator/=(float A)
{
        float a = 1.0f / A;
        
        x *= a;  y *= a;  z *= a;
        
        return *this;
}


inline Vector3& Vector3::operator*= (const Vector3 &A)
{
  x *= A.x;  y *= A.y;  z *= A.z;
  return *this;
}


inline Vector3 operator+ (const Vector3 &A, const Vector3 &B)
{
        return Vector3(A.x + B.x, A.y + B.y, A.z + B.z);
}


inline Vector3 operator- (const Vector3 &A, const Vector3 &B)
{
        return Vector3(A.x - B.x, A.y - B.y, A.z - B.z);
}


inline Vector3 operator* (const Vector3 &A, const Vector3 &B)
{
        return Vector3(A.x * B.x, A.y * B.y, A.z * B.z);
}


inline Vector3 operator* (const Vector3 &A, float B)
{
  return Vector3(A.x * B, A.y * B, A.z * B);
}


inline Vector3 operator* (float A, const Vector3 &B)
{
        return Vector3(B.x * A, B.y * A, B.z * A);
}


inline Vector3 operator/ (const Vector3 &A, const Vector3 &B)
{
        return Vector3(A.x / B.x, A.y / B.y, A.z / B.z);
}


inline Vector3 operator/ (const Vector3 &A, float B)
{
        float b = 1.0f / B;
        return Vector3(A.x * b, A.y * b, A.z * b);
}


inline int operator< (const Vector3 &A, const Vector3 &B)
{
        return A.x < B.x && A.y < B.y && A.z < B.z;
}


inline int operator<= (const Vector3 &A, const Vector3 &B)
{
        return A.x <= B.x && A.y <= B.y && A.z <= B.z;
}


// Might replace floating-point == with comparisons of
// magnitudes, if needed.
inline int operator== (const Vector3 &A, const Vector3 &B)
{ 
        return (A.x == B.x) && (A.y == B.y) && (A.z == B.z);
}


inline Vector3 CrossProd (const Vector3 &A, const Vector3 &B)
{
        return Vector3(A.y * B.z - A.z * B.y,
                           A.z * B.x - A.x * B.z,
                                   A.x * B.y - A.y * B.x);
}


inline void Vector3::Normalize()
{
        float sqrmag = x * x + y * y + z * z;
        
        if (sqrmag > 0.0f)
                (*this) *= 1.0f / sqrtf(sqrmag);
}


// Overload << operator for C++-style output
inline std::ostream& operator<< (std::ostream &s, const Vector3 &A)
{
        return s << "(" << A.x << ", " << A.y << ", " << A.z << ")";
}


// Overload >> operator for C++-style input
inline std::istream& operator>> (std::istream &s, Vector3 &A)
{
        char a;
        // read "(x, y, z)"
        return s >> a >> A.x >> a >> A.y >> a >> A.z >> a;
}


inline int EpsilonEqualV3(const Vector3 &v1, const Vector3 &v2, float thr)
{
        if (fabsf(v1.x-v2.x) > thr)
                return false;
        if (fabsf(v1.y-v2.y) > thr)
                return false;
        if (fabsf(v1.z-v2.z) > thr)
                return false;
        
        return true;
}


inline int EpsilonEqualV3(const Vector3 &v1, const Vector3 &v2)
{
        return EpsilonEqualV3(v1, v2, Limits::Small);
}


}

#endif