//====== Copyright 1996-2005, Valve Corporation, All rights reserved. =======//
//
// Purpose: 
//
// $NoKeywords: $
//
//=============================================================================//

#ifndef VECTOR_H
#define VECTOR_H

#ifdef _WIN32
#pragma once
#endif

#include <math.h>
#include <float.h>

// For vec_t, put this somewhere else?
#include "tier0/basetypes.h"

#if defined( _PS3 )
//#include <ssemath.h>
#include <vectormath/c/vectormath_aos.h>
#include "tier0/platform.h"
#include "mathlib/math_pfns.h"
#endif

#ifndef PLATFORM_PPC // we want our linux with xmm support
// For MMX intrinsics
#include <xmmintrin.h>
#endif

#ifndef ALIGN16_POST
#define ALIGN16_POST
#endif

#define NO_MALLOC_OVERRIDE
#if !defined(NO_MALLOC_OVERRIDE)
#include "tier0/memalloc.h"
#endif // !NO_MALLOC_OVERRIDE
#include "tier0/dbg.h"
#include "tier0/platform.h"
#if !defined( __SPU__ )
#include "tier0/threadtools.h"
#endif
#include "mathlib/vector2d.h"
#include "mathlib/math_pfns.h"
#include "vstdlib/random.h"

// Uncomment this to add extra Asserts to check for NANs, uninitialized vecs, etc.
//#define VECTOR_PARANOIA	1

// Uncomment this to make sure we don't do anything slow with our vectors
//#define VECTOR_NO_SLOW_OPERATIONS 1


// Used to make certain code easier to read.
#define X_INDEX	0
#define Y_INDEX	1
#define Z_INDEX	2


#ifdef VECTOR_PARANOIA
#define CHECK_VALID( _v)	Assert( (_v).IsValid() )
#else
#ifdef GNUC
#define CHECK_VALID( _v)
#else
#define CHECK_VALID( _v)	0
#endif
#endif

#define VecToString(v)	(static_cast<const char *>(CFmtStr("(%f, %f, %f)", (v).x, (v).y, (v).z))) // ** Note: this generates a temporary, don't hold reference!

class VectorByValue;
class QAngle;

//=========================================================
// 3D Vector
//=========================================================
class Vector3D
{
public:
	// Members
	vec_t x, y, z;

	// Construction/destruction:
	Vector3D(void);
	Vector3D(vec_t X, vec_t Y, vec_t Z);

	// Initialization
	void Init(vec_t ix = 0.0f, vec_t iy = 0.0f, vec_t iz = 0.0f);
	// TODO (Ilya): Should there be an init that takes a single float for consistency?

   // Got any nasty NAN's?
	bool IsValid() const;
	bool IsReasonable(float range = 1000000) const;		///< Check for reasonably-sized values (if used as a game world position)
	void Invalidate();

	// array access...
	vec_t operator[](int i) const;
	vec_t& operator[](int i);

	// Base address...
	vec_t* Base();
	vec_t const* Base() const;

	// Cast to Vector2D...
	Vector2D& AsVector2D();
	const Vector2D& AsVector2D() const;

	// Cast to QAngle...
	QAngle& AsQAngle();
	const QAngle& AsQAngle() const;

	// Initialization methods
	void Random(vec_t minVal, vec_t maxVal);
	inline void Zero(); ///< zero out a vector

	// equality
	bool operator==(const Vector3D& v) const;
	bool operator!=(const Vector3D& v) const;

	// arithmetic operations
	FORCEINLINE Vector3D& operator+=(const Vector3D& v);
	FORCEINLINE Vector3D& operator-=(const Vector3D& v);
	FORCEINLINE Vector3D& operator*=(const Vector3D& v);
	FORCEINLINE Vector3D& operator*=(float s);
	FORCEINLINE Vector3D& operator/=(const Vector3D& v);
	FORCEINLINE Vector3D& operator/=(float s);
	FORCEINLINE Vector3D& operator+=(float fl); ///< broadcast add
	FORCEINLINE Vector3D& operator-=(float fl); ///< broadcast sub			

// negate the vector components
	void	Negate();

	// Get the vector's magnitude.
	inline vec_t	Length() const;

	// Get the vector's magnitude squared.
	FORCEINLINE vec_t LengthSqr(void) const
	{
		CHECK_VALID(*this);
		return (x * x + y * y + z * z);
	}

	// Get one over the vector's length
	// via fast hardware approximation
	inline vec_t LengthRecipFast(void) const
	{
		return FastRSqrtFast(LengthSqr());
	}

	// return true if this vector is (0,0,0) within tolerance
	bool IsZero(float tolerance = 0.01f) const
	{
		return (x > -tolerance && x < tolerance&&
			y > -tolerance && y < tolerance&&
			z > -tolerance && z < tolerance);
	}


	// return true if this vector is exactly (0,0,0) -- only fast if vector is coming from memory, not registers
	inline bool IsZeroFast() const RESTRICT
	{
		static_assert(sizeof(vec_t) == sizeof(int));
		return (*(int*)(&x) == 0 &&
			*(int*)(&y) == 0 &&
			*(int*)(&z) == 0);
	}

	vec_t	NormalizeInPlace();								///< Normalize all components
	vec_t	NormalizeInPlaceSafe(const Vector3D& vFallback);///< Normalize all components
	Vector3D	Normalized() const;								///< Return normalized vector
	Vector3D	NormalizedSafe(const Vector3D& vFallback)const;		///< Return normalized vector, falling back to vFallback if the length of this is 0
	bool	IsLengthGreaterThan(float val) const;
	bool	IsLengthLessThan(float val) const;

	// check if a vector is within the box defined by two other vectors
	FORCEINLINE bool WithinAABox(Vector3D const& boxmin, Vector3D const& boxmax);

	// Get the distance from this vector to the other one.
	vec_t	DistTo(const Vector3D& vOther) const;

	// Get the distance from this vector to the other one squared.
	// NJS: note, VC wasn't inlining it correctly in several deeply nested inlines due to being an 'out of line' inline.  
	// may be able to tidy this up after switching to VC7
	FORCEINLINE vec_t DistToSqr(const Vector3D& vOther) const
	{
		Vector3D delta;

		delta.x = x - vOther.x;
		delta.y = y - vOther.y;
		delta.z = z - vOther.z;

		return delta.LengthSqr();
	}

	// Copy
	void	CopyToArray(float* rgfl) const;

	// Multiply, add, and assign to this (ie: *this = a + b * scalar). This
	// is about 12% faster than the actual vector equation (because it's done per-component
	// rather than per-vector).
	void	MulAdd(const Vector3D& a, const Vector3D& b, float scalar);

	// Dot product.
	vec_t	Dot(const Vector3D& vOther) const;

	// assignment
	Vector3D& operator=(const Vector3D& vOther);

	// returns 0, 1, 2 corresponding to the component with the largest absolute value
	inline int LargestComponent() const;
	inline vec_t LargestComponentValue() const;
	inline int SmallestComponent() const;
	inline vec_t SmallestComponentValue() const;

	// 2d
	vec_t	Length2D(void) const;
	vec_t	Length2DSqr(void) const;

	/// get the component of this vector parallel to some other given vector
	inline Vector3D  ProjectOnto(const Vector3D& onto);

	operator VectorByValue& () { return *((VectorByValue*)(this)); }
	operator const VectorByValue& () const { return *((const VectorByValue*)(this)); }

#ifndef VECTOR_NO_SLOW_OPERATIONS
	// copy constructors
//	Vector(const Vector &vOther);

	// arithmetic operations
	Vector3D	operator-(void) const;

	Vector3D	operator+(const Vector3D& v) const;
	Vector3D	operator-(const Vector3D& v) const;
	Vector3D	operator*(const Vector3D& v) const;
	Vector3D	operator/(const Vector3D& v) const;
	Vector3D	operator+(float fl) const;
	Vector3D	operator-(float fl) const;
	Vector3D	operator*(float fl) const;
	Vector3D	operator/(float fl) const;

	// Cross product between two vectors.
	Vector3D	Cross(const Vector3D& vOther) const;

	// Returns a vector with the min or max in X, Y, and Z.
	Vector3D	Min(const Vector3D& vOther) const;
	Vector3D	Max(const Vector3D& vOther) const;

#else

private:
	// No copy constructors allowed if we're in optimal mode
	Vector(const Vector& vOther);
#endif
};

// Zero the object -- necessary for CNetworkVar and possibly other cases.
inline void EnsureValidValue(Vector3D& x) { x.Zero(); }

#define USE_M64S defined( PLATFORM_WINDOWS_PC )



//=========================================================
// 4D Short Vector (aligned on 8-byte boundary)
//=========================================================
class ALIGN8 ShortVector
{
public:

	short x, y, z, w;

	// Initialization
	void Init(short ix = 0, short iy = 0, short iz = 0, short iw = 0);


#if USE_M64S
	__m64& AsM64() { return *(__m64*) & x; }
	const __m64& AsM64() const { return *(const __m64*) & x; }
#endif

	// Setter
	void Set(const ShortVector& vOther);
	void Set(const short ix, const short iy, const short iz, const short iw);

	// array access...
	short operator[](int i) const;
	short& operator[](int i);

	// Base address...
	short* Base();
	short const* Base() const;

	// equality
	bool operator==(const ShortVector& v) const;
	bool operator!=(const ShortVector& v) const;

	// Arithmetic operations
	FORCEINLINE ShortVector& operator+=(const ShortVector& v);
	FORCEINLINE ShortVector& operator-=(const ShortVector& v);
	FORCEINLINE ShortVector& operator*=(const ShortVector& v);
	FORCEINLINE ShortVector& operator*=(float s);
	FORCEINLINE ShortVector& operator/=(const ShortVector& v);
	FORCEINLINE ShortVector& operator/=(float s);
	FORCEINLINE ShortVector operator*(float fl) const;

private:

	// No copy constructors allowed if we're in optimal mode
//	ShortVector(ShortVector const& vOther);

	// No assignment operators either...
//	ShortVector& operator=( ShortVector const& src );

} ALIGN8_POST;






//=========================================================
// 4D Integer Vector
//=========================================================
class IntVector4D
{
public:

	int x, y, z, w;

	// Initialization
	void Init(int ix = 0, int iy = 0, int iz = 0, int iw = 0);

#if USE_M64S
	__m64& AsM64() { return *(__m64*) & x; }
	const __m64& AsM64() const { return *(const __m64*) & x; }
#endif

	// Setter
	void Set(const IntVector4D& vOther);
	void Set(const int ix, const int iy, const int iz, const int iw);

	// array access...
	int operator[](int i) const;
	int& operator[](int i);

	// Base address...
	int* Base();
	int const* Base() const;

	// equality
	bool operator==(const IntVector4D& v) const;
	bool operator!=(const IntVector4D& v) const;

	// Arithmetic operations
	FORCEINLINE IntVector4D& operator+=(const IntVector4D& v);
	FORCEINLINE IntVector4D& operator-=(const IntVector4D& v);
	FORCEINLINE IntVector4D& operator*=(const IntVector4D& v);
	FORCEINLINE IntVector4D& operator*=(float s);
	FORCEINLINE IntVector4D& operator/=(const IntVector4D& v);
	FORCEINLINE IntVector4D& operator/=(float s);
	FORCEINLINE IntVector4D operator*(float fl) const;

private:

	// No copy constructors allowed if we're in optimal mode
	//	IntVector4D(IntVector4D const& vOther);

	// No assignment operators either...
	//	IntVector4D& operator=( IntVector4D const& src );

};



//-----------------------------------------------------------------------------
// Allows us to specifically pass the vector by value when we need to
//-----------------------------------------------------------------------------
class VectorByValue : public Vector3D
{
public:
	// Construction/destruction:
	VectorByValue(void) : Vector3D() {}
	VectorByValue(vec_t X, vec_t Y, vec_t Z) : Vector3D(X, Y, Z) {}
	VectorByValue(const VectorByValue& vOther) { *this = vOther; }
};


//-----------------------------------------------------------------------------
// Utility to simplify table construction. No constructor means can use
// traditional C-style initialization
//-----------------------------------------------------------------------------
class TableVector
{
public:
	vec_t x, y, z;

	operator Vector3D& () { return *((Vector3D*)(this)); }
	operator const Vector3D& () const { return *((const Vector3D*)(this)); }

	// array access...
	inline vec_t& operator[](int i)
	{
		Assert((i >= 0) && (i < 3));
		return ((vec_t*)this)[i];
	}

	inline vec_t operator[](int i) const
	{
		Assert((i >= 0) && (i < 3));
		return ((vec_t*)this)[i];
	}
};


//-----------------------------------------------------------------------------
// Here's where we add all those lovely SSE optimized routines
//-----------------------------------------------------------------------------

class ALIGN16 VectorAligned : public Vector3D
{
public:
	inline VectorAligned(void) {};
	inline VectorAligned(vec_t X, vec_t Y, vec_t Z)
	{
		Init(X, Y, Z);
	}

#ifdef VECTOR_NO_SLOW_OPERATIONS

private:
	// No copy constructors allowed if we're in optimal mode
	VectorAligned(const VectorAligned& vOther);
	VectorAligned(const Vector& vOther);

#else
public:
	explicit VectorAligned(const Vector3D& vOther)
	{
		Init(vOther.x, vOther.y, vOther.z);
	}

	VectorAligned& operator=(const Vector3D& vOther)
	{
		Init(vOther.x, vOther.y, vOther.z);
		return *this;
	}

	VectorAligned& operator=(const VectorAligned& vOther)
	{
		// we know we're aligned, so use simd
		// we can't use the convenient abstract interface coz it gets declared later
#ifdef _X360
		XMStoreVector4A(Base(), XMLoadVector4A(vOther.Base()));
#elif _WIN32
		_mm_store_ps(Base(), _mm_load_ps(vOther.Base()));
#else
		Init(vOther.x, vOther.y, vOther.z);
#endif
		return *this;
	}


#endif
	float w;	// this space is used anyway

#if !defined(NO_MALLOC_OVERRIDE)
	void* operator new[](size_t nSize)
	{
		return MemAlloc_AllocAligned(nSize, 16);
	}

		void* operator new[](size_t nSize, const char* pFileName, int nLine)
	{
		return MemAlloc_AllocAlignedFileLine(nSize, 16, pFileName, nLine);
	}

		void* operator new[](size_t nSize, int /*nBlockUse*/, const char* pFileName, int nLine)
	{
		return MemAlloc_AllocAlignedFileLine(nSize, 16, pFileName, nLine);
	}

		void operator delete[](void* p)
	{
		MemAlloc_FreeAligned(p);
	}

		void operator delete[](void* p, const char* pFileName, int nLine)
	{
		MemAlloc_FreeAligned(p, pFileName, nLine);
	}

		void operator delete[](void* p, int /*nBlockUse*/, const char* pFileName, int nLine)
	{
		MemAlloc_FreeAligned(p, pFileName, nLine);
	}

		// please don't allocate a single quaternion...
		void* operator new   (size_t nSize)
	{
		return MemAlloc_AllocAligned(nSize, 16);
	}
	void* operator new   (size_t nSize, const char* pFileName, int nLine)
	{
		return MemAlloc_AllocAlignedFileLine(nSize, 16, pFileName, nLine);
	}
	void* operator new   (size_t nSize, int /*nBlockUse*/, const char* pFileName, int nLine)
	{
		return MemAlloc_AllocAlignedFileLine(nSize, 16, pFileName, nLine);
	}
	void operator delete (void* p)
	{
		MemAlloc_FreeAligned(p);
	}

	void operator delete (void* p, const char* pFileName, int nLine)
	{
		MemAlloc_FreeAligned(p, pFileName, nLine);
	}

	void operator delete (void* p, int /*nBlockUse*/, const char* pFileName, int nLine)
	{
		MemAlloc_FreeAligned(p, pFileName, nLine);
	}
#endif
} ALIGN16_POST;

//-----------------------------------------------------------------------------
// Vector related operations
//-----------------------------------------------------------------------------

// Vector clear
FORCEINLINE void VectorClear(Vector3D& a);

// Copy
FORCEINLINE void VectorCopy(const Vector3D& src, Vector3D& dst);

// Vector arithmetic
FORCEINLINE void VectorAdd(const Vector3D& a, const Vector3D& b, Vector3D& result);
FORCEINLINE void VectorSubtract(const Vector3D& a, const Vector3D& b, Vector3D& result);
FORCEINLINE void VectorMultiply(const Vector3D& a, vec_t b, Vector3D& result);
FORCEINLINE void VectorMultiply(const Vector3D& a, const Vector3D& b, Vector3D& result);
FORCEINLINE void VectorDivide(const Vector3D& a, vec_t b, Vector3D& result);
FORCEINLINE void VectorDivide(const Vector3D& a, const Vector3D& b, Vector3D& result);
inline void VectorScale(const Vector3D& in, vec_t scale, Vector3D& result);
void VectorMA(const Vector3D& start, float scale, const Vector3D& direction, Vector3D& dest);

// Vector equality with tolerance
bool VectorsAreEqual(const Vector3D& src1, const Vector3D& src2, float tolerance = 0.0f);

#define VectorExpand(v) (v).x, (v).y, (v).z


// Normalization
// FIXME: Can't use quite yet
//vec_t VectorNormalize( Vector& v );

// Length
inline vec_t VectorLength(const Vector3D& v);

// Dot Product
FORCEINLINE vec_t DotProduct(const Vector3D& a, const Vector3D& b);

// Cross product
void CrossProduct(const Vector3D& a, const Vector3D& b, Vector3D& result);

// Store the min or max of each of x, y, and z into the result.
void VectorMin(const Vector3D& a, const Vector3D& b, Vector3D& result);
void VectorMax(const Vector3D& a, const Vector3D& b, Vector3D& result);

// Linearly interpolate between two vectors
void VectorLerp(const Vector3D& src1, const Vector3D& src2, vec_t t, Vector3D& dest);
Vector3D VectorLerp(const Vector3D& src1, const Vector3D& src2, vec_t t);

FORCEINLINE Vector3D ReplicateToVector(float x)
{
	return Vector3D(x, x, x);
}

FORCEINLINE bool PointWithinViewAngle(Vector3D const& vecSrcPosition,
	Vector3D const& vecTargetPosition,
	Vector3D const& vecLookDirection, float flCosHalfFOV)
{
	Vector3D vecDelta = vecTargetPosition - vecSrcPosition;
	float cosDiff = DotProduct(vecLookDirection, vecDelta);

	if (flCosHalfFOV <= 0) // >180
	{
		// signs are different, answer is implicit
		if (cosDiff > 0)
			return true;

		// a/sqrt(b) > c  == a^2 < b * c ^2
		// IFF left and right sides are <= 0
		float flLen2 = vecDelta.LengthSqr();
		return (cosDiff * cosDiff <= flLen2 * flCosHalfFOV * flCosHalfFOV);
	}
	else // flCosHalfFOV > 0
	{
		// signs are different, answer is implicit
		if (cosDiff < 0)
			return false;

		// a/sqrt(b) > c  == a^2 > b * c ^2
		// IFF left and right sides are >= 0
		float flLen2 = vecDelta.LengthSqr();
		return (cosDiff * cosDiff >= flLen2 * flCosHalfFOV * flCosHalfFOV);
	}
}


#ifndef VECTOR_NO_SLOW_OPERATIONS

// Cross product
Vector3D CrossProduct(const Vector3D& a, const Vector3D& b);

// Random vector creation
Vector3D RandomVector(vec_t minVal, vec_t maxVal);

#endif

float RandomVectorInUnitSphere(Vector3D* pVector);
Vector3D RandomVectorInUnitSphere();
Vector3D RandomVectorInUnitSphere(IUniformRandomStream* pRnd);

float RandomVectorInUnitCircle(Vector2D* pVector);

Vector3D RandomVectorOnUnitSphere();
Vector3D RandomVectorOnUnitSphere(IUniformRandomStream* pRnd);


//-----------------------------------------------------------------------------
//
// Inlined Vector methods
//
//-----------------------------------------------------------------------------


//-----------------------------------------------------------------------------
// constructors
//-----------------------------------------------------------------------------
inline Vector3D::Vector3D(void)
{
#ifdef _DEBUG
#ifdef VECTOR_PARANOIA
	// Initialize to NAN to catch errors
	x = y = z = VEC_T_NAN;
#endif
#endif
}

inline Vector3D::Vector3D(vec_t X, vec_t Y, vec_t Z)
{
	x = X; y = Y; z = Z;
	CHECK_VALID(*this);
}

//inline Vector::Vector(const float *pFloat)					
//{
//	Assert( pFloat );
//	x = pFloat[0]; y = pFloat[1]; z = pFloat[2];	
//	CHECK_VALID(*this);
//} 

#if 0
//-----------------------------------------------------------------------------
// copy constructor
//-----------------------------------------------------------------------------

inline Vector::Vector(const Vector& vOther)
{
	CHECK_VALID(vOther);
	x = vOther.x; y = vOther.y; z = vOther.z;
}
#endif

//-----------------------------------------------------------------------------
// initialization
//-----------------------------------------------------------------------------

inline void Vector3D::Init(vec_t ix, vec_t iy, vec_t iz)
{
	x = ix; y = iy; z = iz;
	CHECK_VALID(*this);
}

#if !defined(__SPU__)
inline void Vector3D::Random(vec_t minVal, vec_t maxVal)
{
	x = RandomFloat(minVal, maxVal);
	y = RandomFloat(minVal, maxVal);
	z = RandomFloat(minVal, maxVal);
	CHECK_VALID(*this);
}
#endif

// This should really be a single opcode on the PowerPC (move r0 onto the vec reg)
inline void Vector3D::Zero()
{
	x = y = z = 0.0f;
}

inline void VectorClear(Vector3D& a)
{
	a.x = a.y = a.z = 0.0f;
}

//-----------------------------------------------------------------------------
// assignment
//-----------------------------------------------------------------------------

inline Vector3D& Vector3D::operator=(const Vector3D& vOther)
{
	CHECK_VALID(vOther);
	x = vOther.x; y = vOther.y; z = vOther.z;
	return *this;
}


//-----------------------------------------------------------------------------
// Array access
//-----------------------------------------------------------------------------
inline vec_t& Vector3D::operator[](int i)
{
	Assert((i >= 0) && (i < 3));
	return ((vec_t*)this)[i];
}

inline vec_t Vector3D::operator[](int i) const
{
	Assert((i >= 0) && (i < 3));
	return ((vec_t*)this)[i];
}


//-----------------------------------------------------------------------------
// Base address...
//-----------------------------------------------------------------------------
inline vec_t* Vector3D::Base()
{
	return (vec_t*)this;
}

inline vec_t const* Vector3D::Base() const
{
	return (vec_t const*)this;
}

//-----------------------------------------------------------------------------
// Cast to Vector2D...
//-----------------------------------------------------------------------------

inline Vector2D& Vector3D::AsVector2D()
{
	return *(Vector2D*)this;
}

inline const Vector2D& Vector3D::AsVector2D() const
{
	return *(const Vector2D*)this;
}

//-----------------------------------------------------------------------------
// Cast to QAngle...
//-----------------------------------------------------------------------------

inline QAngle& Vector3D::AsQAngle()
{
	return *(QAngle*)this;
}

inline const QAngle& Vector3D::AsQAngle() const
{
	return *(const QAngle*)this;
}

//-----------------------------------------------------------------------------
// IsValid?
//-----------------------------------------------------------------------------

inline bool Vector3D::IsValid() const
{
	return IsFinite(x) && IsFinite(y) && IsFinite(z);
}

//-----------------------------------------------------------------------------
// IsReasonable?
//-----------------------------------------------------------------------------
inline bool Vector3D::IsReasonable(float range) const
{
	return (Length() < range);
}

//-----------------------------------------------------------------------------
// Invalidate
//-----------------------------------------------------------------------------

inline void Vector3D::Invalidate()
{
	//#ifdef _DEBUG
	//#ifdef VECTOR_PARANOIA
	x = y = z = VEC_T_NAN;
	//#endif
	//#endif
}

//-----------------------------------------------------------------------------
// comparison
//-----------------------------------------------------------------------------

inline bool Vector3D::operator==(const Vector3D& src) const
{
	CHECK_VALID(src);
	CHECK_VALID(*this);
	return (src.x == x) && (src.y == y) && (src.z == z);
}

inline bool Vector3D::operator!=(const Vector3D& src) const
{
	CHECK_VALID(src);
	CHECK_VALID(*this);
	return (src.x != x) || (src.y != y) || (src.z != z);
}


//-----------------------------------------------------------------------------
// Copy
//-----------------------------------------------------------------------------

FORCEINLINE void VectorCopy(const Vector3D& src, Vector3D& dst)
{
	CHECK_VALID(src);
	dst.x = src.x;
	dst.y = src.y;
	dst.z = src.z;
}

inline void	Vector3D::CopyToArray(float* rgfl) const
{
	Assert(rgfl);
	CHECK_VALID(*this);
	rgfl[0] = x, rgfl[1] = y, rgfl[2] = z;
}

//-----------------------------------------------------------------------------
// standard math operations
//-----------------------------------------------------------------------------
// #pragma message("TODO: these should be SSE")

inline void Vector3D::Negate()
{
	CHECK_VALID(*this);
	x = -x; y = -y; z = -z;
}

FORCEINLINE  Vector3D& Vector3D::operator+=(const Vector3D& v)
{
	CHECK_VALID(*this);
	CHECK_VALID(v);
	x += v.x; y += v.y; z += v.z;
	return *this;
}

FORCEINLINE  Vector3D& Vector3D::operator-=(const Vector3D& v)
{
	CHECK_VALID(*this);
	CHECK_VALID(v);
	x -= v.x; y -= v.y; z -= v.z;
	return *this;
}

FORCEINLINE  Vector3D& Vector3D::operator*=(float fl)
{
	x *= fl;
	y *= fl;
	z *= fl;
	CHECK_VALID(*this);
	return *this;
}

FORCEINLINE  Vector3D& Vector3D::operator*=(const Vector3D& v)
{
	CHECK_VALID(v);
	x *= v.x;
	y *= v.y;
	z *= v.z;
	CHECK_VALID(*this);
	return *this;
}

// this ought to be an opcode.
FORCEINLINE Vector3D& Vector3D::operator+=(float fl)
{
	x += fl;
	y += fl;
	z += fl;
	CHECK_VALID(*this);
	return *this;
}

FORCEINLINE Vector3D& Vector3D::operator-=(float fl)
{
	x -= fl;
	y -= fl;
	z -= fl;
	CHECK_VALID(*this);
	return *this;
}



FORCEINLINE  Vector3D& Vector3D::operator/=(float fl)
{
	Assert(fl != 0.0f);
	float oofl = 1.0f / fl;
	x *= oofl;
	y *= oofl;
	z *= oofl;
	CHECK_VALID(*this);
	return *this;
}

FORCEINLINE  Vector3D& Vector3D::operator/=(const Vector3D& v)
{
	CHECK_VALID(v);
	Assert(v.x != 0.0f && v.y != 0.0f && v.z != 0.0f);
	x /= v.x;
	y /= v.y;
	z /= v.z;
	CHECK_VALID(*this);
	return *this;
}


// get the component of this vector parallel to some other given vector
inline Vector3D Vector3D::ProjectOnto(const Vector3D& onto)
{
	return onto * (this->Dot(onto) / (onto.LengthSqr()));
}


//-----------------------------------------------------------------------------
//
// Inlined Short Vector methods
//
//-----------------------------------------------------------------------------


inline void ShortVector::Init(short ix, short iy, short iz, short iw)
{
	x = ix; y = iy; z = iz; w = iw;
}

FORCEINLINE void ShortVector::Set(const ShortVector& vOther)
{
	x = vOther.x;
	y = vOther.y;
	z = vOther.z;
	w = vOther.w;
}

FORCEINLINE void ShortVector::Set(const short ix, const short iy, const short iz, const short iw)
{
	x = ix;
	y = iy;
	z = iz;
	w = iw;
}


//-----------------------------------------------------------------------------
// Array access
//-----------------------------------------------------------------------------
inline short ShortVector::operator[](int i) const
{
	Assert((i >= 0) && (i < 4));
	return ((short*)this)[i];
}

inline short& ShortVector::operator[](int i)
{
	Assert((i >= 0) && (i < 4));
	return ((short*)this)[i];
}

//-----------------------------------------------------------------------------
// Base address...
//-----------------------------------------------------------------------------
inline short* ShortVector::Base()
{
	return (short*)this;
}

inline short const* ShortVector::Base() const
{
	return (short const*)this;
}


//-----------------------------------------------------------------------------
// comparison
//-----------------------------------------------------------------------------

inline bool ShortVector::operator==(const ShortVector& src) const
{
	return (src.x == x) && (src.y == y) && (src.z == z) && (src.w == w);
}

inline bool ShortVector::operator!=(const ShortVector& src) const
{
	return (src.x != x) || (src.y != y) || (src.z != z) || (src.w != w);
}



//-----------------------------------------------------------------------------
// standard math operations
//-----------------------------------------------------------------------------

FORCEINLINE  ShortVector& ShortVector::operator+=(const ShortVector& v)
{
	x += v.x; y += v.y; z += v.z; w += v.w;
	return *this;
}

FORCEINLINE  ShortVector& ShortVector::operator-=(const ShortVector& v)
{
	x -= v.x; y -= v.y; z -= v.z; w -= v.w;
	return *this;
}

FORCEINLINE  ShortVector& ShortVector::operator*=(float fl)
{
	x = (short)(x * fl);
	y = (short)(y * fl);
	z = (short)(z * fl);
	w = (short)(w * fl);
	return *this;
}

FORCEINLINE  ShortVector& ShortVector::operator*=(const ShortVector& v)
{
	x = (short)(x * v.x);
	y = (short)(y * v.y);
	z = (short)(z * v.z);
	w = (short)(w * v.w);
	return *this;
}

FORCEINLINE  ShortVector& ShortVector::operator/=(float fl)
{
	Assert(fl != 0.0f);
	float oofl = 1.0f / fl;
	x = (short)(x * oofl);
	y = (short)(y * oofl);
	z = (short)(z * oofl);
	w = (short)(w * oofl);
	return *this;
}

FORCEINLINE  ShortVector& ShortVector::operator/=(const ShortVector& v)
{
	Assert(v.x != 0 && v.y != 0 && v.z != 0 && v.w != 0);
	x = (short)(x / v.x);
	y = (short)(y / v.y);
	z = (short)(z / v.z);
	w = (short)(w / v.w);
	return *this;
}

FORCEINLINE void ShortVectorMultiply(const ShortVector& src, float fl, ShortVector& res)
{
	Assert(IsFinite(fl));
	res.x = (short)(src.x * fl);
	res.y = (short)(src.y * fl);
	res.z = (short)(src.z * fl);
	res.w = (short)(src.w * fl);
}

FORCEINLINE ShortVector ShortVector::operator*(float fl) const
{
	ShortVector res;
	ShortVectorMultiply(*this, fl, res);
	return res;
}






//-----------------------------------------------------------------------------
//
// Inlined Integer Vector methods
//
//-----------------------------------------------------------------------------


inline void IntVector4D::Init(int ix, int iy, int iz, int iw)
{
	x = ix; y = iy; z = iz; w = iw;
}

FORCEINLINE void IntVector4D::Set(const IntVector4D& vOther)
{
	x = vOther.x;
	y = vOther.y;
	z = vOther.z;
	w = vOther.w;
}

FORCEINLINE void IntVector4D::Set(const int ix, const int iy, const int iz, const int iw)
{
	x = ix;
	y = iy;
	z = iz;
	w = iw;
}


//-----------------------------------------------------------------------------
// Array access
//-----------------------------------------------------------------------------
inline int IntVector4D::operator[](int i) const
{
	Assert((i >= 0) && (i < 4));
	return ((int*)this)[i];
}

inline int& IntVector4D::operator[](int i)
{
	Assert((i >= 0) && (i < 4));
	return ((int*)this)[i];
}

//-----------------------------------------------------------------------------
// Base address...
//-----------------------------------------------------------------------------
inline int* IntVector4D::Base()
{
	return (int*)this;
}

inline int const* IntVector4D::Base() const
{
	return (int const*)this;
}


//-----------------------------------------------------------------------------
// comparison
//-----------------------------------------------------------------------------

inline bool IntVector4D::operator==(const IntVector4D& src) const
{
	return (src.x == x) && (src.y == y) && (src.z == z) && (src.w == w);
}

inline bool IntVector4D::operator!=(const IntVector4D& src) const
{
	return (src.x != x) || (src.y != y) || (src.z != z) || (src.w != w);
}



//-----------------------------------------------------------------------------
// standard math operations
//-----------------------------------------------------------------------------

FORCEINLINE  IntVector4D& IntVector4D::operator+=(const IntVector4D& v)
{
	x += v.x; y += v.y; z += v.z; w += v.w;
	return *this;
}

FORCEINLINE  IntVector4D& IntVector4D::operator-=(const IntVector4D& v)
{
	x -= v.x; y -= v.y; z -= v.z; w -= v.w;
	return *this;
}

FORCEINLINE  IntVector4D& IntVector4D::operator*=(float fl)
{
	x = (int)(x * fl);
	y = (int)(y * fl);
	z = (int)(z * fl);
	w = (int)(w * fl);
	return *this;
}

FORCEINLINE  IntVector4D& IntVector4D::operator*=(const IntVector4D& v)
{
	x = (int)(x * v.x);
	y = (int)(y * v.y);
	z = (int)(z * v.z);
	w = (int)(w * v.w);
	return *this;
}

FORCEINLINE  IntVector4D& IntVector4D::operator/=(float fl)
{
	Assert(fl != 0.0f);
	float oofl = 1.0f / fl;
	x = (int)(x * oofl);
	y = (int)(y * oofl);
	z = (int)(z * oofl);
	w = (int)(w * oofl);
	return *this;
}

FORCEINLINE  IntVector4D& IntVector4D::operator/=(const IntVector4D& v)
{
	Assert(v.x != 0 && v.y != 0 && v.z != 0 && v.w != 0);
	x = (int)(x / v.x);
	y = (int)(y / v.y);
	z = (int)(z / v.z);
	w = (int)(w / v.w);
	return *this;
}

FORCEINLINE void IntVector4DMultiply(const IntVector4D& src, float fl, IntVector4D& res)
{
	Assert(IsFinite(fl));
	res.x = (int)(src.x * fl);
	res.y = (int)(src.y * fl);
	res.z = (int)(src.z * fl);
	res.w = (int)(src.w * fl);
}

FORCEINLINE IntVector4D IntVector4D::operator*(float fl) const
{
	IntVector4D res;
	IntVector4DMultiply(*this, fl, res);
	return res;
}



// =======================


FORCEINLINE void VectorAdd(const Vector3D& a, vec_t b, Vector3D& c)
{
	CHECK_VALID(a);
	CHECK_VALID(b);
	c.x = a.x + b;
	c.y = a.y + b;
	c.z = a.z + b;
}

FORCEINLINE void VectorAdd(const Vector3D& a, const Vector3D& b, Vector3D& c)
{
	CHECK_VALID(a);
	CHECK_VALID(b);
	c.x = a.x + b.x;
	c.y = a.y + b.y;
	c.z = a.z + b.z;
}

FORCEINLINE void VectorSubtract(const Vector3D& a, vec_t b, Vector3D& c)
{
	CHECK_VALID(a);
	CHECK_VALID(b);
	c.x = a.x - b;
	c.y = a.y - b;
	c.z = a.z - b;
}

FORCEINLINE void VectorSubtract(const Vector3D& a, const Vector3D& b, Vector3D& c)
{
	CHECK_VALID(a);
	CHECK_VALID(b);
	c.x = a.x - b.x;
	c.y = a.y - b.y;
	c.z = a.z - b.z;
}

FORCEINLINE void VectorMultiply(const Vector3D& a, vec_t b, Vector3D& c)
{
	CHECK_VALID(a);
	Assert(IsFinite(b));
	c.x = a.x * b;
	c.y = a.y * b;
	c.z = a.z * b;
}

FORCEINLINE void VectorMultiply(const Vector3D& a, const Vector3D& b, Vector3D& c)
{
	CHECK_VALID(a);
	CHECK_VALID(b);
	c.x = a.x * b.x;
	c.y = a.y * b.y;
	c.z = a.z * b.z;
}

// for backwards compatibility
inline void VectorScale(const Vector3D& in, vec_t scale, Vector3D& result)
{
	VectorMultiply(in, scale, result);
}


FORCEINLINE void VectorDivide(const Vector3D& a, vec_t b, Vector3D& c)
{
	CHECK_VALID(a);
	Assert(b != 0.0f);
	vec_t oob = 1.0f / b;
	c.x = a.x * oob;
	c.y = a.y * oob;
	c.z = a.z * oob;
}

FORCEINLINE void VectorDivide(const Vector3D& a, const Vector3D& b, Vector3D& c)
{
	CHECK_VALID(a);
	CHECK_VALID(b);
	Assert((b.x != 0.0f) && (b.y != 0.0f) && (b.z != 0.0f));
	c.x = a.x / b.x;
	c.y = a.y / b.y;
	c.z = a.z / b.z;
}

// FIXME: Remove
// For backwards compatibility
inline void	Vector3D::MulAdd(const Vector3D& a, const Vector3D& b, float scalar)
{
	CHECK_VALID(a);
	CHECK_VALID(b);
	x = a.x + b.x * scalar;
	y = a.y + b.y * scalar;
	z = a.z + b.z * scalar;
}

inline void VectorLerp(const Vector3D& src1, const Vector3D& src2, vec_t t, Vector3D& dest)
{
	CHECK_VALID(src1);
	CHECK_VALID(src2);
	dest.x = src1.x + (src2.x - src1.x) * t;
	dest.y = src1.y + (src2.y - src1.y) * t;
	dest.z = src1.z + (src2.z - src1.z) * t;
}

inline Vector3D VectorLerp(const Vector3D& src1, const Vector3D& src2, vec_t t)
{
	Vector3D result;
	VectorLerp(src1, src2, t, result);
	return result;
}

//-----------------------------------------------------------------------------
// Temporary storage for vector results so const Vector& results can be returned
//-----------------------------------------------------------------------------
#if !defined(__SPU__)
inline Vector3D& AllocTempVector()
{
	static Vector3D s_vecTemp[128];
	static CInterlockedInt s_nIndex;

	int nIndex;
	for (;;)
	{
		int nOldIndex = s_nIndex;
		nIndex = ((nOldIndex + 0x10001) & 0x7F);

		if (s_nIndex.AssignIf(nOldIndex, nIndex))
		{
			break;
		}
		ThreadPause();
	}
	return s_vecTemp[nIndex];
}
#endif


//-----------------------------------------------------------------------------
// dot, cross
//-----------------------------------------------------------------------------
FORCEINLINE vec_t DotProduct(const Vector3D& a, const Vector3D& b)
{
	CHECK_VALID(a);
	CHECK_VALID(b);
	return(a.x * b.x + a.y * b.y + a.z * b.z);
}

// for backwards compatibility
inline vec_t Vector3D::Dot(const Vector3D& vOther) const
{
	CHECK_VALID(vOther);
	return DotProduct(*this, vOther);
}

inline int Vector3D::LargestComponent() const
{
	float flAbsx = fabs(x);
	float flAbsy = fabs(y);
	float flAbsz = fabs(z);
	if (flAbsx > flAbsy)
	{
		if (flAbsx > flAbsz)
			return X_INDEX;
		return Z_INDEX;
	}
	if (flAbsy > flAbsz)
		return Y_INDEX;
	return Z_INDEX;
}

inline int Vector3D::SmallestComponent() const
{
	float flAbsx = fabs(x);
	float flAbsy = fabs(y);
	float flAbsz = fabs(z);
	if (flAbsx < flAbsy)
	{
		if (flAbsx < flAbsz)
			return X_INDEX;
		return Z_INDEX;
	}
	if (flAbsy < flAbsz)
		return Y_INDEX;
	return Z_INDEX;
}


inline float Vector3D::LargestComponentValue() const
{
	float flAbsX = fabs(x);
	float flAbsY = fabs(y);
	float flAbsZ = fabs(z);
	return MAX(MAX(flAbsX, flAbsY), flAbsZ);
}

inline float Vector3D::SmallestComponentValue() const
{
	float flAbsX = fabs(x);
	float flAbsY = fabs(y);
	float flAbsZ = fabs(z);
	return MIN(MIN(flAbsX, flAbsY), flAbsZ);
}


inline void CrossProduct(const Vector3D& a, const Vector3D& b, Vector3D& result)
{
	CHECK_VALID(a);
	CHECK_VALID(b);
	Assert(&a != &result);
	Assert(&b != &result);
	result.x = a.y * b.z - a.z * b.y;
	result.y = a.z * b.x - a.x * b.z;
	result.z = a.x * b.y - a.y * b.x;
}

inline vec_t DotProductAbs(const Vector3D& v0, const Vector3D& v1)
{
	CHECK_VALID(v0);
	CHECK_VALID(v1);
	return FloatMakePositive(v0.x * v1.x) + FloatMakePositive(v0.y * v1.y) + FloatMakePositive(v0.z * v1.z);
}

inline vec_t DotProductAbs(const Vector3D& v0, const float* v1)
{
	return FloatMakePositive(v0.x * v1[0]) + FloatMakePositive(v0.y * v1[1]) + FloatMakePositive(v0.z * v1[2]);
}

//-----------------------------------------------------------------------------
// length
//-----------------------------------------------------------------------------

inline vec_t VectorLength(const Vector3D& v)
{
	CHECK_VALID(v);
	return (vec_t)FastSqrt(v.x * v.x + v.y * v.y + v.z * v.z);
}


inline vec_t Vector3D::Length(void) const
{
	CHECK_VALID(*this);
	return VectorLength(*this);
}


//-----------------------------------------------------------------------------
// Normalization
//-----------------------------------------------------------------------------

/*
// FIXME: Can't use until we're un-macroed in mathlib.h
inline vec_t VectorNormalize( Vector& v )
{
	Assert( v.IsValid() );
	vec_t l = v.Length();
	if (l != 0.0f)
	{
		v /= l;
	}
	else
	{
		// FIXME:
		// Just copying the existing implementation; shouldn't res.z == 0?
		v.x = v.y = 0.0f; v.z = 1.0f;
	}
	return l;
}
*/


// check a point against a box
bool Vector3D::WithinAABox(Vector3D const& boxmin, Vector3D const& boxmax)
{
	return (
		(x >= boxmin.x) && (x <= boxmax.x) &&
		(y >= boxmin.y) && (y <= boxmax.y) &&
		(z >= boxmin.z) && (z <= boxmax.z)
		);
}

//-----------------------------------------------------------------------------
// Get the distance from this vector to the other one 
//-----------------------------------------------------------------------------
inline vec_t Vector3D::DistTo(const Vector3D& vOther) const
{
	Vector3D delta;
	VectorSubtract(*this, vOther, delta);
	return delta.Length();
}


//-----------------------------------------------------------------------------
// Float equality with tolerance
//-----------------------------------------------------------------------------
inline bool FloatsAreEqual(float f1, float f2, float flTolerance)
{
	// Sergiy: the implementation in Source2 is very inefficient, trying to start with a clean slate here, hopefully will reintegrate back to Source2
	const float flAbsToleranceThreshold = 0.000003814697265625; // 2 ^ -FLOAT_EQUALITY_NOISE_CUTOFF, 
	return fabsf(f1 - f2) <= flTolerance * (fabsf(f1) + fabsf(f2)) + flAbsToleranceThreshold;
}


//-----------------------------------------------------------------------------
// Vector equality with percentage tolerance
// are all components within flPercentageTolerance (expressed as a percentage of the larger component, per component)?
// and all components have the same sign
//-----------------------------------------------------------------------------
inline bool VectorsAreWithinPercentageTolerance(const Vector3D& src1, const Vector3D& src2, float flPercentageTolerance)
{
	if (!FloatsAreEqual(src1.x, src2.x, flPercentageTolerance))
		return false;

	if (!FloatsAreEqual(src1.y, src2.y, flPercentageTolerance))
		return false;

	return (FloatsAreEqual(src1.z, src2.z, flPercentageTolerance));
}



//-----------------------------------------------------------------------------
// Vector equality with tolerance
//-----------------------------------------------------------------------------
inline bool VectorsAreEqual(const Vector3D& src1, const Vector3D& src2, float tolerance)
{
	if (FloatMakePositive(src1.x - src2.x) > tolerance)
		return false;
	if (FloatMakePositive(src1.y - src2.y) > tolerance)
		return false;
	return (FloatMakePositive(src1.z - src2.z) <= tolerance);
}


//-----------------------------------------------------------------------------
// Computes the closest point to vecTarget no farther than flMaxDist from vecStart
//-----------------------------------------------------------------------------
inline void ComputeClosestPoint(const Vector3D& vecStart, float flMaxDist, const Vector3D& vecTarget, Vector3D* pResult)
{
	Vector3D vecDelta;
	VectorSubtract(vecTarget, vecStart, vecDelta);
	float flDistSqr = vecDelta.LengthSqr();
	if (flDistSqr <= flMaxDist * flMaxDist)
	{
		*pResult = vecTarget;
	}
	else
	{
		vecDelta /= FastSqrt(flDistSqr);
		VectorMA(vecStart, flMaxDist, vecDelta, *pResult);
	}
}


//-----------------------------------------------------------------------------
// Takes the absolute value of a vector
//-----------------------------------------------------------------------------
inline void VectorAbs(const Vector3D& src, Vector3D& dst)
{
	dst.x = FloatMakePositive(src.x);
	dst.y = FloatMakePositive(src.y);
	dst.z = FloatMakePositive(src.z);
}

inline Vector3D VectorAbs(const Vector3D& src)
{
	return Vector3D(fabsf(src.x), fabsf(src.y), fabsf(src.z));
}


//-----------------------------------------------------------------------------
//
// Slow methods
//
//-----------------------------------------------------------------------------

#ifndef VECTOR_NO_SLOW_OPERATIONS

//-----------------------------------------------------------------------------
// Returns a vector with the min or max in X, Y, and Z.
//-----------------------------------------------------------------------------
inline Vector3D Vector3D::Min(const Vector3D& vOther) const
{
	return Vector3D(x < vOther.x ? x : vOther.x,
		y < vOther.y ? y : vOther.y,
		z < vOther.z ? z : vOther.z);
}

inline Vector3D Vector3D::Max(const Vector3D& vOther) const
{
	return Vector3D(x > vOther.x ? x : vOther.x,
		y > vOther.y ? y : vOther.y,
		z > vOther.z ? z : vOther.z);
}


//-----------------------------------------------------------------------------
// arithmetic operations
//-----------------------------------------------------------------------------

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

inline Vector3D Vector3D::operator+(float fl) const
{
	Vector3D res;
	VectorAdd(*this, fl, res);
	return res;
}

inline Vector3D Vector3D::operator+(const Vector3D& v) const
{
	Vector3D res;
	VectorAdd(*this, v, res);
	return res;
}

inline Vector3D Vector3D::operator-(float fl) const
{
	Vector3D res;
	VectorSubtract(*this, fl, res);
	return res;
}

inline Vector3D Vector3D::operator-(const Vector3D& v) const
{
	Vector3D res;
	VectorSubtract(*this, v, res);
	return res;
}

inline Vector3D Vector3D::operator*(float fl) const
{
	Vector3D res;
	VectorMultiply(*this, fl, res);
	return res;
}

inline Vector3D Vector3D::operator*(const Vector3D& v) const
{
	Vector3D res;
	VectorMultiply(*this, v, res);
	return res;
}

inline Vector3D Vector3D::operator/(float fl) const
{
	Vector3D res;
	VectorDivide(*this, fl, res);
	return res;
}

inline Vector3D Vector3D::operator/(const Vector3D& v) const
{
	Vector3D res;
	VectorDivide(*this, v, res);
	return res;
}

inline Vector3D operator*(float fl, const Vector3D& v)
{
	return v * fl;
}

//-----------------------------------------------------------------------------
// cross product
//-----------------------------------------------------------------------------

inline Vector3D Vector3D::Cross(const Vector3D& vOther) const
{
	Vector3D res;
	CrossProduct(*this, vOther, res);
	return res;
}

//-----------------------------------------------------------------------------
// 2D
//-----------------------------------------------------------------------------

inline vec_t Vector3D::Length2D(void) const
{
	return (vec_t)FastSqrt(x * x + y * y);
}

inline vec_t Vector3D::Length2DSqr(void) const
{
	return (x * x + y * y);
}

inline Vector3D CrossProduct(const Vector3D& a, const Vector3D& b)
{
	return Vector3D(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 VectorMin(const Vector3D& a, const Vector3D& b, Vector3D& result)
{
	result.x = fmin(a.x, b.x);
	result.y = fmin(a.y, b.y);
	result.z = fmin(a.z, b.z);
}

inline void VectorMax(const Vector3D& a, const Vector3D& b, Vector3D& result)
{
	result.x = fmax(a.x, b.x);
	result.y = fmax(a.y, b.y);
	result.z = fmax(a.z, b.z);
}

// and when you want to return the vector rather than cause a LHS with it...
inline Vector3D VectorMin(const Vector3D& a, const Vector3D& b)
{
	return Vector3D(fmin(a.x, b.x), fmin(a.y, b.y), fmin(a.z, b.z));
}

inline Vector3D VectorMax(const Vector3D& a, const Vector3D& b)
{
	return Vector3D(fmax(a.x, b.x), fmax(a.y, b.y), fmax(a.z, b.z));
}

inline float ComputeVolume(const Vector3D& vecMins, const Vector3D& vecMaxs)
{
	Vector3D vecDelta;
	VectorSubtract(vecMaxs, vecMins, vecDelta);
	return DotProduct(vecDelta, vecDelta);
}

#if !defined(__SPU__)
// Get a random vector.
inline Vector3D RandomVector(float minVal, float maxVal)
{
	Vector3D random;
	random.Random(minVal, maxVal);
	return random;
}
#endif

#endif //slow

//-----------------------------------------------------------------------------
// Helper debugging stuff....
//-----------------------------------------------------------------------------
/*
inline bool operator==(float const* f, const Vector3D& v)
{
	// AIIIEEEE!!!!
	Assert(0);
	return false;
}

inline bool operator==(const Vector3D& v, float const* f)
{
	// AIIIEEEE!!!!
	Assert(0);
	return false;
}

inline bool operator!=(float const* f, const Vector3D& v)
{
	// AIIIEEEE!!!!
	Assert(0);
	return false;
}

inline bool operator!=(const Vector3D& v, float const* f)
{
	// AIIIEEEE!!!!
	Assert(0);
	return false;
}
*/

// return a vector perpendicular to another, with smooth variation. The difference between this and
// something like VectorVectors is that there are now discontinuities. _unlike_ VectorVectors,
// you won't get an "u
void VectorPerpendicularToVector(Vector3D const& in, Vector3D* pvecOut);

inline const Vector3D VectorPerpendicularToVector(const Vector3D& in)
{
	Vector3D out;
	VectorPerpendicularToVector(in, &out);
	return out;
}

//-----------------------------------------------------------------------------
// AngularImpulse
//-----------------------------------------------------------------------------
// AngularImpulse are exponential maps (an axis scaled by a "twist" angle in degrees)
typedef Vector3D AngularImpulse;

#ifndef VECTOR_NO_SLOW_OPERATIONS

#if !defined(__SPU__)
inline AngularImpulse RandomAngularImpulse(float minVal, float maxVal)
{
	AngularImpulse	angImp;
	angImp.Random(minVal, maxVal);
	return angImp;
}
#endif

#endif


//-----------------------------------------------------------------------------
// Quaternion
//-----------------------------------------------------------------------------

class RadianEuler;
class DegreeEuler;
class QAngle;

class Quaternion				// same data-layout as engine's vec4_t,
{								//		which is a vec_t[4]
public:
	inline Quaternion(void) {

		// Initialize to NAN to catch errors
#ifdef _DEBUG
#ifdef VECTOR_PARANOIA
		x = y = z = w = VEC_T_NAN;
#endif
#endif
	}
	inline Quaternion(vec_t ix, vec_t iy, vec_t iz, vec_t iw) : x(ix), y(iy), z(iz), w(iw) { }
	inline explicit Quaternion(RadianEuler const& angle);
	inline explicit Quaternion(DegreeEuler const& angle);

	inline void Init(vec_t ix = 0.0f, vec_t iy = 0.0f, vec_t iz = 0.0f, vec_t iw = 0.0f) { x = ix; y = iy; z = iz; w = iw; }
	inline void Init(const Vector3D& vImaginaryPart, float flRealPart) { x = vImaginaryPart.x; y = vImaginaryPart.y; z = vImaginaryPart.z; w = flRealPart; }

	bool IsValid() const;
	void Invalidate();

	bool operator==(const Quaternion& src) const;
	bool operator!=(const Quaternion& src) const;

	inline Quaternion Conjugate() const { return Quaternion(-x, -y, -z, w); }

	// 
	const Vector3D GetForward()const;
	const Vector3D GetLeft()const;
	const Vector3D GetUp()const;

	vec_t* Base() { return (vec_t*)this; }
	const vec_t* Base() const { return (vec_t*)this; }

	// convenience for debugging
	inline void Print() const;

	// Imaginary part
	Vector3D& ImaginaryPart() { return *(Vector3D*)this; }
	const Vector3D& ImaginaryPart() const { return *(Vector3D*)this; }
	float& RealPart() { return w; }
	float RealPart() const { return w; }
	inline QAngle ToQAngle() const;
	inline struct matrix3x4_t ToMatrix() const;

	// array access...
	vec_t operator[](int i) const;
	vec_t& operator[](int i);

	inline Quaternion operator+(void) const { return *this; }
	inline Quaternion operator-(void) const { return Quaternion(-x, -y, -z, -w); }

	vec_t x, y, z, w;
};

// Random Quaternion that is UNIFORMLY distributed over the S^3
// should be good for random generation of orientation for unit tests and for game
// NOTE: Nothing trivial like Quaternion(RandomAngle(0,180)) will do the trick , 
//       one needs to take special care to generate a uniformly distributed quaternion.
const Quaternion RandomQuaternion();
const Quaternion RandomQuaternion();
inline const Quaternion Conjugate(const Quaternion& q)
{
	return Quaternion(-q.x, -q.y, -q.z, q.w);
}



//-----------------------------------------------------------------------------
// Array access
//-----------------------------------------------------------------------------
inline vec_t& Quaternion::operator[](int i)
{
	Assert((i >= 0) && (i < 4));
	return ((vec_t*)this)[i];
}

inline vec_t Quaternion::operator[](int i) const
{
	Assert((i >= 0) && (i < 4));
	return ((vec_t*)this)[i];
}


//-----------------------------------------------------------------------------
// Equality test
//-----------------------------------------------------------------------------
inline bool Quaternion::operator==(const Quaternion& src) const
{
	return (x == src.x) && (y == src.y) && (z == src.z) && (w == src.w);
}

inline bool Quaternion::operator!=(const Quaternion& src) const
{
	return !operator==(src);
}


//-----------------------------------------------------------------------------
// Debugging only
//-----------------------------------------------------------------------------
void Quaternion::Print() const
{
#ifndef _CERT
#if !defined(__SPU__)
	DevMsg(eDLL_T::ENGINE, "q{ %.3fi + %.3fj + %.3fk + %.3f }", x, y, z, w);
#endif
#endif
}




//-----------------------------------------------------------------------------
// Binary operators
//-----------------------------------------------------------------------------
inline Quaternion operator+(const Quaternion& q1, const Quaternion& q2)
{
	return Quaternion(q1.x + q2.x, q1.y + q2.y, q1.z + q2.z, q1.w + q2.w);
}

inline Quaternion operator-(const Quaternion& q1, const Quaternion& q2)
{
	return Quaternion(q1.x - q2.x, q1.y - q2.y, q1.z - q2.z, q1.w - q2.w);
}

inline Quaternion operator*(float s, const Quaternion& q)
{
	return Quaternion(s * q.x, s * q.y, s * q.z, s * q.w);
}

inline Quaternion operator*(const Quaternion& q, float s)
{
	return Quaternion(q.x * s, q.y * s, q.z * s, q.w * s);
}

inline Quaternion operator/(const Quaternion& q, float s)
{
	Assert(s != 0.0f);
	return Quaternion(q.x / s, q.y / s, q.z / s, q.w / s);
}


//-----------------------------------------------------------------------------
// Quaternion equality with tolerance
//-----------------------------------------------------------------------------
inline bool QuaternionsAreEqual(const Quaternion& src1, const Quaternion& src2, float tolerance)
{
	if (FloatMakePositive(src1.x - src2.x) > tolerance)
		return false;
	if (FloatMakePositive(src1.y - src2.y) > tolerance)
		return false;
	if (FloatMakePositive(src1.z - src2.z) > tolerance)
		return false;
	return (FloatMakePositive(src1.w - src2.w) <= tolerance);
}


//-----------------------------------------------------------------------------
// Here's where we add all those lovely SSE optimized routines
//-----------------------------------------------------------------------------
class ALIGN16 QuaternionAligned : public Quaternion
{
public:
	inline QuaternionAligned(void) {};
	inline QuaternionAligned(vec_t X, vec_t Y, vec_t Z, vec_t W)
	{
		Init(X, Y, Z, W);
	}

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

#ifdef VECTOR_NO_SLOW_OPERATIONS

private:
	// No copy constructors allowed if we're in optimal mode
	QuaternionAligned(const QuaternionAligned& vOther);
	QuaternionAligned(const Quaternion& vOther);

#else
public:
	explicit QuaternionAligned(const Quaternion& vOther)
	{
		Init(vOther.x, vOther.y, vOther.z, vOther.w);
	}

	QuaternionAligned& operator=(const Quaternion& vOther)
	{
		Init(vOther.x, vOther.y, vOther.z, vOther.w);
		return *this;
	}

	QuaternionAligned& operator=(const QuaternionAligned& vOther)
	{
		// we know we're aligned, so use simd
		// we can't use the convenient abstract interface coz it gets declared later
#ifdef _X360
		XMStoreVector4A(Base(), XMLoadVector4A(vOther.Base()));
#elif _WIN32
		_mm_store_ps(Base(), _mm_load_ps(vOther.Base()));
#else
		Init(vOther.x, vOther.y, vOther.z, vOther.w);
#endif
		return *this;
	}

#endif

#if !defined(NO_MALLOC_OVERRIDE)
	void* operator new[](size_t nSize)
	{
		return MemAlloc_AllocAligned(nSize, 16);
	}

		void* operator new[](size_t nSize, const char* pFileName, int nLine)
	{
		return MemAlloc_AllocAlignedFileLine(nSize, 16, pFileName, nLine);
	}

		void* operator new[](size_t nSize, int /*nBlockUse*/, const char* pFileName, int nLine)
	{
		return MemAlloc_AllocAlignedFileLine(nSize, 16, pFileName, nLine);
	}

		void operator delete[](void* p)
	{
		MemAlloc_FreeAligned(p);
	}

		void operator delete[](void* p, const char* pFileName, int nLine)
	{
		MemAlloc_FreeAligned(p, pFileName, nLine);
	}

		void operator delete[](void* p, int /*nBlockUse*/, const char* pFileName, int nLine)
	{
		MemAlloc_FreeAligned(p, pFileName, nLine);
	}

		// please don't allocate a single quaternion...
		void* operator new   (size_t nSize)
	{
		return MemAlloc_AllocAligned(nSize, 16);
	}
	void* operator new   (size_t nSize, const char* pFileName, int nLine)
	{
		return MemAlloc_AllocAlignedFileLine(nSize, 16, pFileName, nLine);
	}
	void* operator new   (size_t nSize, int /*nBlockUse*/, const char* pFileName, int nLine)
	{
		return MemAlloc_AllocAlignedFileLine(nSize, 16, pFileName, nLine);
	}
	void operator delete (void* p)
	{
		MemAlloc_FreeAligned(p);
	}

	void operator delete (void* p, const char* pFileName, int nLine)
	{
		MemAlloc_FreeAligned(p, pFileName, nLine);
	}

	void operator delete (void* p, int /*nBlockUse*/, const char* pFileName, int nLine)
	{
		MemAlloc_FreeAligned(p, pFileName, nLine);
	}
#endif
} ALIGN16_POST;


//-----------------------------------------------------------------------------
// Src data hasn't changed, but work data is of a form more friendly for SPU
//-----------------------------------------------------------------------------
#if defined( _PS3 )
//typedef Vector		BoneVector;
typedef VectorAligned		BoneVector;
typedef QuaternionAligned	BoneQuaternion;
typedef QuaternionAligned	BoneQuaternionAligned;
#else
typedef Vector3D				BoneVector;
typedef Quaternion			BoneQuaternion;
typedef QuaternionAligned	BoneQuaternionAligned;
#endif

//-----------------------------------------------------------------------------
// Radian Euler angle aligned to axis (NOT ROLL/PITCH/YAW)
//-----------------------------------------------------------------------------
class QAngle;
#define VEC_DEG2RAD( a ) (a) * (3.14159265358979323846f / 180.0f)
#define VEC_RAD2DEG( a ) (a) * (180.0f / 3.14159265358979323846f)
class RadianEuler
{
public:
	inline RadianEuler(void) { }
	inline RadianEuler(vec_t X, vec_t Y, vec_t Z) { x = X; y = Y; z = Z; }
	inline explicit RadianEuler(Quaternion const& q);
	inline explicit RadianEuler(QAngle const& angles);
	inline explicit RadianEuler(DegreeEuler const& angles);

	// Initialization
	inline void Init(vec_t ix = 0.0f, vec_t iy = 0.0f, vec_t iz = 0.0f) { x = ix; y = iy; z = iz; }

	//	conversion to qangle
	QAngle ToQAngle(void) const;
	bool IsValid() const;
	void Invalidate();

	inline vec_t* Base() { return &x; }
	inline const vec_t* Base() const { return &x; }

	// array access...
	vec_t operator[](int i) const;
	vec_t& operator[](int i);

	vec_t x, y, z;
};


extern void AngleQuaternion(RadianEuler const& angles, Quaternion& qt);
extern void QuaternionAngles(Quaternion const& q, RadianEuler& angles);
inline Quaternion::Quaternion(RadianEuler const& angle)
{
	AngleQuaternion(angle, *this);
}

inline bool Quaternion::IsValid() const
{
	return IsFinite(x) && IsFinite(y) && IsFinite(z) && IsFinite(w);
}


FORCEINLINE float QuaternionLength(const Quaternion& q)
{
	return sqrtf(q.x * q.x + q.y * q.y + q.z * q.z + q.w * q.w);
}

FORCEINLINE bool QuaternionIsNormalized(const Quaternion& q, float  flTolerance = 1e-6f)
{
	float flLen = QuaternionLength(q);
	return (fabs(flLen - 1.0) < flTolerance);
}

inline void Quaternion::Invalidate()
{
	//#ifdef _DEBUG
	//#ifdef VECTOR_PARANOIA
	x = y = z = w = VEC_T_NAN;
	//#endif
	//#endif
}

inline RadianEuler::RadianEuler(Quaternion const& q)
{
	QuaternionAngles(q, *this);
}

inline void VectorCopy(RadianEuler const& src, RadianEuler& dst)
{
	CHECK_VALID(src);
	dst.x = src.x;
	dst.y = src.y;
	dst.z = src.z;
}

inline void VectorScale(RadianEuler const& src, float b, RadianEuler& dst)
{
	CHECK_VALID(src);
	Assert(IsFinite(b));
	dst.x = src.x * b;
	dst.y = src.y * b;
	dst.z = src.z * b;
}

inline bool RadianEuler::IsValid() const
{
	return IsFinite(x) && IsFinite(y) && IsFinite(z);
}

inline void RadianEuler::Invalidate()
{
	//#ifdef _DEBUG
	//#ifdef VECTOR_PARANOIA
	x = y = z = VEC_T_NAN;
	//#endif
	//#endif
}


//-----------------------------------------------------------------------------
// Array access
//-----------------------------------------------------------------------------
inline vec_t& RadianEuler::operator[](int i)
{
	Assert((i >= 0) && (i < 3));
	return ((vec_t*)this)[i];
}

inline vec_t RadianEuler::operator[](int i) const
{
	Assert((i >= 0) && (i < 3));
	return ((vec_t*)this)[i];
}


//-----------------------------------------------------------------------------
// Degree Euler angle aligned to axis (NOT ROLL/PITCH/YAW)
//-----------------------------------------------------------------------------
class DegreeEuler
{
public:
	///\name Initialization 
	//@{
	inline DegreeEuler(void) ///< Create with un-initialized components. If VECTOR_PARANOIA is set, will init with NANS.
	{
		// Initialize to NAN to catch errors
#ifdef VECTOR_PARANOIA
		x = y = z = VEC_T_NAN;
#endif
	}
	inline DegreeEuler(vec_t X, vec_t Y, vec_t Z) { x = X; y = Y; z = Z; }
	inline explicit DegreeEuler(Quaternion const& q);
	inline explicit DegreeEuler(QAngle const& angles);
	inline explicit DegreeEuler(RadianEuler const& angles);

	// Initialization
	inline void Init(vec_t ix = 0.0f, vec_t iy = 0.0f, vec_t iz = 0.0f) { x = ix; y = iy; z = iz; }

	inline QAngle ToQAngle() const;

	//	conversion to qangle
	bool IsValid() const;
	void Invalidate();

	inline vec_t* Base() { return &x; }
	inline const vec_t* Base() const { return &x; }

	// array access...
	vec_t operator[](int i) const;
	vec_t& operator[](int i);

	vec_t x, y, z;
};


//-----------------------------------------------------------------------------
// DegreeEuler equality with tolerance
//-----------------------------------------------------------------------------
inline bool DegreeEulersAreEqual(const DegreeEuler& src1, const DegreeEuler& src2, float tolerance = 0.0f)
{
	if (FloatMakePositive(src1.x - src2.x) > tolerance)
		return false;
	if (FloatMakePositive(src1.y - src2.y) > tolerance)
		return false;
	return (FloatMakePositive(src1.z - src2.z) <= tolerance);
}

/*
extern void AngleQuaternion( DegreeEuler const &angles, Quaternion &qt );
extern void QuaternionAngles( Quaternion const &q, DegreeEuler &angles );
extern void QuaternionVectorsFLU( Quaternion const &q, Vector *pForward, Vector *pLeft, Vector *pUp );
*/

inline Quaternion::Quaternion(DegreeEuler const& angles)
{
	RadianEuler radians(angles);
	AngleQuaternion(radians, *this);
}

inline DegreeEuler::DegreeEuler(RadianEuler const& angles)
{
	Init(VEC_RAD2DEG(angles.x), VEC_RAD2DEG(angles.y), VEC_RAD2DEG(angles.z));
}

inline RadianEuler::RadianEuler(DegreeEuler const& angles)
{
	Init(VEC_DEG2RAD(angles.x), VEC_DEG2RAD(angles.y), VEC_DEG2RAD(angles.z));
}

inline DegreeEuler::DegreeEuler(Quaternion const& q)
{
	RadianEuler radians(q);
	Init(VEC_RAD2DEG(radians.x), VEC_RAD2DEG(radians.y), VEC_RAD2DEG(radians.z));
}

inline bool DegreeEuler::IsValid() const
{
	return IsFinite(x) && IsFinite(y) && IsFinite(z);
}

inline void DegreeEuler::Invalidate()
{
	//#ifdef VECTOR_PARANOIA
	x = y = z = VEC_T_NAN;
	//#endif
}


//-----------------------------------------------------------------------------
// Array access
//-----------------------------------------------------------------------------
inline vec_t& DegreeEuler::operator[](int i)
{
	Assert((i >= 0) && (i < 3));
	return ((vec_t*)this)[i];
}

inline vec_t DegreeEuler::operator[](int i) const
{
	Assert((i >= 0) && (i < 3));
	return ((vec_t*)this)[i];
}



//-----------------------------------------------------------------------------
// Degree Euler QAngle pitch, yaw, roll
//-----------------------------------------------------------------------------
class QAngleByValue;

class QAngle
{
public:
	// Members
	vec_t x, y, z;

	// Construction/destruction
	QAngle(void);
	QAngle(vec_t X, vec_t Y, vec_t Z);
#ifndef _PS3
	//	QAngle(RadianEuler const &angles);	// evil auto type promotion!!!
#endif

	// Allow pass-by-value
	operator QAngleByValue& () { return *((QAngleByValue*)(this)); }
	operator const QAngleByValue& () const { return *((const QAngleByValue*)(this)); }

	// Initialization
	void Init(vec_t ix = 0.0f, vec_t iy = 0.0f, vec_t iz = 0.0f);
	void Random(vec_t minVal, vec_t maxVal);

	// Got any nasty NAN's?
	bool IsValid() const;
	void Invalidate();

	// array access...
	vec_t operator[](int i) const;
	vec_t& operator[](int i);

	// Base address...
	vec_t* Base();
	vec_t const* Base() const;

	// equality
	bool operator==(const QAngle& v) const;
	bool operator!=(const QAngle& v) const;

	// arithmetic operations
	QAngle& operator+=(const QAngle& v);
	QAngle& operator-=(const QAngle& v);
	QAngle& operator*=(float s);
	QAngle& operator/=(float s);

	// Get the vector's magnitude.
	vec_t	Length() const;
	vec_t	LengthSqr() const;

	// negate the QAngle components
	//void	Negate(); 

	// No assignment operators either...
	QAngle& operator=(const QAngle& src);

	void Normalize();
	void NormalizePositive();

	inline struct matrix3x4_t ToMatrix() const;
	inline Quaternion ToQuaternion() const;

#ifndef VECTOR_NO_SLOW_OPERATIONS
	// copy constructors

	// arithmetic operations
	QAngle	operator-(void) const;

	QAngle	operator+(const QAngle& v) const;
	QAngle	operator-(const QAngle& v) const;
	QAngle	operator*(float fl) const;
	QAngle	operator/(float fl) const;
#else

private:
	// No copy constructors allowed if we're in optimal mode
	QAngle(const QAngle& vOther);

#endif
};

// Zero the object -- necessary for CNetworkVar and possibly other cases.
inline void EnsureValidValue(QAngle& x) { x.Init(); }

//-----------------------------------------------------------------------------
// Allows us to specifically pass the vector by value when we need to
//-----------------------------------------------------------------------------
class QAngleByValue : public QAngle
{
public:
	// Construction/destruction:
	QAngleByValue(void) : QAngle() {}
	QAngleByValue(vec_t X, vec_t Y, vec_t Z) : QAngle(X, Y, Z) {}
	QAngleByValue(const QAngleByValue& vOther) { *this = vOther; }
};


inline void VectorAdd(const QAngle& a, const QAngle& b, QAngle& result)
{
	CHECK_VALID(a);
	CHECK_VALID(b);
	result.x = a.x + b.x;
	result.y = a.y + b.y;
	result.z = a.z + b.z;
}

inline void VectorMA(const QAngle& start, float scale, const QAngle& direction, QAngle& dest)
{
	CHECK_VALID(start);
	CHECK_VALID(direction);
	dest.x = start.x + scale * direction.x;
	dest.y = start.y + scale * direction.y;
	dest.z = start.z + scale * direction.z;
}


//-----------------------------------------------------------------------------
// constructors
//-----------------------------------------------------------------------------
inline QAngle::QAngle(void)
{
#ifdef _DEBUG
#ifdef VECTOR_PARANOIA
	// Initialize to NAN to catch errors
	x = y = z = VEC_T_NAN;
#endif
#endif
}

inline QAngle::QAngle(vec_t X, vec_t Y, vec_t Z)
{
	x = X; y = Y; z = Z;
	CHECK_VALID(*this);
}


//-----------------------------------------------------------------------------
// initialization
//-----------------------------------------------------------------------------
inline void QAngle::Init(vec_t ix, vec_t iy, vec_t iz)
{
	x = ix; y = iy; z = iz;
	CHECK_VALID(*this);
}


extern float AngleNormalize(float angle);
extern float AngleNormalizePositive(float angle);

inline void QAngle::Normalize()
{
	x = AngleNormalize(x);
	y = AngleNormalize(y);
	z = AngleNormalize(z);
}

inline void QAngle::NormalizePositive()
{
	x = AngleNormalizePositive(x);
	y = AngleNormalizePositive(y);
	z = AngleNormalizePositive(z);
}


#if !defined(__SPU__)
inline void QAngle::Random(vec_t minVal, vec_t maxVal)
{
	x = RandomFloat(minVal, maxVal);
	y = RandomFloat(minVal, maxVal);
	z = RandomFloat(minVal, maxVal);
	CHECK_VALID(*this);
}
#endif

#ifndef VECTOR_NO_SLOW_OPERATIONS

#if !defined(__SPU__)
inline QAngle RandomAngle(float minVal, float maxVal)
{
	Vector3D random;
	random.Random(minVal, maxVal);
	QAngle ret(random.x, random.y, random.z);
	return ret;
}
#endif

#endif

inline RadianEuler::RadianEuler(QAngle const& angles)
{
	Init(
		angles.z * 3.14159265358979323846f / 180.f,
		angles.x * 3.14159265358979323846f / 180.f,
		angles.y * 3.14159265358979323846f / 180.f);
}

inline DegreeEuler::DegreeEuler(QAngle const& angles)
{
	Init(angles.z, angles.x, angles.y);
}

inline QAngle RadianEuler::ToQAngle(void) const
{
	return QAngle(VEC_RAD2DEG(y), VEC_RAD2DEG(z), VEC_RAD2DEG(x));
}

inline QAngle DegreeEuler::ToQAngle() const
{
	return QAngle(y, z, x);
}


//-----------------------------------------------------------------------------
// assignment
//-----------------------------------------------------------------------------
inline QAngle& QAngle::operator=(const QAngle& vOther)
{
	CHECK_VALID(vOther);
	x = vOther.x; y = vOther.y; z = vOther.z;
	return *this;
}


//-----------------------------------------------------------------------------
// Array access
//-----------------------------------------------------------------------------
inline vec_t& QAngle::operator[](int i)
{
	Assert((i >= 0) && (i < 3));
	return ((vec_t*)this)[i];
}

inline vec_t QAngle::operator[](int i) const
{
	Assert((i >= 0) && (i < 3));
	return ((vec_t*)this)[i];
}


//-----------------------------------------------------------------------------
// Base address...
//-----------------------------------------------------------------------------
inline vec_t* QAngle::Base()
{
	return (vec_t*)this;
}

inline vec_t const* QAngle::Base() const
{
	return (vec_t const*)this;
}


//-----------------------------------------------------------------------------
// IsValid?
//-----------------------------------------------------------------------------
inline bool QAngle::IsValid() const
{
	return IsFinite(x) && IsFinite(y) && IsFinite(z);
}

//-----------------------------------------------------------------------------
// Invalidate
//-----------------------------------------------------------------------------

inline void QAngle::Invalidate()
{
	//#ifdef _DEBUG
	//#ifdef VECTOR_PARANOIA
	x = y = z = VEC_T_NAN;
	//#endif
	//#endif
}

//-----------------------------------------------------------------------------
// comparison
//-----------------------------------------------------------------------------
inline bool QAngle::operator==(const QAngle& src) const
{
	CHECK_VALID(src);
	CHECK_VALID(*this);
	return (src.x == x) && (src.y == y) && (src.z == z);
}

inline bool QAngle::operator!=(const QAngle& src) const
{
	CHECK_VALID(src);
	CHECK_VALID(*this);
	return (src.x != x) || (src.y != y) || (src.z != z);
}


//-----------------------------------------------------------------------------
// Copy
//-----------------------------------------------------------------------------
inline void VectorCopy(const QAngle& src, QAngle& dst)
{
	CHECK_VALID(src);
	dst.x = src.x;
	dst.y = src.y;
	dst.z = src.z;
}


//-----------------------------------------------------------------------------
// standard math operations
//-----------------------------------------------------------------------------
inline QAngle& QAngle::operator+=(const QAngle& v)
{
	CHECK_VALID(*this);
	CHECK_VALID(v);
	x += v.x; y += v.y; z += v.z;
	return *this;
}

inline QAngle& QAngle::operator-=(const QAngle& v)
{
	CHECK_VALID(*this);
	CHECK_VALID(v);
	x -= v.x; y -= v.y; z -= v.z;
	return *this;
}

inline QAngle& QAngle::operator*=(float fl)
{
	x *= fl;
	y *= fl;
	z *= fl;
	CHECK_VALID(*this);
	return *this;
}

inline QAngle& QAngle::operator/=(float fl)
{
	Assert(fl != 0.0f);
	float oofl = 1.0f / fl;
	x *= oofl;
	y *= oofl;
	z *= oofl;
	CHECK_VALID(*this);
	return *this;
}


//-----------------------------------------------------------------------------
// length
//-----------------------------------------------------------------------------
inline vec_t QAngle::Length() const
{
	CHECK_VALID(*this);
	return (vec_t)FastSqrt(LengthSqr());
}


inline vec_t QAngle::LengthSqr() const
{
	CHECK_VALID(*this);
	return x * x + y * y + z * z;
}


//-----------------------------------------------------------------------------
// Vector equality with tolerance
//-----------------------------------------------------------------------------
inline bool QAnglesAreEqual(const QAngle& src1, const QAngle& src2, float tolerance = 0.0f)
{
	if (FloatMakePositive(src1.x - src2.x) > tolerance)
		return false;
	if (FloatMakePositive(src1.y - src2.y) > tolerance)
		return false;
	return (FloatMakePositive(src1.z - src2.z) <= tolerance);
}


//-----------------------------------------------------------------------------
// arithmetic operations (SLOW!!)
//-----------------------------------------------------------------------------
#ifndef VECTOR_NO_SLOW_OPERATIONS

inline QAngle QAngle::operator-(void) const
{
	QAngle ret(-x, -y, -z);
	return ret;
}

inline QAngle QAngle::operator+(const QAngle& v) const
{
	QAngle res;
	res.x = x + v.x;
	res.y = y + v.y;
	res.z = z + v.z;
	return res;
}

inline QAngle QAngle::operator-(const QAngle& v) const
{
	QAngle res;
	res.x = x - v.x;
	res.y = y - v.y;
	res.z = z - v.z;
	return res;
}

inline QAngle QAngle::operator*(float fl) const
{
	QAngle res;
	res.x = x * fl;
	res.y = y * fl;
	res.z = z * fl;
	return res;
}

inline QAngle QAngle::operator/(float fl) const
{
	QAngle res;
	res.x = x / fl;
	res.y = y / fl;
	res.z = z / fl;
	return res;
}

inline QAngle operator*(float fl, const QAngle& v)
{
	QAngle ret(v * fl);
	return ret;
}

#endif // VECTOR_NO_SLOW_OPERATIONS


//-----------------------------------------------------------------------------
// NOTE: These are not completely correct.  The representations are not equivalent
// unless the QAngle represents a rotational impulse along a coordinate axis (x,y,z)
inline void QAngleToAngularImpulse(const QAngle& angles, AngularImpulse& impulse)
{
	impulse.x = angles.z;
	impulse.y = angles.x;
	impulse.z = angles.y;
}

inline void AngularImpulseToQAngle(const AngularImpulse& impulse, QAngle& angles)
{
	angles.x = impulse.y;
	angles.y = impulse.z;
	angles.z = impulse.x;
}

inline QAngle Quaternion::ToQAngle() const
{
	extern void QuaternionAngles(const Quaternion & q, QAngle & angles);

	QAngle anglesOut;
	QuaternionAngles(*this, anglesOut);
	return anglesOut;
}

#if !defined( _X360 ) && !defined( _PS3 )

FORCEINLINE vec_t InvRSquared(const float* v)
{
	return 1.0f / MAX(1.0f, v[0] * v[0] + v[1] * v[1] + v[2] * v[2]);
}

FORCEINLINE vec_t InvRSquared(const Vector3D& v)
{
	return InvRSquared(v.Base());
}

#else

// call directly
#if defined(__SPU__)
FORCEINLINE float _VMX_InvRSquared(Vector& v)
#else
FORCEINLINE float _VMX_InvRSquared(const Vector& v)
#endif
{
#if !defined (_PS3)
	XMVECTOR xmV = XMVector3ReciprocalLength(XMLoadVector3(v.Base()));
	xmV = XMVector3Dot(xmV, xmV);
	return xmV.x;
#else	//!_PS3
	vector_float_union vRet;
	vec_float4 v0, v1, vIn, vOut;
	vector unsigned char permMask;
	v0 = vec_ld(0, v.Base());
	permMask = vec_lvsl(0, v.Base());
	v1 = vec_ld(11, v.Base());
	vIn = vec_perm(v0, v1, permMask);
	vOut = vec_madd(vIn, vIn, _VEC_ZEROF);
	vec_float4 vTmp = vec_sld(vIn, vIn, 4);
	vec_float4 vTmp2 = vec_sld(vIn, vIn, 8);
	vOut = vec_madd(vTmp, vTmp, vOut);
	vOut = vec_madd(vTmp2, vTmp2, vOut);
	vOut = vec_re(vec_add(vOut, _VEC_EPSILONF));
	vec_st(vOut, 0, &vRet.vf);
	float ret = vRet.f[0];
	return ret;
#endif	//!_PS3
}

#define InvRSquared(x) _VMX_InvRSquared(x)

#endif // _X360

#if !defined( _X360 ) && !defined( _PS3 )

// FIXME: Change this back to a #define once we get rid of the vec_t version
float VectorNormalize(Vector3D& v);

// FIXME: Obsolete version of VectorNormalize, once we remove all the friggin float*s
FORCEINLINE float VectorNormalize(float* v)
{
	return VectorNormalize(*(reinterpret_cast<Vector3D*>(v)));
}

#else
#if !defined( _PS3 )
// modified version of Microsoft's XMVector3Length
// microsoft's version will return INF for very small vectors
// e.g. 	Vector vTest(7.98555446e-20,-6.85012984e-21,0); VectorNormalize( vTest );
// so we clamp to epsilon instead of checking for zero
XMFINLINE XMVECTOR XMVector3Length_Fixed
(
	FXMVECTOR V
)
{
	// Returns a QNaN on infinite vectors.
	static CONST XMVECTOR g_fl4SmallVectorEpsilon = { 1e-24f,1e-24f,1e-24f,1e-24f };

	XMVECTOR D;
	XMVECTOR Rsq;
	XMVECTOR Rcp;
	XMVECTOR Zero;
	XMVECTOR RT;
	XMVECTOR Result;
	XMVECTOR Length;
	XMVECTOR H;

	H = __vspltisw(1);
	D = __vmsum3fp(V, V);
	H = __vcfsx(H, 1);
	Rsq = __vrsqrtefp(D);
	RT = __vmulfp(D, H);
	Rcp = __vmulfp(Rsq, Rsq);
	H = __vnmsubfp(RT, Rcp, H);
	Rsq = __vmaddfp(Rsq, H, Rsq);
	Zero = __vspltisw(0);
	Result = __vcmpgefp(g_fl4SmallVectorEpsilon, D);
	Length = __vmulfp(D, Rsq);
	Result = __vsel(Length, Zero, Result);

	return Result;
}
#endif

// call directly
FORCEINLINE float _VMX_VectorNormalize(Vector& vec)
{
#if !defined _PS3
	float mag = XMVector3Length_Fixed(XMLoadVector3(vec.Base())).x;
	float den = 1.f / (mag + FLT_EPSILON);
	vec.x *= den;
	vec.y *= den;
	vec.z *= den;
	return mag;
#else	// !_PS3
	vec_float4 vIn;
	vec_float4 v0, v1;
	vector unsigned char permMask;
	v0 = vec_ld(0, vec.Base());
	permMask = vec_lvsl(0, vec.Base());
	v1 = vec_ld(11, vec.Base());
	vIn = vec_perm(v0, v1, permMask);
	float mag = vmathV3Length((VmathVector3*)&vIn);
	float den = 1.f / (mag + FLT_EPSILON);
	vec.x *= den;
	vec.y *= den;
	vec.z *= den;
	return mag;
#endif	// !_PS3
}
// FIXME: Change this back to a #define once we get rid of the vec_t version
FORCEINLINE float VectorNormalize(Vector& v)
{
	return _VMX_VectorNormalize(v);
}
// FIXME: Obsolete version of VectorNormalize, once we remove all the friggin float*s
FORCEINLINE float VectorNormalize(float* pV)
{
	return _VMX_VectorNormalize(*(reinterpret_cast<Vector*>(pV)));
}

#endif // _X360

#if !defined( _X360 ) && !defined( _PS3 )
FORCEINLINE void VectorNormalizeFast(Vector3D& vec)
{
	float ool = FastRSqrt(FLT_EPSILON + vec.x * vec.x + vec.y * vec.y + vec.z * vec.z);

	vec.x *= ool;
	vec.y *= ool;
	vec.z *= ool;
}
#else

// call directly
FORCEINLINE void VectorNormalizeFast(Vector& vec)
{
#if !defined (_PS3)
	XMVECTOR xmV = XMVector3LengthEst(XMLoadVector3(vec.Base()));
	float den = 1.f / (xmV.x + FLT_EPSILON);
	vec.x *= den;
	vec.y *= den;
	vec.z *= den;
#else	// !_PS3
	vector_float_union vVec;

	vec_float4 vIn, vOut, vOOLen, vDot;

	// load
	vec_float4 v0, v1;
	vector unsigned char permMask;
	v0 = vec_ld(0, vec.Base());
	permMask = vec_lvsl(0, vec.Base());
	v1 = vec_ld(11, vec.Base());
	vIn = vec_perm(v0, v1, permMask);

	// vec.vec
	vOut = vec_madd(vIn, vIn, _VEC_ZEROF);
	vec_float4 vTmp = vec_sld(vIn, vIn, 4);
	vec_float4 vTmp2 = vec_sld(vIn, vIn, 8);
	vOut = vec_madd(vTmp, vTmp, vOut);
	vOut = vec_madd(vTmp2, vTmp2, vOut);

	// splat dot to all 
	vDot = vec_splat(vOut, 0);

	vOOLen = vec_rsqrte(vec_add(vDot, _VEC_EPSILONF));

	// vec * 1.0/sqrt(vec.vec)
	vOut = vec_madd(vIn, vOOLen, _VEC_ZEROF);

	// store
	vec_st(vOut, 0, &vVec.vf);

	// store vec
	vec.x = vVec.f[0];
	vec.y = vVec.f[1];
	vec.z = vVec.f[2];

#endif	// !_PS3
}

#endif // _X360

inline vec_t Vector3D::NormalizeInPlace()
{
	return VectorNormalize(*this);
}

inline vec_t Vector3D::NormalizeInPlaceSafe(const Vector3D& vFallback)
{
	float flLength = VectorNormalize(*this);
	if (flLength == 0.0f)
	{
		*this = vFallback;
	}
	return flLength;
}

inline Vector3D Vector3D::Normalized() const
{
	Vector3D norm = *this;
	VectorNormalize(norm);
	return norm;
}


inline Vector3D Vector3D::NormalizedSafe(const Vector3D& vFallback)const
{
	Vector3D vNorm = *this;
	float flLength = VectorNormalize(vNorm);
	return (flLength != 0.0f) ? vNorm : vFallback;
}


inline bool Vector3D::IsLengthGreaterThan(float val) const
{
	return LengthSqr() > val * val;
}

inline bool Vector3D::IsLengthLessThan(float val) const
{
	return LengthSqr() < val * val;
}


inline const Vector3D ScaleVector(const Vector3D& a, const Vector3D& b)
{
	return Vector3D(a.x * b.x, a.y * b.y, a.z * b.z);
}



inline const Quaternion Exp(const Vector3D& v)
{
	float theta = v.Length();
	if (theta < 0.001f)
	{
		// limit case, cos(theta)       ~= 1 - theta^2/2 + theta^4/24
		//             sin(theta)/theta ~= 1 - theta^2/6 + theta^4/120
		float theta2_2 = theta * theta * 0.5f, theta4_24 = theta2_2 * theta2_2 * (1.0f / 6.0f);
		float k = 1.0f - theta2_2 * (1.0f / 3.0f) + theta4_24 * 0.05f;
		return Quaternion(k * v.x, k * v.y, k * v.z, 1 - theta2_2 + theta4_24);
	}
	else
	{
		float k = sinf(theta) / theta;
		return Quaternion(k * v.x, k * v.y, k * v.z, cosf(theta));
	}
}


inline const Vector3D QuaternionLog(const Quaternion& q)
{
	Vector3D axis = q.ImaginaryPart();
	float sinTheta = axis.Length(), factor;
	if (sinTheta > 0.001f)
	{
		// there's some substantial rotation; if w < 0, it's an over-180-degree rotation (in real space)
		float theta = asinf(MIN(sinTheta, 1.0f));
		factor = (q.w < 0.0f ? M_PI_F - theta : theta) / sinTheta;
	}
	else
	{
		// ArcSin[x]/x = 1 + x^2/6 + x^4 * 3/40 + o( x^5 )
		float sinTheta2 = sinTheta * sinTheta;
		float sinTheta4 = sinTheta2 * sinTheta2;
		factor = (1 + sinTheta2 * (1.0f / 6.0f) + sinTheta4 * (3.0f / 40.0f));
		if (q.w < 0)
		{
			factor = -factor; // because the axis of rotation is not defined, we'll just consider this rotation to be close enough to identity
		}
	}
	return axis * factor;
}



inline float Snap(float a, float flSnap)
{
	return floorf(a / flSnap + 0.5f) * flSnap;
}

inline  const Vector3D Snap(const Vector3D& a, float flSnap)
{
	return Vector3D(Snap(a.x, flSnap), Snap(a.y, flSnap), Snap(a.z, flSnap));
}

#endif