1 #ifndef __RIMATH_H

     2 #define __RIMATH_H

     3 

     4 /*------------------------------------------------------------------------

     5  *

     6  * OpenVG 1.1 Reference Implementation

     7  * -----------------------------------

     8  *

     9  * Copyright (c) 2007 The Khronos Group Inc.

    10  *

    11  * Permission is hereby granted, free of charge, to any person obtaining a

    12  * copy of this software and /or associated documentation files

    13  * (the "Materials "), to deal in the Materials without restriction,

    14  * including without limitation the rights to use, copy, modify, merge,

    15  * publish, distribute, sublicense, and/or sell copies of the Materials,

    16  * and to permit persons to whom the Materials are furnished to do so,

    17  * subject to the following conditions: 

    18  *

    19  * The above copyright notice and this permission notice shall be included 

    20  * in all copies or substantial portions of the Materials. 

    21  *

    22  * THE MATERIALS ARE PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,

    23  * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF

    24  * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.

    25  * IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,

    26  * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR

    27  * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE MATERIALS OR

    28  * THE USE OR OTHER DEALINGS IN THE MATERIALS.

    29  *

    30  *//**

    31  * \file

    32  * \brief	Math functions, Vector and Matrix classes.

    33  * \note	

    34  *//*-------------------------------------------------------------------*/

    35 

    36 #ifndef __RIDEFS_H

    37 #include "riDefs.h"

    38 #endif

    39 

    40 #include <math.h>

    41 

    42 namespace OpenVGRI

    43 {

    44 

    45 /*-------------------------------------------------------------------*//*!

    46 * \brief	

    47 * \param	

    48 * \return	

    49 * \note		

    50 *//*-------------------------------------------------------------------*/

    51 

    52 RI_INLINE int		RI_ISNAN(float a)

    53 {

    54 	RIfloatInt p;

    55 	p.f = a;

    56 	unsigned int exponent = (p.i>>23) & 0xff;

    57 	unsigned int mantissa = p.i & 0x7fffff;

    58 	if(exponent == 255 && mantissa)

    59 		return 1;

    60 	return 0;

    61 }

    62 

    63 #if (RI_MANTISSA_BITS > 23)

    64 #error RI_MANTISSA_BITS is greater than 23

    65 #elif (RI_EXPONENT_BITS > 8)

    66 #error RI_EXPONENT_BITS is greater than 8

    67 #elif (RI_MANTISSA_BITS != 23) || (RI_EXPONENT_BITS != 8)

    68 

    69 class RIfloat

    70 {

    71 public:

    72 	RIfloat() : v(0.0f)						{ removeBits(); }

    73 	RIfloat(float a) : v(a)					{ removeBits(); }

    74 	RIfloat(double a) : v((float)a)			{ removeBits(); }

    75 	RIfloat(int a) : v((float)a)			{ removeBits(); }

    76 	RIfloat(unsigned int a) : v((float)a)	{ removeBits(); }

    77 	RIfloat&	operator=(const RIfloat &a)	{ v = a.v; removeBits(); return *this; }

    78 	RIfloat&	operator+=(const RIfloat &a){ v += a.v; removeBits(); return *this; }

    79 	RIfloat&	operator-=(const RIfloat &a){ v -= a.v; removeBits(); return *this; }

    80 	RIfloat&	operator*=(const RIfloat &a){ v *= a.v; removeBits(); return *this; }

    81 	RIfloat&	operator/=(const RIfloat &a){ v /= a.v; removeBits(); return *this; }

    82 	RIfloat		operator-() const			{ return -v; }

    83 	operator float() const					{ return v; }

    84 	operator double() const					{ return (double)v; }

    85 	operator int() const					{ return (int)v; }

    86 

    87 	friend RIfloat	operator+(const RIfloat &a, const RIfloat &b);

    88 	friend RIfloat	operator+(float a, const RIfloat &b);

    89 	friend RIfloat	operator+(const RIfloat &a, float b);

    90 	friend RIfloat	operator-(const RIfloat &a, const RIfloat &b);

    91 	friend RIfloat	operator-(float a, const RIfloat &b);

    92 	friend RIfloat	operator-(const RIfloat &a, float b);

    93 	friend RIfloat	operator*(const RIfloat &a, const RIfloat &b);

    94 	friend RIfloat	operator*(float a, const RIfloat &b);

    95 	friend RIfloat	operator*(const RIfloat &a, float b);

    96 	friend RIfloat	operator/(const RIfloat &a, const RIfloat &b);

    97 	friend RIfloat	operator/(float a, const RIfloat &b);

    98 	friend RIfloat	operator/(const RIfloat &a, float b);

    99 

   100 	friend bool		operator<(const RIfloat &a, const RIfloat &b);

   101 	friend bool		operator<(float a, const RIfloat &b);

   102 	friend bool		operator<(const RIfloat &a, float b);

   103 	friend bool		operator>(const RIfloat &a, const RIfloat &b);

   104 	friend bool		operator>(float a, const RIfloat &b);

   105 	friend bool		operator>(const RIfloat &a, float b);

   106 	friend bool		operator<=(const RIfloat &a, const RIfloat &b);

   107 	friend bool		operator<=(float a, const RIfloat &b);

   108 	friend bool		operator<=(const RIfloat &a, float b);

   109 	friend bool		operator>=(const RIfloat &a, const RIfloat &b);

   110 	friend bool		operator>=(float a, const RIfloat &b);

   111 	friend bool		operator>=(const RIfloat &a, float b);

   112 	friend bool		operator==(const RIfloat &a, const RIfloat &b);

   113 	friend bool		operator==(float a, const RIfloat &b);

   114 	friend bool		operator==(const RIfloat &a, float b);

   115 	friend bool		operator!=(const RIfloat &a, const RIfloat &b);

   116 	friend bool		operator!=(float a, const RIfloat &b);

   117 	friend bool		operator!=(const RIfloat &a, float b);

   118 private:

   119 	void	removeBits()

   120 	{

   121 		RIfloatInt p;

   122 		p.f = v;

   123 		unsigned int exponent = (p.i>>23) & 0xff;

   124 		if(exponent == 0 || exponent == 255)

   125 			return;	//zero, denormal, infinite, or NaN

   126 

   127 		p.i &= ~((1<<(23-RI_MANTISSA_BITS))-1);

   128 

   129 #if (RI_EXPONENT_BITS != 8)

   130 		if (exponent > 127 + (1 << (RI_EXPONENT_BITS-1)))

   131 			exponent = 127 + (1 << (RI_EXPONENT_BITS-1));

   132 

   133 		if (exponent < 127 + 1 - (1 << (RI_EXPONENT_BITS-1)))

   134 			exponent = 127 + 1 - (1 << (RI_EXPONENT_BITS-1));

   135 

   136 		p.i &= ~(0xff<<23);

   137 		p.i |= exponent<<23;

   138 #endif

   139 		v = p.f;

   140 	}

   141 

   142 	float		v;

   143 };

   144 

   145 RI_INLINE RIfloat operator+(const RIfloat &a, const RIfloat &b)	{ return RIfloat(a.v+b.v); }

   146 RI_INLINE RIfloat operator+(float a, const RIfloat &b)			{ return RIfloat(a+b.v); }

   147 RI_INLINE RIfloat operator+(const RIfloat &a, float b)			{ return RIfloat(a.v+b); }

   148 RI_INLINE RIfloat operator-(const RIfloat &a, const RIfloat &b)	{ return RIfloat(a.v-b.v); }

   149 RI_INLINE RIfloat operator-(float a, const RIfloat &b)			{ return RIfloat(a-b.v); }

   150 RI_INLINE RIfloat operator-(const RIfloat &a, float b)			{ return RIfloat(a.v-b); }

   151 RI_INLINE RIfloat operator*(const RIfloat &a, const RIfloat &b)	{ return RIfloat(a.v*b.v); }

   152 RI_INLINE RIfloat operator*(float a, const RIfloat &b)			{ return RIfloat(a*b.v); }

   153 RI_INLINE RIfloat operator*(const RIfloat &a, float b)			{ return RIfloat(a.v*b); }

   154 RI_INLINE RIfloat operator/(const RIfloat &a, const RIfloat &b)	{ return RIfloat(a.v/b.v); }

   155 RI_INLINE RIfloat operator/(float a, const RIfloat &b)			{ return RIfloat(a/b.v); }

   156 RI_INLINE RIfloat operator/(const RIfloat &a, float b)			{ return RIfloat(a.v/b); }

   157 

   158 RI_INLINE bool operator<(const RIfloat &a, const RIfloat &b)	{ return a.v < b.v ? true : false; }

   159 RI_INLINE bool operator<(float a, const RIfloat &b)				{ return a < b.v ? true : false; }

   160 RI_INLINE bool operator<(const RIfloat &a, float b)				{ return a.v < b ? true : false; }

   161 RI_INLINE bool operator>(const RIfloat &a, const RIfloat &b)	{ return a.v > b.v ? true : false; }

   162 RI_INLINE bool operator>(float a, const RIfloat &b)				{ return a > b.v ? true : false; }

   163 RI_INLINE bool operator>(const RIfloat &a, float b)				{ return a.v > b ? true : false; }

   164 RI_INLINE bool operator<=(const RIfloat &a, const RIfloat &b)	{ return a.v <= b.v ? true : false; }

   165 RI_INLINE bool operator<=(float a, const RIfloat &b)			{ return a <= b.v ? true : false; }

   166 RI_INLINE bool operator<=(const RIfloat &a, float b)			{ return a.v <= b ? true : false; }

   167 RI_INLINE bool operator>=(const RIfloat &a, const RIfloat &b)	{ return a.v >= b.v ? true : false; }

   168 RI_INLINE bool operator>=(float a, const RIfloat &b)			{ return a >= b.v ? true : false; }

   169 RI_INLINE bool operator>=(const RIfloat &a, float b)			{ return a.v >= b ? true : false; }

   170 RI_INLINE bool operator==(const RIfloat &a, const RIfloat &b)	{ return a.v == b.v ? true : false; }

   171 RI_INLINE bool operator==(float a, const RIfloat &b)			{ return a == b.v ? true : false; }

   172 RI_INLINE bool operator==(const RIfloat &a, float b)			{ return a.v == b ? true : false; }

   173 RI_INLINE bool operator!=(const RIfloat &a, const RIfloat &b)	{ return a.v != b.v ? true : false; }

   174 RI_INLINE bool operator!=(float a, const RIfloat &b)			{ return a != b.v ? true : false; }

   175 RI_INLINE bool operator!=(const RIfloat &a, float b)			{ return a.v != b ? true : false; }

   176 

   177 #else

   178 typedef float RIfloat;

   179 #endif

   180 

   181 #define	PI						3.141592654f

   182 

   183 RI_INLINE RIfloat	RI_MAX(RIfloat a, RIfloat b)				{ return (a > b) ? a : b; }

   184 RI_INLINE RIfloat	RI_MIN(RIfloat a, RIfloat b)				{ return (a < b) ? a : b; }

   185 RI_INLINE RIfloat	RI_CLAMP(RIfloat a, RIfloat l, RIfloat h)	{ if(RI_ISNAN(a)) return l; RI_ASSERT(l <= h); return (a < l) ? l : (a > h) ? h : a; }

   186 RI_INLINE void		RI_SWAP(RIfloat &a, RIfloat &b)				{ RIfloat tmp = a; a = b; b = tmp; }

   187 RI_INLINE RIfloat	RI_ABS(RIfloat a)							{ return (a < 0.0f) ? -a : a; }

   188 RI_INLINE RIfloat	RI_SQR(RIfloat a)							{ return a * a; }

   189 RI_INLINE RIfloat	RI_DEG_TO_RAD(RIfloat a)					{ return a * PI / 180.0f; }

   190 RI_INLINE RIfloat	RI_RAD_TO_DEG(RIfloat a)					{ return a * 180.0f/ PI; }

   191 RI_INLINE RIfloat	RI_MOD(RIfloat a, RIfloat b)				{ if(RI_ISNAN(a) || RI_ISNAN(b)) return 0.0f; RI_ASSERT(b >= 0.0f); if(b == 0.0f) return 0.0f; RIfloat f = (RIfloat)fmod(a, b); if(f < 0.0f) f += b; RI_ASSERT(f >= 0.0f && f <= b); return f; }

   192 

   193 RI_INLINE int		RI_INT_MAX(int a, int b)			{ return (a > b) ? a : b; }

   194 RI_INLINE int		RI_INT_MIN(int a, int b)			{ return (a < b) ? a : b; }

   195 RI_INLINE void		RI_INT_SWAP(int &a, int &b)			{ int tmp = a; a = b; b = tmp; }

   196 RI_INLINE int		RI_INT_MOD(int a, int b)			{ RI_ASSERT(b >= 0); if(!b) return 0; int i = a % b; if(i < 0) i += b; RI_ASSERT(i >= 0 && i < b); return i; }

   197 RI_INLINE int		RI_INT_ADDSATURATE(int a, int b)	{ RI_ASSERT(b >= 0); int r = a + b; return (r >= a) ? r : RI_INT32_MAX; }

   198 

   199 class Matrix3x3;

   200 class Vector2;

   201 class Vector3;

   202 

   203 //==============================================================================================

   204 

   205 //MatrixRxC, R = number of rows, C = number of columns

   206 //indexing: matrix[row][column]

   207 //Matrix3x3 inline functions cannot be inside the class because Vector3 is not defined yet when Matrix3x3 is defined

   208 

   209 class Matrix3x3

   210 {

   211 public:

   212 	RI_INLINE					Matrix3x3		();						//initialized to identity

   213 	RI_INLINE					Matrix3x3		( const Matrix3x3& m );

   214 	RI_INLINE					Matrix3x3		( RIfloat m00, RIfloat m01, RIfloat m02, RIfloat m10, RIfloat m11, RIfloat m12, RIfloat m20, RIfloat m21, RIfloat m22 );

   215 	RI_INLINE					~Matrix3x3		();

   216 	RI_INLINE Matrix3x3&		operator=		( const Matrix3x3& m );

   217 	RI_INLINE Vector3&			operator[]		( int i );				//returns a row vector

   218 	RI_INLINE const Vector3&	operator[]		( int i ) const;

   219 	RI_INLINE void				set				( RIfloat m00, RIfloat m01, RIfloat m02, RIfloat m10, RIfloat m11, RIfloat m12, RIfloat m20, RIfloat m21, RIfloat m22 );

   220 	RI_INLINE const Vector3		getRow			( int i ) const;

   221 	RI_INLINE const Vector3		getColumn		( int i ) const;

   222 	RI_INLINE void				setRow			( int i, const Vector3& v );

   223 	RI_INLINE void				setColumn		( int i, const Vector3& v );

   224 	RI_INLINE void				operator*=		( const Matrix3x3& m );

   225 	RI_INLINE void				operator*=		( RIfloat f );

   226 	RI_INLINE void				operator+=		( const Matrix3x3& m );

   227 	RI_INLINE void				operator-=		( const Matrix3x3& m );

   228 	RI_INLINE const Matrix3x3	operator-		() const;

   229 	RI_INLINE void				identity		();

   230 	RI_INLINE void				transpose		();

   231 	bool						invert			();	//if the matrix is singular, returns false and leaves it unmodified

   232 	RI_INLINE RIfloat				det				() const;

   233 	RI_INLINE bool				isAffine		() const;

   234 

   235 private:

   236 	RIfloat						matrix[3][3];

   237 };

   238 

   239 //==============================================================================================

   240 

   241 class Vector2

   242 {

   243 public:

   244 	RI_INLINE					Vector2			() : x(0.0f), y(0.0f)					{}

   245 	RI_INLINE					Vector2			( const Vector2& v ) : x(v.x), y(v.y)	{}

   246 	RI_INLINE					Vector2			( RIfloat fx, RIfloat fy ) : x(fx), y(fy)	{}

   247 	RI_INLINE					~Vector2		()								{}

   248 	RI_INLINE Vector2&			operator=		( const Vector2& v )			{ x = v.x; y = v.y; return *this; }

   249 	RI_INLINE RIfloat&			operator[]		( int i )						{ RI_ASSERT(i>=0&&i<2); return (&x)[i]; }

   250 	RI_INLINE const RIfloat&	operator[]		( int i ) const					{ RI_ASSERT(i>=0&&i<2); return (&x)[i]; }

   251 	RI_INLINE void				set				( RIfloat fx, RIfloat fy )			{ x = fx; y = fy; }

   252 	RI_INLINE void				operator*=		( RIfloat f )						{ x *= f; y *= f; }

   253 	RI_INLINE void				operator+=		( const Vector2& v )			{ x += v.x; y += v.y; }

   254 	RI_INLINE void				operator-=		( const Vector2& v )			{ x -= v.x; y -= v.y; }

   255 	RI_INLINE const Vector2		operator-		() const						{ return Vector2(-x,-y); }

   256 	//if the vector is zero, returns false and leaves it unmodified

   257 	RI_INLINE bool				normalize		()								{ double l = (double)x*(double)x+(double)y*(double)y; if( l == 0.0 ) return false; l = 1.0 / sqrt(l); x = (RIfloat)((double)x * l); y = (RIfloat)((double)y * l); return true; }

   258 	RI_INLINE RIfloat			length			() const						{ return (RIfloat)sqrt((double)x*(double)x+(double)y*(double)y); }

   259 	RI_INLINE void				scale			( const Vector2& v )			{ x *= v.x; y *= v.y; }	//component-wise scale

   260 	RI_INLINE void				negate			()								{ x = -x; y = -y; }

   261 

   262 	RIfloat						x,y;

   263 };

   264 

   265 //==============================================================================================

   266 

   267 class Vector3

   268 {

   269 public:

   270 	RI_INLINE					Vector3			() : x(0.0f), y(0.0f), z(0.0f)							{}

   271 	RI_INLINE					Vector3			( const Vector3& v ) : x(v.x), y(v.y), z(v.z)			{}

   272 	RI_INLINE					Vector3			( RIfloat fx, RIfloat fy, RIfloat fz ) : x(fx), y(fy), z(fz)	{}

   273 	RI_INLINE					~Vector3		()								{}

   274 	RI_INLINE Vector3&			operator=		( const Vector3& v )			{ x = v.x; y = v.y; z = v.z; return *this; }

   275 	RI_INLINE RIfloat&			operator[]		( int i )						{ RI_ASSERT(i>=0&&i<3); return (&x)[i]; }

   276 	RI_INLINE const RIfloat&	operator[]		( int i ) const					{ RI_ASSERT(i>=0&&i<3); return (&x)[i]; }

   277 	RI_INLINE void				set				( RIfloat fx, RIfloat fy, RIfloat fz ){ x = fx; y = fy; z = fz; }

   278 	RI_INLINE void				operator*=		( RIfloat f )						{ x *= f; y *= f; z *= f; }

   279 	RI_INLINE void				operator+=		( const Vector3& v )			{ x += v.x; y += v.y; z += v.z; }

   280 	RI_INLINE void				operator-=		( const Vector3& v )			{ x -= v.x; y -= v.y; z -= v.z; }

   281 	RI_INLINE const Vector3		operator-		() const						{ return Vector3(-x,-y,-z); }

   282 	//if the vector is zero, returns false and leaves it unmodified

   283 	RI_INLINE bool				normalize		()								{ double l = (double)x*(double)x+(double)y*(double)y+(double)z*(double)z; if( l == 0.0 ) return false; l = 1.0 / sqrt(l); x = (RIfloat)((double)x * l); y = (RIfloat)((double)y * l); z = (RIfloat)((double)z * l); return true; }

   284 	RI_INLINE RIfloat			length			() const						{ return (RIfloat)sqrt((double)x*(double)x+(double)y*(double)y+(double)z*(double)z); }

   285 	RI_INLINE void				scale			( const Vector3& v )			{ x *= v.x; y *= v.y; z *= v.z; }	//component-wise scale

   286 	RI_INLINE void				negate			()								{ x = -x; y = -y; z = -z; }

   287 

   288 	RIfloat						x,y,z;

   289 };

   290 

   291 //==============================================================================================

   292 

   293 //Vector2 global functions

   294 RI_INLINE bool			operator==	( const Vector2& v1, const Vector2& v2 )	{ return (v1.x == v2.x) && (v1.y == v2.y); }

   295 RI_INLINE bool			operator!=	( const Vector2& v1, const Vector2& v2 )	{ return (v1.x != v2.x) || (v1.y != v2.y); }

   296 RI_INLINE bool			isEqual		( const Vector2& v1, const Vector2& v2, RIfloat epsilon )	{ return RI_SQR(v2.x-v1.x) + RI_SQR(v2.y-v1.y) <= epsilon*epsilon; }

   297 RI_INLINE bool			isZero		( const Vector2& v )						{ return (v.x == 0.0f) && (v.y == 0.0f); }

   298 RI_INLINE const Vector2	operator*	( RIfloat f, const Vector2& v )				{ return Vector2(v.x*f,v.y*f); }

   299 RI_INLINE const Vector2	operator*	( const Vector2& v, RIfloat f )				{ return Vector2(v.x*f,v.y*f); }

   300 RI_INLINE const Vector2	operator+	( const Vector2& v1, const Vector2& v2 )	{ return Vector2(v1.x+v2.x, v1.y+v2.y); }

   301 RI_INLINE const Vector2	operator-	( const Vector2& v1, const Vector2& v2 )	{ return Vector2(v1.x-v2.x, v1.y-v2.y); }

   302 RI_INLINE RIfloat		dot			( const Vector2& v1, const Vector2& v2 )	{ return v1.x*v2.x+v1.y*v2.y; }

   303 //if v is a zero vector, returns a zero vector

   304 RI_INLINE const Vector2	normalize	( const Vector2& v )						{ double l = (double)v.x*(double)v.x+(double)v.y*(double)v.y; if( l != 0.0 ) l = 1.0 / sqrt(l); return Vector2((RIfloat)((double)v.x * l), (RIfloat)((double)v.y * l)); }

   305 //if onThis is a zero vector, returns a zero vector

   306 RI_INLINE const Vector2	project		( const Vector2& v, const Vector2& onThis ) { RIfloat l = dot(onThis,onThis); if( l != 0.0f ) l = dot(v, onThis)/l; return onThis * l; }

   307 RI_INLINE const Vector2	lerp		( const Vector2& v1, const Vector2& v2, RIfloat ratio )	{ return v1 + ratio * (v2 - v1); }

   308 RI_INLINE const Vector2	scale		( const Vector2& v1, const Vector2& v2 )	{ return Vector2(v1.x*v2.x, v1.y*v2.y); }

   309 //matrix * column vector. The input vector2 is implicitly expanded to (x,y,1)

   310 RI_INLINE const Vector2 affineTransform( const Matrix3x3& m, const Vector2& v )	{ RI_ASSERT(m.isAffine()); return Vector2(v.x * m[0][0] + v.y * m[0][1] + m[0][2], v.x * m[1][0] + v.y * m[1][1] + m[1][2]); }

   311 //matrix * column vector. The input vector2 is implicitly expanded to (x,y,0)

   312 RI_INLINE const Vector2 affineTangentTransform(const Matrix3x3& m, const Vector2& v)	{ RI_ASSERT(m.isAffine()); return Vector2(v.x * m[0][0] + v.y * m[0][1], v.x * m[1][0] + v.y * m[1][1]); }

   313 RI_INLINE const Vector2 perpendicularCW(const Vector2& v)						{ return Vector2(v.y, -v.x); }

   314 RI_INLINE const Vector2 perpendicularCCW(const Vector2& v)						{ return Vector2(-v.y, v.x); }

   315 RI_INLINE const Vector2 perpendicular(const Vector2& v, bool cw)				{ if(cw) return Vector2(v.y, -v.x); return Vector2(-v.y, v.x); }

   316 

   317 //==============================================================================================

   318 

   319 //Vector3 global functions

   320 RI_INLINE bool			operator==	( const Vector3& v1, const Vector3& v2 )	{ return (v1.x == v2.x) && (v1.y == v2.y) && (v1.z == v2.z); }

   321 RI_INLINE bool			operator!=	( const Vector3& v1, const Vector3& v2 )	{ return (v1.x != v2.x) || (v1.y != v2.y) || (v1.z != v2.z); }

   322 RI_INLINE bool			isEqual		( const Vector3& v1, const Vector3& v2, RIfloat epsilon )	{ return RI_SQR(v2.x-v1.x) + RI_SQR(v2.y-v1.y) + RI_SQR(v2.z-v1.z) <= epsilon*epsilon; }

   323 RI_INLINE const Vector3	operator*	( RIfloat f, const Vector3& v )				{ return Vector3(v.x*f,v.y*f,v.z*f); }

   324 RI_INLINE const Vector3	operator*	( const Vector3& v, RIfloat f )				{ return Vector3(v.x*f,v.y*f,v.z*f); }

   325 RI_INLINE const Vector3	operator+	( const Vector3& v1, const Vector3& v2 )	{ return Vector3(v1.x+v2.x, v1.y+v2.y, v1.z+v2.z); }

   326 RI_INLINE const Vector3	operator-	( const Vector3& v1, const Vector3& v2 )	{ return Vector3(v1.x-v2.x, v1.y-v2.y, v1.z-v2.z); }

   327 RI_INLINE RIfloat		dot			( const Vector3& v1, const Vector3& v2 )	{ return v1.x*v2.x+v1.y*v2.y+v1.z*v2.z; }

   328 RI_INLINE const Vector3	cross		( const Vector3& v1, const Vector3& v2 )	{ return Vector3( v1.y*v2.z-v1.z*v2.y, v1.z*v2.x-v1.x*v2.z, v1.x*v2.y-v1.y*v2.x ); }

   329 //if v is a zero vector, returns a zero vector

   330 RI_INLINE const Vector3	normalize	( const Vector3& v )						{ double l = (double)v.x*(double)v.x+(double)v.y*(double)v.y+(double)v.z*(double)v.z; if( l != 0.0 ) l = 1.0 / sqrt(l); return Vector3((RIfloat)((double)v.x * l), (RIfloat)((double)v.y * l), (RIfloat)((double)v.z * l)); }

   331 RI_INLINE const Vector3	lerp		( const Vector3& v1, const Vector3& v2, RIfloat ratio )	{ return v1 + ratio * (v2 - v1); }

   332 RI_INLINE const Vector3	scale		( const Vector3& v1, const Vector3& v2 )	{ return Vector3(v1.x*v2.x, v1.y*v2.y, v1.z*v2.z); }

   333 

   334 //==============================================================================================

   335 

   336 //matrix * column vector

   337 RI_INLINE const Vector3	operator*	( const Matrix3x3& m, const Vector3& v)		{ return Vector3( v.x*m[0][0]+v.y*m[0][1]+v.z*m[0][2], v.x*m[1][0]+v.y*m[1][1]+v.z*m[1][2], v.x*m[2][0]+v.y*m[2][1]+v.z*m[2][2] ); }

   338 

   339 //==============================================================================================

   340 

   341 //Matrix3x3 global functions

   342 RI_INLINE bool				operator==	( const Matrix3x3& m1, const Matrix3x3& m2 )	{ for(int i=0;i<3;i++) for(int j=0;j<3;j++) if( m1[i][j] != m2[i][j] ) return false; return true; }

   343 RI_INLINE bool				operator!=	( const Matrix3x3& m1, const Matrix3x3& m2 )	{ return !(m1 == m2); }

   344 RI_INLINE const Matrix3x3	operator*	( const Matrix3x3& m1, const Matrix3x3& m2 )	{ Matrix3x3 t; for(int i=0;i<3;i++) for(int j=0;j<3;j++) t[i][j] = m1[i][0] * m2[0][j] + m1[i][1] * m2[1][j] + m1[i][2] * m2[2][j]; return t; }

   345 RI_INLINE const Matrix3x3	operator*	( RIfloat f, const Matrix3x3& m )					{ Matrix3x3 t(m); t *= f; return t; }

   346 RI_INLINE const Matrix3x3	operator*	( const Matrix3x3& m, RIfloat f )					{ Matrix3x3 t(m); t *= f; return t; }

   347 RI_INLINE const Matrix3x3	operator+	( const Matrix3x3& m1, const Matrix3x3& m2 )	{ Matrix3x3 t(m1); t += m2; return t; }

   348 RI_INLINE const Matrix3x3	operator-	( const Matrix3x3& m1, const Matrix3x3& m2 )	{ Matrix3x3 t(m1); t -= m2; return t; }

   349 RI_INLINE const Matrix3x3	transpose	( const Matrix3x3& m )							{ Matrix3x3 t(m); t.transpose(); return t; }

   350 // if the matrix is singular, returns it unmodified

   351 RI_INLINE const Matrix3x3	invert		( const Matrix3x3& m )							{ Matrix3x3 t(m); t.invert(); return t; }

   352 

   353 //==============================================================================================

   354 

   355 //Matrix3x3 inline functions (cannot be inside the class because Vector3 is not defined yet when Matrix3x3 is defined)

   356 RI_INLINE					Matrix3x3::Matrix3x3	()									{ identity(); }

   357 RI_INLINE					Matrix3x3::Matrix3x3	( const Matrix3x3& m )				{ *this = m; }

   358 RI_INLINE					Matrix3x3::Matrix3x3	( RIfloat m00, RIfloat m01, RIfloat m02, RIfloat m10, RIfloat m11, RIfloat m12, RIfloat m20, RIfloat m21, RIfloat m22 )	{ set(m00,m01,m02,m10,m11,m12,m20,m21,m22); }

   359 RI_INLINE					Matrix3x3::~Matrix3x3	()									{}

   360 RI_INLINE Matrix3x3&		Matrix3x3::operator=	( const Matrix3x3& m )				{ for(int i=0;i<3;i++) for(int j=0;j<3;j++) matrix[i][j] = m.matrix[i][j]; return *this; }

   361 RI_INLINE Vector3&			Matrix3x3::operator[]	( int i )							{ RI_ASSERT(i>=0&&i<3); return (Vector3&)matrix[i][0]; }

   362 RI_INLINE const Vector3&	Matrix3x3::operator[]	( int i ) const						{ RI_ASSERT(i>=0&&i<3); return (const Vector3&)matrix[i][0]; }

   363 RI_INLINE void				Matrix3x3::set			( RIfloat m00, RIfloat m01, RIfloat m02, RIfloat m10, RIfloat m11, RIfloat m12, RIfloat m20, RIfloat m21, RIfloat m22 ) { matrix[0][0] = m00; matrix[0][1] = m01; matrix[0][2] = m02; matrix[1][0] = m10; matrix[1][1] = m11; matrix[1][2] = m12; matrix[2][0] = m20; matrix[2][1] = m21; matrix[2][2] = m22; }

   364 RI_INLINE const Vector3		Matrix3x3::getRow		( int i ) const						{ RI_ASSERT(i>=0&&i<3); return Vector3(matrix[i][0], matrix[i][1], matrix[i][2]); }

   365 RI_INLINE const Vector3		Matrix3x3::getColumn	( int i ) const						{ RI_ASSERT(i>=0&&i<3); return Vector3(matrix[0][i], matrix[1][i], matrix[2][i]); }

   366 RI_INLINE void				Matrix3x3::setRow		( int i, const Vector3& v )			{ RI_ASSERT(i>=0&&i<3); matrix[i][0] = v.x; matrix[i][1] = v.y; matrix[i][2] = v.z; }

   367 RI_INLINE void				Matrix3x3::setColumn	( int i, const Vector3& v )			{ RI_ASSERT(i>=0&&i<3); matrix[0][i] = v.x; matrix[1][i] = v.y; matrix[2][i] = v.z; }

   368 RI_INLINE void				Matrix3x3::operator*=	( const Matrix3x3& m )				{ *this = *this * m; }

   369 RI_INLINE void				Matrix3x3::operator*=	( RIfloat f )							{ for(int i=0;i<3;i++) for(int j=0;j<3;j++) matrix[i][j] *= f; }

   370 RI_INLINE void				Matrix3x3::operator+=	( const Matrix3x3& m )				{ for(int i=0;i<3;i++) for(int j=0;j<3;j++) matrix[i][j] += m.matrix[i][j]; }

   371 RI_INLINE void				Matrix3x3::operator-=	( const Matrix3x3& m )				{ for(int i=0;i<3;i++) for(int j=0;j<3;j++) matrix[i][j] -= m.matrix[i][j]; }

   372 RI_INLINE const Matrix3x3	Matrix3x3::operator-	() const							{ return Matrix3x3( -matrix[0][0],-matrix[0][1],-matrix[0][2], -matrix[1][0],-matrix[1][1],-matrix[1][2], -matrix[2][0],-matrix[2][1],-matrix[2][2]); }

   373 RI_INLINE void				Matrix3x3::identity		()									{ for(int i=0;i<3;i++) for(int j=0;j<3;j++) matrix[i][j] = (i == j) ? 1.0f : 0.0f; }

   374 RI_INLINE void				Matrix3x3::transpose	()									{ RI_SWAP(matrix[1][0], matrix[0][1]); RI_SWAP(matrix[2][0], matrix[0][2]); RI_SWAP(matrix[2][1], matrix[1][2]); }

   375 RI_INLINE RIfloat			Matrix3x3::det			() const							{ return matrix[0][0] * (matrix[1][1]*matrix[2][2] - matrix[2][1]*matrix[1][2]) + matrix[0][1] * (matrix[2][0]*matrix[1][2] - matrix[1][0]*matrix[2][2]) + matrix[0][2] * (matrix[1][0]*matrix[2][1] - matrix[2][0]*matrix[1][1]); }

   376 RI_INLINE bool				Matrix3x3::isAffine		() const							{ if(matrix[2][0] == 0.0f && matrix[2][1] == 0.0f && matrix[2][2] == 1.0f) return true; return false; }

   377 

   378 //==============================================================================================

   379 

   380 }	//namespace OpenVGRI

   381 

   382 #endif /* __RIMATH_H */