//////////////////////////////////////////////////////////////////////////////
//
// Asura_PS3_VectorMaths.h
//
// Copyright (c) 2005 Rebellion. All rights reserved.
//
// Description:
//
//     aligned vector types for PS3 (PPU and SPU)
//
//////////////////////////////////////////////////////////////////////////////

// B E G I N   F I L E   W R A P P E R ///////////////////////////////////////

#ifdef ASURA_USE_PRAGMA_ONCE
    #pragma once
#endif

#ifndef ___ASURA_PS3_VECTORMATHS_H___
#define ___ASURA_PS3_VECTORMATHS_H___ 1

////////////////////////////////////////////////////////////////////////////////
// Include files...
////////////////////////////////////////////////////////////////////////////////

// Asura Platform-independent
#include "Core/Asura_HalfFloat.h"
#include "Debugging/Asura_Logfile.h"

////////////////////////////////////////////////////////////////////////////////
// Forward class declarations...
////////////////////////////////////////////////////////////////////////////////

class Asura_AMatrix3x3;

////////////////////////////////////////////////////////////////////////////////
// Constant Defines...
////////////////////////////////////////////////////////////////////////////////

#define PERM_X 0x00010203
#define PERM_Y 0x04050607
#define PERM_Z 0x08090a0b
#define PERM_W 0x0c0d0e0f
#define PERM_A 0x10111213
#define PERM_B 0x14151617
#define PERM_C 0x18191a1b
#define PERM_D 0x1c1d1e1f
#define PERM_XYZA (vec_uchar16)(vec_uint4){ PERM_X, PERM_Y, PERM_Z, PERM_A }
#define PERM_ZXYW (vec_uchar16)(vec_uint4){ PERM_Z, PERM_X, PERM_Y, PERM_W }
#define PERM_YZXW (vec_uchar16)(vec_uint4){ PERM_Y, PERM_Z, PERM_X, PERM_W }
#define PERM_YZAB (vec_uchar16)(vec_uint4){ PERM_Y, PERM_Z, PERM_A, PERM_B }
#define PERM_ZABC (vec_uchar16)(vec_uint4){ PERM_Z, PERM_A, PERM_B, PERM_C }
#define PERM_XYAW (vec_uchar16)(vec_uint4){ PERM_X, PERM_Y, PERM_A, PERM_W }
#define PERM_XAZW (vec_uchar16)(vec_uint4){ PERM_X, PERM_A, PERM_Z, PERM_W }

#define PERM_ZBWX ((vec_uchar16)(vec_uint4){ PERM_Z, PERM_B, PERM_W, PERM_X })
#define PERM_XCYX ((vec_uchar16)(vec_uint4){ PERM_X, PERM_C, PERM_Y, PERM_X })
#define PERM_XYAB ((vec_uchar16)(vec_uint4){ PERM_X, PERM_Y, PERM_A, PERM_B })
#define PERM_ZWCD ((vec_uchar16)(vec_uint4){ PERM_Z, PERM_W, PERM_C, PERM_D })
#define PERM_XZBX ((vec_uchar16)(vec_uint4){ PERM_X, PERM_Z, PERM_B, PERM_X })     
#define PERM_CXXX ((vec_uchar16)(vec_uint4){ PERM_C, PERM_X, PERM_X, PERM_X })
#define PERM_YAXX ((vec_uchar16)(vec_uint4){ PERM_Y, PERM_A, PERM_X, PERM_X })
#define PERM_XAZC ((vec_uchar16)(vec_uint4){ PERM_X, PERM_A, PERM_Z, PERM_C })
#define PERM_YXWZ ((vec_uchar16)(vec_uint4){ PERM_Y, PERM_X, PERM_W, PERM_Z })
#define PERM_YBWD ((vec_uchar16)(vec_uint4){ PERM_Y, PERM_B, PERM_W, PERM_D })
#define PERM_XYCX ((vec_uchar16)(vec_uint4){ PERM_X, PERM_Y, PERM_C, PERM_X })
#define PERM_YCXY ((vec_uchar16)(vec_uint4){ PERM_Y, PERM_C, PERM_X, PERM_Y })
#define PERM_CXYC ((vec_uchar16)(vec_uint4){ PERM_C, PERM_X, PERM_Y, PERM_C })
#define PERM_ZAYX ((vec_uchar16)(vec_uint4){ PERM_Z, PERM_A, PERM_Y, PERM_X })
#define PERM_BZXX ((vec_uchar16)(vec_uint4){ PERM_B, PERM_Z, PERM_X, PERM_X })
#define PERM_XZYA ((vec_uchar16)(vec_uint4){ PERM_X, PERM_Z, PERM_Y, PERM_A })
#define PERM_ZXXB ((vec_uchar16)(vec_uint4){ PERM_Z, PERM_X, PERM_X, PERM_B })
#define PERM_YXXC ((vec_uchar16)(vec_uint4){ PERM_Y, PERM_X, PERM_X, PERM_C })
#define PERM_BBYX ((vec_uchar16)(vec_uint4){ PERM_B, PERM_B, PERM_Y, PERM_X })

////////////////////////////////////////////////////////////////////////////////
// Class declaration...
////////////////////////////////////////////////////////////////////////////////

// Temporary define for toggling platform implementation. 
// Platforms with their own implementation should define this 
#if defined(ASURA_USE_PLATFORM_VECTORS)

// Asura_AVector3...
#define Asura_AVector3_InParam		const Asura_AVector3&
#define Asura_AVector3_OutParam		Asura_AVector3&
#define Asura_AVector3_InOutParam	Asura_AVector3&

// Asura_AMatrix3x3...
#define Asura_AMatrix3x3_InParam	const Asura_AMatrix3x3&
#define Asura_AMatrix3x3_OutParam	Asura_AMatrix3x3&
#define Asura_AMatrix3x3_InOutParam	Asura_AMatrix3x3&

// Asura_AVector4...
#define Asura_AVector4_InParam		const Asura_AVector4&
#define Asura_AVector4_OutParam		Asura_AVector4&
#define Asura_AVector4_InOutParam	Asura_AVector4&

// Asura_AMatrix4x4...
#define Asura_AMatrix4x4_InParam	const Asura_AMatrix4x4&
#define Asura_AMatrix4x4_OutParam	Asura_AMatrix4x4&
#define Asura_AMatrix4x4_InOutParam	Asura_AMatrix4x4&

class Asura_AVector3
{

public:
				Asura_AVector3() {}
	inline 		Asura_AVector3( float fX, float fY, float fZ );
	inline 		Asura_AVector3( const vector float xVector3 );
	inline		Asura_AVector3( const Asura_Vector_3& xVec3 );

	inline void	operator = ( const Asura_Vector_3& xVec3 );
	inline void	operator = ( const Asura_AVector3& xVec3 );

	inline void		LoadAligned( const float* pfFloats );
	inline void		LoadUnaligned( const float*	pfFloats );

	inline void		StoreAligned( float* pfFloats )	const;
	inline void		StoreUnaligned(	float* pfFloats	) const;

	inline void		StoreAlignedHalf( Asura_HalfFloat* pxFloats ) const;

	// assignment operators
	inline void	operator += ( Asura_AVector3_InParam xV );
	inline void	operator -= ( Asura_AVector3_InParam xV );
	inline void	operator *= ( const float f );
	inline void	operator /= ( const float f );
	inline void	operator *= ( Asura_AVector3_InParam xV );

	// unary operators
	inline Asura_AVector3	operator -  () const;

	// binary operators
	inline Asura_AVector3	operator  + ( Asura_AVector3_InParam xV ) const;
	inline Asura_AVector3	operator  - ( Asura_AVector3_InParam xV ) const;
	inline Asura_AVector3	operator  * ( const float f ) const;
	inline Asura_AVector3	operator  / ( const float f ) const;
	inline Asura_AVector3	operator  ^ ( Asura_AVector3_InParam v ) const;
	inline float			operator  * ( Asura_AVector3_InParam v ) const;
	inline bool				operator == ( Asura_AVector3_InParam v ) const;

	inline bool ApproximatelyEquals	    ( Asura_AVector3_InParam v, const float epsilon = 1.0e-6f) const;
	inline bool				operator != ( Asura_AVector3_InParam v) const;
	inline bool				operator  > ( Asura_AVector3_InParam v) const;
	inline bool				operator  < ( Asura_AVector3_InParam v) const;
	inline bool				operator >= ( Asura_AVector3_InParam v) const;
	inline bool				operator <= ( Asura_AVector3_InParam v) const;

	inline void		Max(Asura_AVector3_InParam xV1, Asura_AVector3_InParam xV2);
	inline void		Min(Asura_AVector3_InParam xV1, Asura_AVector3_InParam xV2);

	inline float	MagnitudeSquared(void) const;
	inline float	Magnitude(void) const;

	inline void		Normalise(void);
	inline float	NormaliseReturningMagnitude(void);
	inline void		ApproximateNormalise(void);
	inline float	ApproximateNormaliseReturningMagnitude(void);

	inline void		Zero(void) ;

	inline void		SelectGE( Asura_AVector3_InParam xLHS, Asura_AVector3_InParam xRHS, Asura_AVector3_InParam xTrue, Asura_AVector3_InParam xFalse );
	inline void		SelectGT( Asura_AVector3_InParam xLHS, Asura_AVector3_InParam xRHS, Asura_AVector3_InParam xTrue, Asura_AVector3_InParam xFalse );
	inline void		SelectLE( Asura_AVector3_InParam xLHS, Asura_AVector3_InParam xRHS, Asura_AVector3_InParam xTrue, Asura_AVector3_InParam xFalse );
	inline void		SelectLT( Asura_AVector3_InParam xLHS, Asura_AVector3_InParam xRHS, Asura_AVector3_InParam xTrue, Asura_AVector3_InParam xFalse );
	inline void		SelectEQ( Asura_AVector3_InParam xLHS, Asura_AVector3_InParam xRHS, Asura_AVector3_InParam xTrue, Asura_AVector3_InParam xFalse );
	inline void		SelectNE( Asura_AVector3_InParam xLHS, Asura_AVector3_InParam xRHS, Asura_AVector3_InParam xTrue, Asura_AVector3_InParam xFalse );

	inline float	GetX(void)			const	{ return m_xVector3[0]; }
	inline float	GetY(void)			const	{ return m_xVector3[1]; }
	inline float	GetZ(void)			const	{ return m_xVector3[2]; }
	inline float	GetElement(u_int u) const	{ return m_xVector3[u]; }
    inline const float* GetElements( void ) const { return &m_xVector3[0]; }

	inline void		SetX(float fX)				{ m_xVector3[0] = fX; }
	inline void		SetY(float fY)				{ m_xVector3[1] = fY; }
	inline void		SetZ(float fZ)				{ m_xVector3[2] = fZ; }
	inline void		SetElement(u_int u, float f){ m_xVector3[u] = f; }

	inline void		Lerp(Asura_AVector3_InParam xTo, const float fWeight);
	inline void		Lerp(Asura_AVector3_InParam xFrom, Asura_AVector3_InParam xTo, const float fWeight);

	inline void		Output() const;

private:
	vector float m_xVector3;

	friend class Asura_AMatrix3x3;
	friend class Asura_AVector4;
	friend class Asura_AMatrix4x4;
};

inline Asura_AVector3::Asura_AVector3( float fX, float fY, float fZ )
{
	m_xVector3 = ( vector float ){ fX, fY, fZ, 0.0f };
}

inline Asura_AVector3::Asura_AVector3 ( const Asura_Vector_3& xVec3 )
{
	m_xVector3 = (vector float){ xVec3.x, xVec3.y, xVec3.z, 0.0f };
}

inline Asura_AVector3::Asura_AVector3( const vector float xVector3 )
{
	m_xVector3 = xVector3;
}

inline void Asura_AVector3::operator = ( const Asura_Vector_3& xVec3 ) 
{
	m_xVector3 = (vector float){ xVec3.x, xVec3.y, xVec3.z, 0.0f };

}

inline void Asura_AVector3::operator = ( const Asura_AVector3& xVec3 )
{
	m_xVector3 = xVec3.m_xVector3;
}

inline void	Asura_AVector3::LoadAligned( const float* pfFloats )
{
	m_xVector3 = vec_ld( 0, const_cast<float*>(pfFloats) ); //BW stupid sony headers

	// this loads four floats from pfFloats, so we need to set the last
	// one to zero
    m_xVector3[3] = 0.0f;
}

inline void	Asura_AVector3::LoadUnaligned( const float*	pfFloats )
{
	m_xVector3 = (vector float){ pfFloats[0], pfFloats[1], pfFloats[2], 0.0f };

}

inline void	Asura_AVector3::StoreAligned( float* pfFloats )	const 
{
	// vec_st will replace four values in *pfFloats, which may not be appropriate,
	// but should be okay for aligned data
	// for the moment just use unaligned style
	// (MP - 08.01.2008) TODO: May be better to use vec_ste here.
	pfFloats[0] = m_xVector3[0];
	pfFloats[1] = m_xVector3[1];
	pfFloats[2] = m_xVector3[2];
	//vec_st( m_xVector3, 0, pfFloats );
}

inline void	Asura_AVector3::StoreUnaligned(	float* pfFloats	) const	
{
	pfFloats[0] = m_xVector3[0];
	pfFloats[1] = m_xVector3[1];
	pfFloats[2] = m_xVector3[2];
}

inline void Asura_AVector3::StoreAlignedHalf( Asura_HalfFloat* pxHalves ) const
{
	pxHalves[0] = m_xVector3[0];
	pxHalves[1] = m_xVector3[1];
	pxHalves[2] = m_xVector3[2];
}

inline void Asura_AVector3::operator += ( Asura_AVector3_InParam xV )
{
	m_xVector3 = vec_add( m_xVector3, xV.m_xVector3 );
}

inline void Asura_AVector3::operator -= ( Asura_AVector3_InParam xV )
{
	m_xVector3 = vec_sub( m_xVector3, xV.m_xVector3 );
}

inline void Asura_AVector3::operator *= ( const float f )
{
	vector float	xScale	= { f, f, f, 0.0f };

	m_xVector3 = vec_madd( m_xVector3, xScale, (vector float)(0.0f) );
}

inline void Asura_AVector3::operator /= ( const float f )
{
	float fInv = 1.0f / f;
	vector float	xScale	= { fInv, fInv, fInv, 0.0f };

	m_xVector3 = vec_madd( m_xVector3, xScale, (vector float)(0.0f) );
}

inline void Asura_AVector3::operator *= (	Asura_AVector3_InParam xV ) 
{
	m_xVector3 = vec_madd( m_xVector3, xV.m_xVector3, (vector float)(0.0f) );
}

inline Asura_AVector3 Asura_AVector3::operator - () const
{
	return Asura_AVector3( negatef4(m_xVector3) );
}

inline Asura_AVector3 Asura_AVector3::operator + ( Asura_AVector3_InParam xV) const
{
	return Asura_AVector3( vec_add( m_xVector3, xV.m_xVector3 ) );
}

inline Asura_AVector3 Asura_AVector3::operator - ( Asura_AVector3_InParam xV) const
{
	return Asura_AVector3( vec_sub( m_xVector3, xV.m_xVector3 ) );
}

inline Asura_AVector3 Asura_AVector3::operator * ( const float f) const
{
	vector float	xScale	= { f, f, f, 0.0f };

	return Asura_AVector3( vec_madd( m_xVector3, xScale, (vector float)(0.0f) ) );
}

inline Asura_AVector3 Asura_AVector3::operator / ( float f) const
{
	float fInv = 1.0f / f;
	vector float	xScale	= { fInv, fInv, fInv, 0.0f };

	return Asura_AVector3( vec_madd( m_xVector3, xScale, (vector float)(0.0f) ) );
}                                               

inline Asura_AVector3 Asura_AVector3::operator ^ ( Asura_AVector3_InParam v) const
{
	// doing the cross product with vector operations is perhaps slightly slower in the speed tests
	// than casting to floats and using floating point arithmetic instead, but the inlining should 
	// mean that in real code the compiler can interleave e.g. the load instructions to achieve
	// speed improvements
	// cross product ( 1.y*2.z - 1.z*2.y , 1.z*2.x - 1.x*2.z , 1.x*2.y - 1.y * 2.x )
	// produce two vectors (1.y*2.z,1.z*2.x,1.x*2.y) and (1.z*2.y,1.x*2.z,1.y*2.x) then subtract
	vector unsigned char xLeftShift = {0x04,0x05,0x06,0x07,0x08,0x09,0x0A,0x0B,0x00,0x01,0x02,0x03,0x0C,0x0D,0x0E,0x0F};
	vector float xVec1 = vec_perm(  m_xVector3, (vector float)(0.0f), xLeftShift );
	vector unsigned char xRightShift = {0x08,0x09,0x0A,0x0B,0x00,0x01,0x02,0x03,0x04,0x05,0x06,0x07,0x0C,0x0D,0x0E,0x0F};
	vector float xVec2 = vec_madd( xVec1, vec_perm(v.m_xVector3, (vector float)(0.0f), xRightShift), (vector float)(0.0f) );
	vector float xVec4 = vec_perm( m_xVector3, (vector float)(0.0f), xRightShift );
	vector float xVec3 = vec_perm( v.m_xVector3, (vector float)(0.0f), xLeftShift);
	vector float result = vec_nmsub( xVec4, xVec3, xVec2 );
	return Asura_AVector3( result );
}

inline float Asura_AVector3::operator * ( Asura_AVector3_InParam v) const
{
	// for short operations like this it is quicker to use components than vector operations
	float *pf = (float*)&m_xVector3;
	float *pf2 = (float*)&v.m_xVector3;
	return ( pf[0]*pf2[0] + pf[1]*pf2[1] + pf[2]*pf2[2] );
}

inline bool Asura_AVector3::operator == ( Asura_AVector3_InParam v) const
{
	return vec_all_eq( m_xVector3, v.m_xVector3 );
}

inline bool  Asura_AVector3::ApproximatelyEquals( Asura_AVector3_InParam v, const float epsilon) const
{
	const vector float delta = vec_sub( m_xVector3, v.m_xVector3 );
	const float* pf = &delta[0];
	return ( (pf[0]*pf[0]+pf[1]*pf[1]+pf[2]*pf[2]) < epsilon );
}

inline bool Asura_AVector3::operator != ( Asura_AVector3_InParam v) const
{
	return vec_any_ne( m_xVector3, v.m_xVector3 );
}

inline bool	Asura_AVector3::operator > ( Asura_AVector3_InParam v) const
{
	// greater than and less than functions are quite tricky because there is a
	// fourth component in the vector which is always 0, so vector comparisons will 
	// not return correct answers
	// instead, we have to use a float style comparison
	const float* pf = &m_xVector3[0];
	const float* pf2 = &v.m_xVector3[0];
	return ( pf[0]>pf2[0] && pf[1]>pf2[1] && pf[2]>pf2[2] );
}

inline bool	Asura_AVector3::operator < ( Asura_AVector3_InParam v) const
{
	// see note for > operator above
	const float* pf = &m_xVector3[0];
	const float* pf2 = &v.m_xVector3[0];
	return ( pf[0]<pf2[0] && pf[1]<pf2[1] && pf[2]<pf2[2] );
}

inline bool	Asura_AVector3::operator >=	( Asura_AVector3_InParam v) const
{
	return vec_all_ge( m_xVector3, v.m_xVector3 );
}

inline bool	Asura_AVector3::operator <=	( Asura_AVector3_InParam v) const
{
	return vec_all_le( m_xVector3, v.m_xVector3 );
}

inline void Asura_AVector3::Max(Asura_AVector3_InParam xV1, Asura_AVector3_InParam xV2)
{
	m_xVector3 = vec_max(xV1.m_xVector3,xV2.m_xVector3);
}

inline void Asura_AVector3::Min(Asura_AVector3_InParam xV1, Asura_AVector3_InParam xV2)
{
	m_xVector3 = vec_min(xV1.m_xVector3,xV2.m_xVector3);
}

inline float Asura_AVector3::MagnitudeSquared(void) const
{   
	const float *pf = (const float*) &m_xVector3;
	return pf[0]*pf[0]+pf[1]*pf[1]+pf[2]*pf[2];
}

inline float Asura_AVector3::Magnitude(void) const
{   
	const float *pf = &m_xVector3[0];
	return Asura_Maths::Sqrt( pf[0]*pf[0]+pf[1]*pf[1]+pf[2]*pf[2] );
}

inline void Asura_AVector3::Normalise(void)
{
	// the PS3 includes the SIMD instruction rsqrtf4 which will calculate the reciprocal square root for all four components of the vector
	// at once, enabling ~ 3x faster normalisation. However, it seems to only be accurate up to around 7d.p. which causes it to fail
	// some of the unit tests == operators. For the moment it is safer to leave this function using the reference implementation, until
	// it is determined whether the errors will break any code
#ifdef ASURA_PS3_FASTNORMALISE
	vector float result;
	result = vec_madd( m_xVector3, m_xVector3, (vec_float4)(0.0f) );
	result = vec_madd( vec_sld( m_xVector3, m_xVector3, 4 ), vec_sld( m_xVector3, m_xVector3, 4 ), result );
	result = vec_madd( vec_sld( m_xVector3, m_xVector3, 8 ), vec_sld( m_xVector3, m_xVector3, 8 ), result );
	result = vec_splat( result, 0 );
	result = rsqrtf4( result );
	m_xVector3 = vec_madd(m_xVector3,result,(vector float)(0.0f) );
#else
	float fMag = Asura_Maths::Sqrt( MagnitudeSquared() );
	if( fMag != 0.0f )
	{
		float fInvMag = 1.0f/fMag;
		float* pf = &m_xVector3[0];
		pf[0] *= fInvMag;
		pf[1] *= fInvMag;
		pf[2] *= fInvMag;
	}
#endif
}

inline float Asura_AVector3::NormaliseReturningMagnitude(void)
{
#ifdef ASURA_PS3_FASTNORMALISE
	vector float result;
	result = vec_madd( m_xVector3, m_xVector3, (vec_float4)(0.0f) );
	result = vec_madd( vec_sld( m_xVector3, m_xVector3, 4 ), vec_sld( m_xVector3, m_xVector3, 4 ), result );
	result = vec_madd( vec_sld( m_xVector3, m_xVector3, 8 ), vec_sld( m_xVector3, m_xVector3, 8 ), result );
	result = vec_splat( result, 0 );
	result = rsqrtf4( result );
	m_xVector3 = vec_madd(m_xVector3,result,(vector float)(0.0f) );
	return 1.0f/(*(reinterpret_cast<float*>(&result) ) );
#else
	float fMag = Asura_Maths::Sqrt( MagnitudeSquared() );
	if( fMag != 0.0f )
	{
		float fInvMag = 1.0f/fMag;
		float* pf = &m_xVector3[0];
		pf[0] *= fInvMag;
		pf[1] *= fInvMag;
		pf[2] *= fInvMag;

		return fMag;
	}
	return 0.0f;
#endif
}

inline void Asura_AVector3::ApproximateNormalise(void)
{
	vector float result;
	result = vec_madd( m_xVector3, m_xVector3, (vec_float4)(0.0f) );
	result = vec_madd( vec_sld( m_xVector3, m_xVector3, 4 ), vec_sld( m_xVector3, m_xVector3, 4 ), result );
	result = vec_madd( vec_sld( m_xVector3, m_xVector3, 8 ), vec_sld( m_xVector3, m_xVector3, 8 ), result );
	result = vec_splat( result, 0 );
	result = vec_rsqrte( result );
	m_xVector3 = vec_madd(m_xVector3,result,(vector float)(0.0f) );
}       

inline float Asura_AVector3::ApproximateNormaliseReturningMagnitude(void)
{
	vector float result;
	result = vec_madd( m_xVector3, m_xVector3, (vec_float4)(0.0f) );
	result = vec_madd( vec_sld( m_xVector3, m_xVector3, 4 ), vec_sld( m_xVector3, m_xVector3, 4 ), result );
	result = vec_madd( vec_sld( m_xVector3, m_xVector3, 8 ), vec_sld( m_xVector3, m_xVector3, 8 ), result );
	result = vec_splat( result, 0 );
	result = vec_rsqrte( result );
	m_xVector3 = vec_madd(m_xVector3,result,(vector float)(0.0f) );
	return 1.0f/(*(reinterpret_cast<float*>(&result) ) );
}

inline void Asura_AVector3::Zero(void) 
{
	m_xVector3 = (vector float)(0.0f);
}

inline void	Asura_AVector3::SelectGE( Asura_AVector3_InParam xLHS, Asura_AVector3_InParam xRHS, Asura_AVector3_InParam xTrue, Asura_AVector3_InParam xFalse )
{
	const vec_bint4 xSelection = vec_cmpge( xLHS.m_xVector3, xRHS.m_xVector3 );
	m_xVector3 = vec_sel( xFalse.m_xVector3, xTrue.m_xVector3, xSelection );
}

inline void	Asura_AVector3::SelectGT( Asura_AVector3_InParam xLHS, Asura_AVector3_InParam xRHS, Asura_AVector3_InParam xTrue, Asura_AVector3_InParam xFalse )
{
	const vec_bint4 xSelection = vec_cmpgt( xLHS.m_xVector3, xRHS.m_xVector3 );
	m_xVector3 = vec_sel( xFalse.m_xVector3, xTrue.m_xVector3, xSelection );
}

inline void	Asura_AVector3::SelectLE( Asura_AVector3_InParam xLHS, Asura_AVector3_InParam xRHS, Asura_AVector3_InParam xTrue, Asura_AVector3_InParam xFalse )
{
	const vec_bint4 xSelection = vec_cmple( xLHS.m_xVector3, xRHS.m_xVector3 );
	m_xVector3 = vec_sel( xFalse.m_xVector3, xTrue.m_xVector3, xSelection );
}

inline void	Asura_AVector3::SelectLT( Asura_AVector3_InParam xLHS, Asura_AVector3_InParam xRHS, Asura_AVector3_InParam xTrue, Asura_AVector3_InParam xFalse )
{
	const vec_bint4 xSelection = vec_cmplt( xLHS.m_xVector3, xRHS.m_xVector3 );
	m_xVector3 = vec_sel( xFalse.m_xVector3, xTrue.m_xVector3, xSelection );
}

inline void	Asura_AVector3::SelectEQ( Asura_AVector3_InParam xLHS, Asura_AVector3_InParam xRHS, Asura_AVector3_InParam xTrue, Asura_AVector3_InParam xFalse )
{
	const vec_bint4 xSelection = vec_cmpeq( xLHS.m_xVector3, xRHS.m_xVector3 );
	m_xVector3 = vec_sel( xFalse.m_xVector3, xTrue.m_xVector3, xSelection );
}

inline void	Asura_AVector3::SelectNE( Asura_AVector3_InParam xLHS, Asura_AVector3_InParam xRHS, Asura_AVector3_InParam xTrue, Asura_AVector3_InParam xFalse )
{
	const vec_bint4 xSelection = vec_cmpeq( xLHS.m_xVector3, xRHS.m_xVector3 );
	m_xVector3 = vec_sel( xTrue.m_xVector3, xFalse.m_xVector3, xSelection );
}

inline void Asura_AVector3::Lerp(Asura_AVector3_InParam xTo, const float fWeight)
{
	Asura_AVector3 xTemp( vec_sub(xTo.m_xVector3, m_xVector3) );
	xTemp*=fWeight;
	m_xVector3 = vec_add( xTemp.m_xVector3, m_xVector3 );
}

inline void Asura_AVector3::Lerp(Asura_AVector3_InParam xFrom, Asura_AVector3_InParam xTo, const float fWeight)
{
	Asura_AVector3 xTemp( vec_sub(xTo.m_xVector3, xFrom.m_xVector3) );
	xTemp*=fWeight;
	m_xVector3 = vec_add( xTemp.m_xVector3, xFrom.m_xVector3 );
}

inline void Asura_AVector3::Output() const
{
	Asura_OutputToDebugger_Retail( ("Int 3Vector %f %f %f\n", GetX(), GetY(), GetZ() ) );
}

inline void Asura_Vector_3::operator = ( const Asura_AVector3& xV4)
{
	x = xV4.GetX();
	y = xV4.GetY();
	z = xV4.GetZ();
}

// explicit constructor
inline Asura_Vector_3::Asura_Vector_3( const Asura_AVector3& xV4 )
{
	x = xV4.GetX();
	y = xV4.GetY();
	z = xV4.GetZ();
}

// Note that in Asura_Matrix_3x3 the data is stored as 3 rows. For efficiency when using 
// vector instructions it's better for us to store 3 columns. Since the old 3x3 matrix
// use e11 to e33 I'm keeping the convention of starting named columns at 1 rather than 0.
class Asura_AMatrix3x3
{
	Asura_AVector3 m_xCol0;
	Asura_AVector3 m_xCol1;
	Asura_AVector3 m_xCol2;

	friend class Asura_AVector3;

public:

	inline Asura_AMatrix3x3( void ) { }
	inline Asura_AMatrix3x3(	float e11, float e12, float e13,
								float e21, float e22, float e23,
								float e31, float e32, float e33 );
	inline Asura_AMatrix3x3(	const Asura_AVector3& xCol0, const Asura_AVector3& xCol1, const Asura_AVector3& xCol2 );
	inline Asura_AMatrix3x3( const Asura_Matrix_3x3& xMat3 ) { FromMat3( xMat3 ); }

	inline void				ToMat3(Asura_Matrix_3x3& xMat3) const;
	inline void				FromMat3(const Asura_Matrix_3x3& xMat3);

	inline Asura_AVector3	GetColumn0(void) const { return m_xCol0; }
	inline Asura_AVector3	GetColumn1(void) const { return m_xCol1; }
	inline Asura_AVector3	GetColumn2(void) const { return m_xCol2; }
	inline Asura_AVector3	GetColumn(u_int u) const { return *(&m_xCol0 + u); }
	inline float			GetElement(u_int uRow, u_int uCol) const { return GetColumn(uCol).GetElement(uRow); }
	
	const float* GetElements( void ) const { return m_xCol0.GetElements(); } //pointer to the data with padding for alignment, i.e. {11,12,13,<pad>,21,22,23,<pad>,31,32,33,<pad>}

	inline void				SetColumn0(Asura_AVector3_InParam xCol0) { m_xCol0 = xCol0; }
	inline void				SetColumn1(Asura_AVector3_InParam xCol1) { m_xCol1 = xCol1; }
	inline void				SetColumn2(Asura_AVector3_InParam xCol2) { m_xCol2 = xCol2; }

	// Tilde operator to make a skew-symmetric matrix from a vector
	inline void				Tilde( Asura_AVector3_InParam xVector );
	// Try to make the matrix orthogonal by forcing the axes represented by the matrix
	// columns to be a) unit-length and b) perpendicular.
	inline void				Orthogonalise( void );
	// matrix-matrix per component addition - used in integration
	inline void				Add( Asura_AMatrix3x3_InParam xM );
	inline Asura_AVector3	Apply( Asura_AVector3_InParam xInput) const;
	inline void				Apply( Asura_AMatrix3x3_InParam xInput, Asura_AMatrix3x3_OutParam xOutput ) const;
	inline Asura_AVector3	ApplyTranspose( Asura_AVector3_InParam xVec ) const;
	inline void				Transpose(void);
	inline void				Identity(void);
	inline void				Zero(void);

	inline Asura_AMatrix3x3	operator  * ( const float s) const;
	inline bool				operator == ( Asura_AMatrix3x3_InParam xMat) const;
	inline bool				operator != ( Asura_AMatrix3x3_InParam xMat) const;
};

inline Asura_AMatrix3x3::Asura_AMatrix3x3(	float e11, float e12, float e13,
											float e21, float e22, float e23,
											float e31, float e32, float e33 )
{
	m_xCol0 = Asura_AVector3(e11, e21, e31);
	m_xCol1 = Asura_AVector3(e12, e22, e32);
	m_xCol2 = Asura_AVector3(e13, e23, e33);
}

inline Asura_AMatrix3x3::Asura_AMatrix3x3( const Asura_AVector3& xCol0, const Asura_AVector3& xCol1, const Asura_AVector3& xCol2 )
{
	m_xCol0 = xCol0;
	m_xCol1 = xCol1;
	m_xCol2 = xCol2;
}

inline void Asura_AMatrix3x3::ToMat3(Asura_Matrix_3x3& xMat3) const
{
	xMat3.e11 = GetElement(0,0);
	xMat3.e12 = GetElement(0,1);
	xMat3.e13 = GetElement(0,2);
	xMat3.e21 = GetElement(1,0);
	xMat3.e22 = GetElement(1,1);
	xMat3.e23 = GetElement(1,2);
	xMat3.e31 = GetElement(2,0);
	xMat3.e32 = GetElement(2,1);
	xMat3.e33 = GetElement(2,2);
}

inline void Asura_AMatrix3x3::FromMat3(const Asura_Matrix_3x3& xMat3)
{
	m_xCol0 = Asura_AVector3(xMat3.e11, xMat3.e21, xMat3.e31);
	m_xCol1 = Asura_AVector3(xMat3.e12, xMat3.e22, xMat3.e32);
	m_xCol2 = Asura_AVector3(xMat3.e13, xMat3.e23, xMat3.e33);
}

inline void Asura_AMatrix3x3::Tilde( Asura_AVector3_InParam xVector )
{
	m_xCol0 = Asura_AVector3( 0.0f,				+xVector.GetZ(),	-xVector.GetY() );
	m_xCol1 = Asura_AVector3( -xVector.GetZ(),	0.0f,				+xVector.GetX() );
	m_xCol2 = Asura_AVector3( +xVector.GetY(),	-xVector.GetX(),	0.0f			);
}

inline void Asura_AMatrix3x3::Orthogonalise( void )
{
	// normalise "down" column
	m_xCol1.Normalise();

	// make "forwards" column perpendicular to this, and normalise
	m_xCol2 -= m_xCol1 * ( m_xCol2 * m_xCol1 );
	m_xCol2.Normalise();

	// and now the "right" column must be perpendicular to both of these: x = + y^z
	m_xCol0 = m_xCol1 ^ m_xCol2;
}

inline void Asura_AMatrix3x3::Add( Asura_AMatrix3x3_InParam xM ) 
{
	m_xCol0 += xM.m_xCol0;
	m_xCol1 += xM.m_xCol1;
	m_xCol2 += xM.m_xCol2;
}

inline Asura_AVector3 Asura_AMatrix3x3::Apply( Asura_AVector3_InParam xInput) const
{
	vector float xResult;
	vector float xX, xY, xZ;
	xX = vec_splat( xInput.m_xVector3, 0 );
	xY = vec_splat( xInput.m_xVector3, 1 );
	xZ = vec_splat( xInput.m_xVector3, 2 );
	xResult = vec_madd( GetColumn0().m_xVector3, xX, (vector float)(0.0f) );
	xResult = vec_madd( GetColumn1().m_xVector3, xY, xResult );
	xResult = vec_madd( GetColumn2().m_xVector3, xZ, xResult );
	return Asura_AVector3( xResult );
}

inline void Asura_AMatrix3x3::Apply( Asura_AMatrix3x3_InParam xInput, Asura_AMatrix3x3_OutParam xOutput ) const
{
	vector float xInCol0 = xInput.m_xCol0.m_xVector3;
	vector float xOutCol0 = vec_madd( m_xCol0.m_xVector3, vec_splat( xInCol0, 0 ), (vector float)(0.0f) );
	xOutCol0 = vec_madd( m_xCol1.m_xVector3, vec_splat( xInCol0, 1 ), xOutCol0 );
	xOutCol0 = vec_madd( m_xCol2.m_xVector3, vec_splat( xInCol0, 2 ), xOutCol0 );

	vector float xInCol1 = xInput.m_xCol1.m_xVector3;
	vector float xOutCol1 = vec_madd( m_xCol0.m_xVector3, vec_splat( xInCol1, 0 ), (vector float)(0.0f) );
	xOutCol1 = vec_madd( m_xCol1.m_xVector3, vec_splat( xInCol1, 1 ), xOutCol1 );
	xOutCol1 = vec_madd( m_xCol2.m_xVector3, vec_splat( xInCol1, 2 ), xOutCol1 );

	vector float xInCol2 = xInput.m_xCol2.m_xVector3;
	vector float xOutCol2 = vec_madd( m_xCol0.m_xVector3, vec_splat( xInCol2, 0 ), (vector float)(0.0f) );
	xOutCol2 = vec_madd( m_xCol1.m_xVector3, vec_splat( xInCol2, 1 ), xOutCol2 );
	xOutCol2 = vec_madd( m_xCol2.m_xVector3, vec_splat( xInCol2, 2 ), xOutCol2 );

	xOutput.m_xCol0 = xOutCol0;
	xOutput.m_xCol1 = xOutCol1;
	xOutput.m_xCol2 = xOutCol2;
}

inline void Asura_AMatrix3x3::Transpose(void)
{
	float f = m_xCol0.GetY();
	m_xCol0.SetY( m_xCol1.GetX() );
	m_xCol1.SetX( f );

	f = m_xCol0.GetZ();
	m_xCol0.SetZ( m_xCol2.GetX() );
	m_xCol2.SetX( f );

	f = m_xCol1.GetZ();
	m_xCol1.SetZ( m_xCol2.GetY() );
	m_xCol2.SetY( f );
}

inline Asura_AVector3 Asura_AMatrix3x3::ApplyTranspose( Asura_AVector3_InParam xVec ) const
{
	// transpose turns rows into columns, so instead of transposing and multiplying by rows,
	// we can just multiply by columns, which is quicker
	return(Asura_AVector3((xVec * m_xCol0), (xVec * m_xCol1), (xVec * m_xCol2)));
}

inline void Asura_AMatrix3x3::Identity(void)
{
	m_xCol0 = Asura_AVector3( 1.0f, 0.0f, 0.0f );
	m_xCol1 = Asura_AVector3( 0.0f, 1.0f, 0.0f );
	m_xCol2 = Asura_AVector3( 0.0f, 0.0f, 1.0f );
}

inline Asura_AMatrix3x3 Asura_AMatrix3x3::operator * ( const float s) const
{
	return Asura_AMatrix3x3(m_xCol0 * s, m_xCol1 * s, m_xCol2 * s);
}

inline void Asura_AMatrix3x3::Zero(void)
{
	m_xCol0.Zero();
	m_xCol1.Zero();
	m_xCol2.Zero();
}

inline bool Asura_AMatrix3x3::operator == ( Asura_AMatrix3x3_InParam xMat) const
{
	return	m_xCol0 == xMat.m_xCol0 &&
		m_xCol1 == xMat.m_xCol1 &&
		m_xCol2 == xMat.m_xCol2;
}
inline bool Asura_AMatrix3x3::operator != ( Asura_AMatrix3x3_InParam xMat) const
{
	return ! (*this==xMat);
}

inline void Asura_Quat::ToMatrix(Asura_AMatrix3x3_OutParam m) const
{
	/* Implementation from Graphics Gems IV */

	float s = 2.0f;///Magnitude();      

	float xs = x*s;
	float ys = y*s;
	float zs = z*s;

	float wx = w*xs;
	float wy = w*ys;
	float wz = w*zs;

	float xx = x*xs;
	float xy = x*ys;
	float xz = x*zs;

	float yy = y*ys;
	float yz = y*zs;
	float zz = z*zs;

	m.SetColumn0(Asura_AVector3( 1.0f - (yy+zz), xy + wz, xz - wy ));
	m.SetColumn1(Asura_AVector3( xy - wz, 1.0f - (xx + zz), yz + wx ));
	m.SetColumn2(Asura_AVector3( xz + wy, yz - wx, 1.0f - (xx+yy) ));
}

inline void Asura_Quat::FromMatrix(Asura_AMatrix3x3_InParam xM4)
{
	Asura_Matrix_3x3 xM3;
	xM4.ToMat3(xM3);
	FromMatrix(xM3);
}

inline Asura_AVector3 Asura_Quat::RotateVector(Asura_AVector3_InParam p) const
{
	//-----------------------------------------------------
	// 
	// This function rotates a vector p by the quaternion q
	// The formula is:
	//                  p' = (q)(p)(q*)
	//
	// but we calculate it using a vector formula which is an
	// expansion of this that Kev worked out.
	// for q = (v,w)
	// v' = 2w(v^p) + (2(v.p)v + (w^2-v^2)p
	//
	//--------------------------------------------------------

	const Asura_AVector3 v(x,y,z); // the vector part of q
	Asura_AVector3 section1, section2, section3;

	float vMagSq = v.MagnitudeSquared();

	section1 = v^p;
	section1 *= (2*w);

	section2 = v;
	section2 *= (2*(v*p));

	section3 = p;
	section3 *= (w*w - vMagSq);

	return section1+section2+section3;
}

// creates a quaternion from an angle and axis of rotation
inline void Asura_Quat::FromRotation(Asura_Angle angle, Asura_AVector3_InParam axis)
{
	Asura_Assert(Asura_Maths::Abs(axis.MagnitudeSquared() - 1.0f) < 0.0001f, "axis should be normalised");

	angle *= 0.5f;
	float s;
	Asura_Maths::SinCos(angle, s, w);
	x = axis.GetX() * s;
	y = axis.GetY() * s;
	z = axis.GetZ() * s;
}

class Asura_SoA_Float
{
	float m_afValues[4];

public:

	Asura_SoA_Float() {}
	inline Asura_SoA_Float( float fSplat );
	inline Asura_SoA_Float( float f0, float f1, float f2, float f3 );
	inline Asura_SoA_Float( const Asura_SoA_Float& xRHS );

	inline Asura_SoA_Float&	operator  = ( const Asura_SoA_Float& xRHS );
	inline Asura_SoA_Float	operator  + ( const Asura_SoA_Float& xRHS ) const;
	inline Asura_SoA_Float	operator  * ( const Asura_SoA_Float& xRHS ) const;
	inline bool				operator == ( const Asura_SoA_Float& xRHS ) const;
	inline bool				operator != ( const Asura_SoA_Float& xRHS ) const;

	friend class Asura_SoA_Vector_3;
};

inline Asura_SoA_Float::Asura_SoA_Float( float fSplat )
{
	m_afValues[0] = fSplat;
	m_afValues[1] = fSplat;
	m_afValues[2] = fSplat;
	m_afValues[3] = fSplat;
}

inline Asura_SoA_Float::Asura_SoA_Float( float f0, float f1, float f2, float f3 )
{
	m_afValues[0] = f0;
	m_afValues[1] = f1;
	m_afValues[2] = f2;
	m_afValues[3] = f3;
}

inline Asura_SoA_Float& Asura_SoA_Float::operator = ( const Asura_SoA_Float& xRHS )
{
	m_afValues[0] = xRHS.m_afValues[0];
	m_afValues[1] = xRHS.m_afValues[1];
	m_afValues[2] = xRHS.m_afValues[2];
	m_afValues[3] = xRHS.m_afValues[3];
	return *this;
}

inline Asura_SoA_Float::Asura_SoA_Float( const Asura_SoA_Float& xRHS )
{
	m_afValues[0] = xRHS.m_afValues[0];
	m_afValues[1] = xRHS.m_afValues[1];
	m_afValues[2] = xRHS.m_afValues[2];
	m_afValues[3] = xRHS.m_afValues[3];
}

inline Asura_SoA_Float Asura_SoA_Float::operator + ( const Asura_SoA_Float& xRHS ) const
{
	return Asura_SoA_Float( m_afValues[0] + xRHS.m_afValues[0],
		m_afValues[1] + xRHS.m_afValues[1],
		m_afValues[2] + xRHS.m_afValues[2],
		m_afValues[3] + xRHS.m_afValues[3] );
}

inline Asura_SoA_Float Asura_SoA_Float::operator * ( const Asura_SoA_Float& xRHS ) const
{
	return Asura_SoA_Float( m_afValues[0] * xRHS.m_afValues[0],
		m_afValues[1] * xRHS.m_afValues[1],
		m_afValues[2] * xRHS.m_afValues[2],
		m_afValues[3] * xRHS.m_afValues[3] );
}

inline bool Asura_SoA_Float::operator == ( const Asura_SoA_Float& xRHS ) const
{
	return	m_afValues[0] == xRHS.m_afValues[0] && 
		m_afValues[1] == xRHS.m_afValues[1] && 
		m_afValues[2] == xRHS.m_afValues[2] && 
		m_afValues[3] == xRHS.m_afValues[3];
}

inline bool Asura_SoA_Float::operator != ( const Asura_SoA_Float& xRHS ) const
{
	return	m_afValues[0] != xRHS.m_afValues[0] || 
		m_afValues[1] != xRHS.m_afValues[1] || 
		m_afValues[2] != xRHS.m_afValues[2] || 
		m_afValues[3] != xRHS.m_afValues[3];
}

class Asura_SoA_Vector_3
{
	Asura_SoA_Float m_xX;
	Asura_SoA_Float m_xY;
	Asura_SoA_Float m_xZ;

public:

	Asura_SoA_Vector_3() {}
	inline Asura_SoA_Vector_3( const Asura_SoA_Float& xX, const Asura_SoA_Float& xY, const Asura_SoA_Float& xZ );
	inline Asura_SoA_Vector_3( const Asura_SoA_Vector_3& xRHS);

	// create new SoA floats from members
	inline Asura_SoA_Float GetXXXX() const { return m_xX; }
	inline Asura_SoA_Float GetYYYY() const { return m_xY; }
	inline Asura_SoA_Float GetZZZZ() const { return m_xZ; }


	inline void		operator = ( const Asura_SoA_Vector_3& xRHS );
	inline Asura_SoA_Float			operator * ( const Asura_SoA_Vector_3& xRHS ) const;// dot product
	inline const Asura_SoA_Vector_3	operator * ( const Asura_AMatrix3x3& xMat ) const; // transform by matrix
};

inline Asura_SoA_Vector_3::Asura_SoA_Vector_3( const Asura_SoA_Float& xX, const Asura_SoA_Float& xY, const Asura_SoA_Float& xZ ) 
{
	m_xX = xX;
	m_xY = xY;
	m_xZ = xZ;
}

inline Asura_SoA_Vector_3::Asura_SoA_Vector_3( const Asura_SoA_Vector_3& xRHS)
{
	m_xX = xRHS.m_xX;
	m_xY = xRHS.m_xY;
	m_xZ = xRHS.m_xZ;
}

inline void Asura_SoA_Vector_3::operator = ( const Asura_SoA_Vector_3& xRHS )
{
	m_xX = xRHS.m_xX;
	m_xY = xRHS.m_xY;
	m_xZ = xRHS.m_xZ;
}

inline Asura_SoA_Float Asura_SoA_Vector_3::operator * ( const Asura_SoA_Vector_3& xRHS ) const
{
	//		return	m_xX * xRHS.m_xX + 
	//				m_xY * xRHS.m_xY + 
	//				m_xZ * xRHS.m_xZ;
	float fOut0, fOut1, fOut2, fOut3;

	float fX0 = m_xX.m_afValues[0];
	float fX1 = xRHS.m_xX.m_afValues[0];
	float fY0 = m_xY.m_afValues[0];
	float fY1 = xRHS.m_xY.m_afValues[0];
	float fZ0 = m_xZ.m_afValues[0];
	float fZ1 = xRHS.m_xZ.m_afValues[0];
	float fX2 = m_xX.m_afValues[1];
	float fX3 = xRHS.m_xX.m_afValues[1];
	float fY2 = m_xY.m_afValues[1];
	float fY3 = xRHS.m_xY.m_afValues[1];
	float fZ2 = m_xZ.m_afValues[1];
	float fZ3 = xRHS.m_xZ.m_afValues[1];

	fOut0 = fX0 * fX1 + fY0 * fY1 + fZ0 * fZ1;
	fOut1 = fX2 * fX3 + fY2 * fY3 + fZ2 * fZ3;

	fX0 = m_xX.m_afValues[2];
	fX1 = xRHS.m_xX.m_afValues[2];
	fY0 = m_xY.m_afValues[2];
	fY1 = xRHS.m_xY.m_afValues[2];
	fZ0 = m_xZ.m_afValues[2];
	fZ1 = xRHS.m_xZ.m_afValues[2];
	fX2 = m_xX.m_afValues[3];
	fX3 = xRHS.m_xX.m_afValues[3];
	fY2 = m_xY.m_afValues[3];
	fY3 = xRHS.m_xY.m_afValues[3];
	fZ2 = m_xZ.m_afValues[3];
	fZ3 = xRHS.m_xZ.m_afValues[3];

	fOut2 = fX0 * fX1 + fY0 * fY1 + fZ0 * fZ1;
	fOut3 = fX2 * fX3 + fY2 * fY3 + fZ2 * fZ3;

	return Asura_SoA_Float( fOut0, fOut1, fOut2, fOut3 );
}

inline const Asura_SoA_Vector_3 Asura_SoA_Vector_3::operator * ( const Asura_AMatrix3x3& xMat ) const
{
	const Asura_AVector3 xCol0 = xMat.GetColumn0();
	// splat M11 to all elements of xM11, etc.
	Asura_SoA_Float xM11(xCol0.GetX());
	Asura_SoA_Float xM21(xCol0.GetY());
	Asura_SoA_Float xM31(xCol0.GetZ());

	Asura_SoA_Float xOutX = m_xX * xM11;
	Asura_SoA_Float xOutY = m_xX * xM21;
	Asura_SoA_Float xOutZ = m_xX * xM31;

	const Asura_AVector3 xCol1 = xMat.GetColumn1();
	Asura_SoA_Float xM12(xCol1.GetX());
	Asura_SoA_Float xM22(xCol1.GetY());
	Asura_SoA_Float xM32(xCol1.GetZ());

	xOutX = m_xY * xM12 + xOutX;
	xOutY = m_xY * xM22 + xOutY;
	xOutZ = m_xY * xM32 + xOutZ;

	const Asura_AVector3 xCol2 = xMat.GetColumn2();
	Asura_SoA_Float xM13(xCol2.GetX());
	Asura_SoA_Float xM23(xCol2.GetY());
	Asura_SoA_Float xM33(xCol2.GetZ());

	xOutX = m_xZ * xM13 + xOutX;
	xOutY = m_xZ * xM23 + xOutY;
	xOutZ = m_xZ * xM33 + xOutZ;

	return Asura_SoA_Vector_3( xOutX, xOutY, xOutZ );
}

//////////////////////////////////////////////////////////////////////////
// Asura_AVector4
//////////////////////////////////////////////////////////////////////////
class Asura_AVector4
{
public:
	// Constructors...
	Asura_AVector4(void) {}
	inline Asura_AVector4(const float fX, const float fY, const float fZ, const float fW);
	inline Asura_AVector4(const Asura_AVector4& xOther);
	inline Asura_AVector4(const Asura_AVector3& xOther, const float fW);
	inline Asura_AVector4(const Asura_Vector_3& xOther, const float fW);
	inline Asura_AVector4( const vector float xVector4 );

	// Zero functions...
	inline void ZeroW0(void);
	inline void ZeroW1(void);

	// Selection functions...
	inline void SelectGE(Asura_AVector4_InParam xLHS, Asura_AVector4_InParam xRHS, Asura_AVector4_InParam xTrue, Asura_AVector4_InParam xFalse);
	inline void SelectGT(Asura_AVector4_InParam xLHS, Asura_AVector4_InParam xRHS, Asura_AVector4_InParam xTrue, Asura_AVector4_InParam xFalse);
	inline void SelectLE(Asura_AVector4_InParam xLHS, Asura_AVector4_InParam xRHS, Asura_AVector4_InParam xTrue, Asura_AVector4_InParam xFalse);
	inline void SelectLT(Asura_AVector4_InParam xLHS, Asura_AVector4_InParam xRHS, Asura_AVector4_InParam xTrue, Asura_AVector4_InParam xFalse);
	inline void SelectEQ(Asura_AVector4_InParam xLHS, Asura_AVector4_InParam xRHS, Asura_AVector4_InParam xTrue, Asura_AVector4_InParam xFalse);
	inline void SelectNE(Asura_AVector4_InParam xLHS, Asura_AVector4_InParam xRHS, Asura_AVector4_InParam xTrue, Asura_AVector4_InParam xFalse);

	// Load/store functions...
	inline void LoadAligned(const float* pfFloats);
	inline void LoadUnaligned(const float* pfFloats);
	inline void StoreAligned(float* pfFloats) const;
	inline void StoreUnaligned(float* pfFloats) const;
	inline void StoreAlignedHalf(Asura_HalfFloat* pxHalves) const;

	// Manipulation functions...
	inline Asura_AVector4	ComponentwiseProduct(Asura_AVector4_InParam xOther) const;
	inline float	DotProduct4(Asura_AVector4_InParam xOther) const;
	inline float	DotProduct3(Asura_AVector4_InParam xOther) const;
	inline float	MagnitudeSquared(void) const;
	inline float	Magnitude(void) const;
	inline void		Normalise(void);
	inline float	NormaliseReturningMagnitude(void);
	inline void		ApproximateNormalise(void);
	inline float	ApproximateNormaliseReturningMagnitude(void);

	// Assignment operators...
	inline void	operator  = ( Asura_AVector4_InParam xOther );
	inline void	operator += ( Asura_AVector4_InParam xOther );
	inline void	operator -= ( Asura_AVector4_InParam xOther );
	inline void	operator *= ( const float fScale );
	inline void	operator /= ( const float fScale );

	// Unary operators...
	inline Asura_AVector4	operator - ( void ) const;

	//Binary operators
	inline Asura_AVector4	operator + ( Asura_AVector4_InParam xOther ) const;
	inline Asura_AVector4	operator - ( Asura_AVector4_InParam xOther ) const;
	inline Asura_AVector4	operator ^ ( Asura_AVector4_InParam xOther ) const;
	inline float			operator * ( Asura_AVector4_InParam xOther ) const;
	inline Asura_AVector4	operator * ( const float fScale ) const;
	inline Asura_AVector4	operator / ( const float fScale ) const;

	// Comparison operators...
	bool					operator == (Asura_AVector4_InParam xOther) const;
	bool		ApproximatelyEquals		(Asura_AVector4_InParam xOther, const float fEpsilon = 1.0e-6f) const;
	bool					operator != (Asura_AVector4_InParam xOther) const;

	// Minimum/maximum functions...
	inline void Max(Asura_AVector4_InParam xV1, Asura_AVector4_InParam xV2) { m_xVector4 = vec_max( xV1.m_xVector4, xV2.m_xVector4 ); }
	inline void Min(Asura_AVector4_InParam xV1, Asura_AVector4_InParam xV2) { m_xVector4 = vec_min( xV1.m_xVector4, xV2.m_xVector4 ); }

	// Properties...
	inline float	GetX(void) const { return m_xVector4[0]; }
	inline float	GetY(void) const { return m_xVector4[1]; }
	inline float	GetZ(void) const { return m_xVector4[2]; }
	inline float	GetW(void) const { return m_xVector4[3]; }
	inline float	GetElement(const u_int uIndex) const;
	inline const float* GetElements( void ) const { return &m_xVector4[0]; }
	inline vector float GetVector() { return m_xVector4; }

	inline void		SetX(const float fX) { m_xVector4[0] = fX; }
	inline void		SetY(const float fY) { m_xVector4[1] = fY; }
	inline void		SetZ(const float fZ) { m_xVector4[2] = fZ; }
	inline void		SetW(const float fW) { m_xVector4[3] = fW; }
	inline void		SetElement(const u_int uIndex, const float fValue);
	inline void		SetElements(const float fX, const float fY, const float fZ, const float fW);
	inline void		SetVector(const vector float xVector) { m_xVector4 = xVector; }

private:
	vector float m_xVector4;

	friend class Asura_AMatrix4x4;
};

inline Asura_AVector4::Asura_AVector4(const float fX, const float fY, const float fZ, const float fW)
{
	m_xVector4 = (vector float){ fX, fY, fZ, fW };
}

inline Asura_AVector4::Asura_AVector4(const Asura_AVector4& xOther)
: m_xVector4( xOther.m_xVector4 )
{
}

inline Asura_AVector4::Asura_AVector4(const Asura_AVector3& xOther, const float fW)
: m_xVector4( xOther.m_xVector3 )
{
	m_xVector4[3] = fW;
}

inline Asura_AVector4::Asura_AVector4(const Asura_Vector_3& xOther, const float fW)
{
	m_xVector4 = (vector float){ xOther.GetX(), xOther.GetY(), xOther.GetZ(), fW };
}

inline Asura_AVector4::Asura_AVector4( const vector float xVector4 )
: m_xVector4( xVector4 )
{
}

inline void Asura_AVector4::ZeroW0(void)
{
	m_xVector4 = (vector float)(0.0f);
}

inline void Asura_AVector4::ZeroW1(void)
{
	m_xVector4 = (vector float){ 0.0f, 0.0f, 0.0f, 1.0f };
}

inline void Asura_AVector4::SelectGE(Asura_AVector4_InParam xLHS, Asura_AVector4_InParam xRHS, Asura_AVector4_InParam xTrue, Asura_AVector4_InParam xFalse)
{
	const vec_bint4 xSelection = vec_cmpge( xLHS.m_xVector4, xRHS.m_xVector4 );
	m_xVector4 = vec_sel( xFalse.m_xVector4, xTrue.m_xVector4, xSelection );
}

inline void Asura_AVector4::SelectGT(Asura_AVector4_InParam xLHS, Asura_AVector4_InParam xRHS, Asura_AVector4_InParam xTrue, Asura_AVector4_InParam xFalse)
{
	const vec_bint4 xSelection = vec_cmpgt( xLHS.m_xVector4, xRHS.m_xVector4 );
	m_xVector4 = vec_sel( xFalse.m_xVector4, xTrue.m_xVector4, xSelection );
}

inline void Asura_AVector4::SelectLE(Asura_AVector4_InParam xLHS, Asura_AVector4_InParam xRHS, Asura_AVector4_InParam xTrue, Asura_AVector4_InParam xFalse)
{
	const vec_bint4 xSelection = vec_cmple( xLHS.m_xVector4, xRHS.m_xVector4 );
	m_xVector4 = vec_sel( xFalse.m_xVector4, xTrue.m_xVector4, xSelection );
}

inline void Asura_AVector4::SelectLT(Asura_AVector4_InParam xLHS, Asura_AVector4_InParam xRHS, Asura_AVector4_InParam xTrue, Asura_AVector4_InParam xFalse)
{
	const vec_bint4 xSelection = vec_cmplt( xLHS.m_xVector4, xRHS.m_xVector4 );
	m_xVector4 = vec_sel( xFalse.m_xVector4, xTrue.m_xVector4, xSelection );
}

inline void Asura_AVector4::SelectEQ(Asura_AVector4_InParam xLHS, Asura_AVector4_InParam xRHS, Asura_AVector4_InParam xTrue, Asura_AVector4_InParam xFalse)
{
	const vec_bint4 xSelection = vec_cmpeq( xLHS.m_xVector4, xRHS.m_xVector4 );
	m_xVector4 = vec_sel( xFalse.m_xVector4, xTrue.m_xVector4, xSelection );
}

inline void Asura_AVector4::SelectNE(Asura_AVector4_InParam xLHS, Asura_AVector4_InParam xRHS, Asura_AVector4_InParam xTrue, Asura_AVector4_InParam xFalse)
{
	const vec_bint4 xSelection = vec_cmpeq( xLHS.m_xVector4, xRHS.m_xVector4 );
	m_xVector4 = vec_sel( xTrue.m_xVector4, xFalse.m_xVector4, xSelection );
}

inline void Asura_AVector4::LoadAligned(const float* pfFloats)
{
	m_xVector4 = vec_ld( 0, const_cast<float*>(pfFloats) );
}

inline void Asura_AVector4::LoadUnaligned(const float* pfFloats)
{
	m_xVector4 = (vector float){ pfFloats[0], pfFloats[1], pfFloats[2], pfFloats[3] };
}

inline void Asura_AVector4::StoreAligned(float* pfFloats) const
{
	vec_st( m_xVector4, 0, pfFloats );
}

inline void Asura_AVector4::StoreUnaligned(float* pfFloats) const
{
	pfFloats[0] = m_xVector4[0];
	pfFloats[1] = m_xVector4[1];
	pfFloats[2] = m_xVector4[2];	
	pfFloats[3] = m_xVector4[3];
}

inline void Asura_AVector4::StoreAlignedHalf(Asura_HalfFloat* pxHalves) const
{
	pxHalves[0] = m_xVector4[0];
	pxHalves[1] = m_xVector4[1];
	pxHalves[2] = m_xVector4[2];
	pxHalves[3] = m_xVector4[3];
}

inline Asura_AVector4 Asura_AVector4::ComponentwiseProduct(Asura_AVector4_InParam xOther) const
{
	return Asura_AVector4( vec_madd(m_xVector4, xOther.m_xVector4, (vector float)(0.0f) ) );
}

inline float Asura_AVector4::DotProduct4(Asura_AVector4_InParam xOther) const
{

// Alternative implementation using VMX intrinsic operators
// Tests show no better performance compared to the non ev
#if 0
	vector float xV1  = vec_madd(m_xVector4, xOther.m_xVector4, (vector float)(0.0f)); // result of multiplication = [A , B, C, D]
	vector float xV2 = vec_sld(  xV1, xV1, 4 ); // B , C , D , A 
	xV1 = vec_add(xV1, xV2); // A + B , B + C, C + D, A + D 
	xV2 =  vec_sld(  xV1, xV1, 8 ); // C + D, A + D, A + B, B + C 
	xV1 = vec_add( xV1, xV2); // A + B + C + D  
	return vec_extract(xV1, 0);
#endif

	const float* pfTemp = &m_xVector4[0];
	const float* pfOther = &xOther.m_xVector4[0];
	return ( pfTemp[0]*pfOther[0] + pfTemp[1]*pfOther[1] + pfTemp[2]*pfOther[2] + pfTemp[3]*pfOther[3] ) ;
}

inline float Asura_AVector4::DotProduct3(Asura_AVector4_InParam xOther) const
{
	const float* pfTemp = &m_xVector4[0];
	const float* pfOther = &xOther.m_xVector4[0];
	return ( pfTemp[0]*pfOther[0] + pfTemp[1]*pfOther[1] + pfTemp[2]*pfOther[2] ) ;
}

inline float Asura_AVector4::MagnitudeSquared(void) const
{
	const float* pfTemp = &m_xVector4[0];
	return ( pfTemp[0]*pfTemp[0] + pfTemp[1]*pfTemp[1] + pfTemp[2]*pfTemp[2] + pfTemp[3]*pfTemp[3] ) ;
}

inline float Asura_AVector4::Magnitude(void) const
{
	return(Asura_Maths::Sqrt(MagnitudeSquared()));
}

inline void Asura_AVector4::Normalise(void)
{
	const float	fMagnitude	= Magnitude();

	if( fMagnitude != 0.0f )
	{
		float fInvMag = 1.0f/fMagnitude;
		*this	*= fInvMag;
	}
}

inline float Asura_AVector4::NormaliseReturningMagnitude(void)
{
	const float	fMagnitude	= Magnitude();	

	if(fMagnitude == 0.0f)
	{
		return(0.0f);
	}
	else
	{
		float fInvMag = 1.0f/fMagnitude;
		*this	*= fInvMag;
		return (fMagnitude);
	}
}

inline void Asura_AVector4::ApproximateNormalise(void)
{
	vector float xResult;
	xResult = vec_madd( m_xVector4, m_xVector4, (vector float)(0.0f) );
	xResult = vec_madd( vec_sld( m_xVector4, m_xVector4, 4 ), vec_sld( m_xVector4, m_xVector4, 4 ), xResult );
	xResult = vec_madd( vec_sld( m_xVector4, m_xVector4, 8 ), vec_sld( m_xVector4, m_xVector4, 8 ), xResult );
	xResult = vec_madd( vec_sld( m_xVector4, m_xVector4, 12), vec_sld( m_xVector4, m_xVector4, 12 ), xResult );
	xResult = vec_splat( xResult, 0 );
	xResult = vec_rsqrte( xResult );
	m_xVector4 = vec_madd(m_xVector4,xResult, (vector float)(0.0f) );
}

inline float Asura_AVector4::ApproximateNormaliseReturningMagnitude(void)
{
	vector float xResult;
	xResult = vec_madd( m_xVector4, m_xVector4, (vector float)(0.0f) );
	xResult = vec_madd( vec_sld( m_xVector4, m_xVector4, 4 ), vec_sld( m_xVector4, m_xVector4, 4 ), xResult );
	xResult = vec_madd( vec_sld( m_xVector4, m_xVector4, 8 ), vec_sld( m_xVector4, m_xVector4, 8 ), xResult );
	xResult = vec_madd( vec_sld( m_xVector4, m_xVector4, 12), vec_sld( m_xVector4, m_xVector4, 12 ), xResult );
	xResult = vec_splat( xResult, 0 );
	xResult = vec_rsqrte( xResult );
	m_xVector4 = vec_madd(m_xVector4,xResult, (vector float)(0.0f) );
	return 1.0f / *(reinterpret_cast<float*>(&xResult) );
}

inline void Asura_AVector4::operator = (Asura_AVector4_InParam xOther)
{
	m_xVector4 = xOther.m_xVector4;
}

inline void Asura_AVector4::operator += (Asura_AVector4_InParam xOther)
{
	m_xVector4 = vec_add(m_xVector4,xOther.m_xVector4);
}

inline void Asura_AVector4::operator -= (Asura_AVector4_InParam xOther)
{
	m_xVector4 = vec_sub(m_xVector4,xOther.m_xVector4);
}

inline void Asura_AVector4::operator *= (const float fScale)
{
	vector float xScale = { fScale, fScale, fScale, fScale };
	m_xVector4 = vec_madd( m_xVector4, xScale, (vector float)(0.0f) );
}

inline void Asura_AVector4::operator /= (const float fScale)
{
	float fInvScale = 1.0f/fScale;
	vector float xInvScale = { fInvScale, fInvScale, fInvScale, fInvScale };
	m_xVector4 = vec_madd( m_xVector4, xInvScale, (vector float)(0.0f) );
}

inline Asura_AVector4 Asura_AVector4::operator - (void) const
{
	return Asura_AVector4( negatef4(m_xVector4) );
}

inline Asura_AVector4 Asura_AVector4::operator + (Asura_AVector4_InParam xOther) const
{
	return Asura_AVector4( vec_add( m_xVector4, xOther.m_xVector4 ) );
}

inline Asura_AVector4 Asura_AVector4::operator - (Asura_AVector4_InParam xOther) const
{
	return Asura_AVector4( vec_sub( m_xVector4, xOther.m_xVector4 ) );
}

inline Asura_AVector4 Asura_AVector4::operator ^ (Asura_AVector4_InParam xOther) const
{
	const float* pf = &m_xVector4[0];
	const float* pfOther = &xOther.m_xVector4[0];
	return(Asura_AVector4((pf[1] * pfOther[2]) - (pf[2] * pfOther[1]), (pf[2] * pfOther[0]) - (pf[0] * pfOther[2]), (pf[0] * pfOther[1]) - (pf[1] * pfOther[0]), 0.0f));
}

inline float Asura_AVector4::operator * (Asura_AVector4_InParam xOther) const
{
	return DotProduct4(xOther);
}

inline Asura_AVector4 Asura_AVector4::operator * (const float fScale) const
{
	vector float xScale = { fScale, fScale, fScale, fScale };
	return Asura_AVector4( vec_madd( m_xVector4, xScale, (vector float)(0.0f) ) );
}

inline Asura_AVector4 Asura_AVector4::operator / (const float fScale) const
{
	float fInv = 1.0f / fScale;
	vector float xInv = { fInv, fInv, fInv, fInv };
	return Asura_AVector4( vec_madd( m_xVector4, xInv, (vector float)(0.0f) ) );
}

inline bool Asura_AVector4::operator == (Asura_AVector4_InParam xOther) const
{
	return vec_all_eq( m_xVector4, xOther.m_xVector4 );
}

inline bool Asura_AVector4::ApproximatelyEquals(Asura_AVector4_InParam xOther, const float fEpsilon ) const
{
	Asura_AVector4 xDifference( *this - xOther );

	return( xDifference.MagnitudeSquared() < fEpsilon );
}

inline bool Asura_AVector4::operator != (Asura_AVector4_InParam xOther) const
{
	return !(vec_all_eq( m_xVector4, xOther.m_xVector4 ) );
}

inline float Asura_AVector4::GetElement(const u_int uIndex) const
{
	return m_xVector4[uIndex];
}

inline void Asura_AVector4::SetElement(const u_int uIndex, const float fValue)
{
	m_xVector4[uIndex] = fValue;
}

inline void Asura_AVector4::SetElements(const float fX, const float fY, const float fZ, const float fW)
{
	m_xVector4 = (vector float){ fX, fY, fZ, fW };
}

//////////////////////////////////////////////////////////////////////////
// Asura_AMatrix4x4
//////////////////////////////////////////////////////////////////////////
class Asura_AMatrix4x4
{
public:
	// Constructors/initializers
	Asura_AMatrix4x4(void) {}
	Asura_AMatrix4x4(const float fE11, const float fE12, const float fE13, const float fE14,
		const float fE21, const float fE22, const float fE23, const float fE24,
		const float fE31, const float fE32, const float fE33, const float fE34,
		const float fE41, const float fE42, const float fE43, const float fE44)
	{
		m_axColumns[0].SetElements(fE11, fE21, fE31, fE41);
		m_axColumns[1].SetElements(fE12, fE22, fE32, fE42);
		m_axColumns[2].SetElements(fE13, fE23, fE33, fE43);
		m_axColumns[3].SetElements(fE14, fE24, fE34, fE44);
	}
	Asura_AMatrix4x4(const Asura_AMatrix4x4& xOther)
	{
		m_axColumns[0] = xOther.m_axColumns[0];
		m_axColumns[1] = xOther.m_axColumns[1];
		m_axColumns[2] = xOther.m_axColumns[2];
		m_axColumns[3] = xOther.m_axColumns[3];
	}
	Asura_AMatrix4x4(const Asura_Matrix_3x3& xOrientation, const Asura_Vector_3& xTranslation)
	{
		m_axColumns[0].SetElements(xOrientation.e11,		xOrientation.e21,		xOrientation.e31,		0.0f);
		m_axColumns[1].SetElements(xOrientation.e12,		xOrientation.e22,		xOrientation.e32,		0.0f);
		m_axColumns[2].SetElements(xOrientation.e13,		xOrientation.e23,		xOrientation.e33,		0.0f);
		m_axColumns[3].SetElements(xTranslation.GetX(),	xTranslation.GetY(),	xTranslation.GetZ(),	1.0f);
	}
	Asura_AMatrix4x4(const Asura_Quat& xOrientation, const Asura_Vector_3& xTranslation)
	{
		const float	fS	= 2.0f;	//xRotation.Magnitude();

		const float	fXS	= xOrientation.x * fS;
		const float	fYS	= xOrientation.y * fS;
		const float	fZS	= xOrientation.z * fS;

		const float	fWX	= xOrientation.w * fXS;
		const float	fWY	= xOrientation.w * fYS;
		const float	fWZ	= xOrientation.w * fZS;

		const float	fXX	= xOrientation.x * fXS;
		const float	fXY	= xOrientation.x * fYS;
		const float	fXZ	= xOrientation.x * fZS;

		const float	fYY	= xOrientation.y * fYS;
		const float	fYZ	= xOrientation.y * fZS;
		const float	fZZ	= xOrientation.z * fZS;

		// Implementation from Graphics Gems IV...
		m_axColumns[0].SetElements(1.0f - (fYY + fZZ),	fXY + fWZ,			fXZ - fWY,			0.0f);
		m_axColumns[1].SetElements(fXY - fWZ,			1.0f - (fXX + fZZ),	fYZ + fWX,			0.0f);
		m_axColumns[2].SetElements(fXZ + fWY,			fYZ - fWX,			1.0f - (fXX + fYY),	0.0f);
		m_axColumns[3].SetElements(xTranslation.GetX(),	xTranslation.GetY(),	xTranslation.GetZ(),	1.0f);
	}

	// Constant functions...
	void Identity(void)
	{
		m_axColumns[0].SetElements(1.0f, 0.0f, 0.0f, 0.0f);
		m_axColumns[1].SetElements(0.0f, 1.0f, 0.0f, 0.0f);
		m_axColumns[2].SetElements(0.0f, 0.0f, 1.0f, 0.0f);
		m_axColumns[3].ZeroW1();
	}
	void Zero(void)
	{
		m_axColumns[0].ZeroW0();
		m_axColumns[1].ZeroW0();
		m_axColumns[2].ZeroW0();
		m_axColumns[3].ZeroW0();
	}

	// Orientation/translation functions...
	void ToMatrix(Asura_Matrix_3x3& xOrientation) const
	{
		xOrientation.e11	= m_axColumns[0].GetX();
		xOrientation.e21	= m_axColumns[0].GetY();
		xOrientation.e31	= m_axColumns[0].GetZ();

		xOrientation.e12	= m_axColumns[1].GetX();
		xOrientation.e22	= m_axColumns[1].GetY();
		xOrientation.e32	= m_axColumns[1].GetZ();

		xOrientation.e13	= m_axColumns[2].GetX();
		xOrientation.e23	= m_axColumns[2].GetY();
		xOrientation.e33	= m_axColumns[2].GetZ();
	}
	void FromMatrix(const Asura_Matrix_3x3& xOrientation)
	{
		m_axColumns[0].SetElements(xOrientation.e11, xOrientation.e21, xOrientation.e31, 0.0f);
		m_axColumns[1].SetElements(xOrientation.e12, xOrientation.e22, xOrientation.e32, 0.0f);
		m_axColumns[2].SetElements(xOrientation.e13, xOrientation.e23, xOrientation.e33, 0.0f);
		m_axColumns[3].ZeroW1();
	}

	Asura_Vector_3 GetTranslation() const
	{
		return Asura_Vector_3(m_axColumns[3].GetX(),m_axColumns[3].GetY(),m_axColumns[3].GetZ());
	}

	void SetTranslation(const Asura_Vector_3& xTranslation)
	{
		m_axColumns[3]	= Asura_AVector4(xTranslation, 1.0f);
	}

	Asura_AVector3 GetATranslation() const
	{
		return Asura_AVector3(m_axColumns[3].GetX(),m_axColumns[3].GetY(),m_axColumns[3].GetZ());
	}

	void SetTranslation(Asura_AVector3_InParam xTranslation)
	{
		m_axColumns[3].m_xVector4 = xTranslation.m_xVector3;
		m_axColumns[3].SetElement( 3, 1.0f );
	}

	Asura_AVector4 GetA4Translation() const
	{
		return m_axColumns[3];
	}
	void SetTranslation(Asura_AVector4_InParam xTranslation)
	{
		m_axColumns[3].m_xVector4 = xTranslation.m_xVector4;
	}

	// Manipulation functions...
	Asura_AVector4 Apply(Asura_AVector4_InParam xInput) const
	{
		vector float xTemp0, xTemp1, xResult;
		vector float xX, xY, xZ, xW;
		vec_float4 xZero = (vector float)(0.0f);
		xX = vec_splat( xInput.m_xVector4, 0 );
		xY = vec_splat( xInput.m_xVector4, 1 );
		xZ = vec_splat( xInput.m_xVector4, 2 );
		xW = vec_splat( xInput.m_xVector4, 3 );
		xTemp0 = vec_madd( m_axColumns[0].m_xVector4, xX, xZero );
		xTemp1 = vec_madd( m_axColumns[1].m_xVector4, xY, xZero );
		xTemp0 = vec_madd( m_axColumns[2].m_xVector4, xZ, xTemp0 );
	    xTemp1 = vec_madd( m_axColumns[3].m_xVector4, xW, xTemp1 );
	    xResult = vec_add( xTemp0, xTemp1 );
		Asura_AVector4 xOutput;
		xOutput.m_xVector4 = xResult;
		return xOutput;
	}
	void ApplyToArray(const Asura_AVector4* pxInput, const u_int uNumElements, Asura_AVector4* pxOutput) const
	{
		vector float xTemp0, xTemp1, xResult;
		vector float xX, xY, xZ, xW;
		vec_float4 xZero = (vector float)(0.0f);
		for( u_int u = 0; u < uNumElements; u++, pxInput++, pxOutput++ )
		{
			xX = vec_splat( pxInput->m_xVector4, 0 );
			xY = vec_splat( pxInput->m_xVector4, 1 );
			xZ = vec_splat( pxInput->m_xVector4, 2 );
			xW = vec_splat( pxInput->m_xVector4, 3 );
			xTemp0 = vec_madd( m_axColumns[0].m_xVector4, xX, xZero );
			xTemp1 = vec_madd( m_axColumns[1].m_xVector4, xY, xZero );
			xTemp0 = vec_madd( m_axColumns[2].m_xVector4, xZ, xTemp0 );
			xTemp1 = vec_madd( m_axColumns[3].m_xVector4, xW, xTemp1 );
			xResult = vec_add( xTemp0, xTemp1 );
			pxOutput->m_xVector4 = xResult;
		}
	}

	void Apply(Asura_AMatrix4x4_InParam xInput, Asura_AMatrix4x4_OutParam xOutput) const
	{
		// calculate the results into a temporary output matrix in case our output matrix is the
		// same matrix as ourself
		Asura_AMatrix4x4 xTempOutput;
		
		xTempOutput.m_axColumns[0] = this->Apply( xInput.m_axColumns[0] );
		xTempOutput.m_axColumns[1] = this->Apply( xInput.m_axColumns[1] );
		xTempOutput.m_axColumns[2] = this->Apply( xInput.m_axColumns[2] );
		xTempOutput.m_axColumns[3] = this->Apply( xInput.m_axColumns[3] );
		
		xOutput = Asura_AMatrix4x4( xTempOutput );
	}

	void ApplyToArray(const Asura_AMatrix4x4* pxInput, const u_int uNumElements, Asura_AMatrix4x4* pxOutput) const
	{
		for( u_int u = 0; u < uNumElements; u++, pxInput++, pxOutput++ )
		{
			pxOutput->m_axColumns[0] = this->Apply( pxInput->m_axColumns[0] );
			pxOutput->m_axColumns[1] = this->Apply( pxInput->m_axColumns[1] );
			pxOutput->m_axColumns[2] = this->Apply( pxInput->m_axColumns[2] );
			pxOutput->m_axColumns[3] = this->Apply( pxInput->m_axColumns[3] );
		}
	}

	float xCofactor(const u_int uRow, const u_int uCol) const
	{
		static const u_int	aauIndices[4][3]	= { {1, 2, 3}, {0, 2, 3}, {0, 1, 3}, {0, 1, 2} };
		const u_int*		puRowsToUse			= aauIndices[uRow];
		const u_int*		puColsToUse			= aauIndices[uCol];

		// Calculate the determinant of the sub-matrix produced by ignoring row uRow and column uCol
		const float	fA	= GetElement(puRowsToUse[0], puColsToUse[0]);
		const float	fB	= GetElement(puRowsToUse[0], puColsToUse[1]);
		const float	fC	= GetElement(puRowsToUse[0], puColsToUse[2]);

		const float	fD	= GetElement(puRowsToUse[1], puColsToUse[0]);
		const float	fE	= GetElement(puRowsToUse[1], puColsToUse[1]);
		const float	fF	= GetElement(puRowsToUse[1], puColsToUse[2]);

		const float	fG	= GetElement(puRowsToUse[2], puColsToUse[0]);
		const float	fH	= GetElement(puRowsToUse[2], puColsToUse[1]);
		const float	fI	= GetElement(puRowsToUse[2], puColsToUse[2]);

		return( ((fA * fE * fI) + (fD * fH * fC) + (fG * fB * fF)) -
			((fA * fH * fF) + (fD * fB * fI) + (fG * fE * fC)) );
	}
	float Determinant(void) const
	{
		/* function implementation based on code from STIDC SDK:           */
		/* --------------------------------------------------------------  */
		/* PLEASE DO NOT MODIFY THIS SECTION                               */
		/* This prolog section is automatically generated.                 */
		/*                                                                 */
		/* (C)Copyright                                                    */
		/* Sony Computer Entertainment, Inc.,                              */
		/* Toshiba Corporation,                                            */
		/* International Business Machines Corporation,                    */
		/* 2001,2002.                                                      */
		/* S/T/I Confidential Information                                  */
		/* --------------------------------------------------------------  */
		vector float xIn0, xIn1, xIn2, xIn3;
		vector float xTemp0, xTemp1, xTemp2, xTemp3;
		vector float xCofactor0;
		vector float xT0, xT1, xT2, xT3;
		vector float xT12, xT23;
		vector float xT1r, xT2r;
		vector float xT12r, xT23r;
		vector float xT1r3, xT1r3r;
		vector float xZero = (vector float)(0.0);
		vector float xResult;
		xIn0 = m_axColumns[0].m_xVector4;
		xIn1 = m_axColumns[1].m_xVector4;
		xIn2 = m_axColumns[2].m_xVector4;
		xIn3 = m_axColumns[3].m_xVector4;
		/* Perform transform of the input matrix of the form:
		 *    A B C D
		 *    E F G H
		 *    I J K L
		 *    M N O P
		 *
		 * The pseudo transpose of the input matrix is trans:
		 *    A E I M
		 *    J N B F
		 *    C G K O
		 *    L P D H
		 */
		xTemp0 = vec_perm(xIn0, xIn1, PERM_XAZC );//(PERM_XAZC);	/* A E C G */
		xTemp1 = vec_perm(xIn2, xIn3, PERM_XAZC );	/* I M K O */
		xTemp2 = vec_perm(xIn0, xIn1, PERM_YBWD);	/* B F D H */
		xTemp3 = vec_perm(xIn2, xIn3, PERM_YBWD);	/* J N L P */
		xT0 = vec_perm(xTemp0, xTemp1, PERM_XYAB);	/* A E I M */
		xT1 = vec_perm(xTemp3, xTemp2, PERM_XYAB);	/* J N B F */
		xT2 = vec_perm(xTemp0, xTemp1, PERM_ZWCD);	/* C G K O */
		xT3 = vec_perm(xTemp3, xTemp2, PERM_ZWCD);	/* L P D H */
		/* Generate a xCofactoractor matrix. The computed xCofactoractors reside in
		 * xCofactor0, xCofactor1, xCofactor2, xCofactor3.
		 */
		xT23 = vec_madd(xT2, xT3, xZero);		/* CL GP KD OH */
		xT23 = vec_perm(xT23, xT23, PERM_YXWZ);	/* GP CL OH KD */
		xCofactor0 = vec_nmsub(xT1, xT23, xZero);		/* -(JGP NCL FOH BKD) */
		xT23r = vec_sld(xT23, xT23, 8);			/* OH KD GP CL */
		xCofactor0 = vec_madd(xT1, xT23r, xCofactor0);		/* JOH NKD BGP FCL + xCofactor0 */
		xT12 = vec_madd(xT1, xT2, xZero);		/* JC NG BK FO */
		xT12 = vec_perm(xT12, xT12, PERM_YXWZ);	/* NG JC FO BK */
		xCofactor0 = vec_madd(xT3, xT12, xCofactor0);		/* LNG PJC DFO HBK + xCofactor0 */
		xT12r = vec_sld(xT12, xT12, 8);			/* FO BK NG JC */
		xCofactor0 = vec_nmsub(xT3, xT12r, xCofactor0);		/* xCofactor0 - LFO PBK DNG HJC */
		xT1r = vec_sld(xT1, xT1, 8);			/* B F J N */
		xT2r = vec_sld(xT2, xT2, 8);			/* K O C G */
		xT1r3 = vec_madd(xT1r, xT3, xZero);		/* BL FP JD NH */
		xT1r3 = vec_perm(xT1r3, xT1r3, PERM_YXWZ);	/* FP BL NH JD */
		xCofactor0 = vec_madd(xT2r, xT1r3, xCofactor0);		/* KFP OBL CNH GJD + xCofactor0 */
		xT1r3r = vec_sld(xT1r3, xT1r3, 8);		/* NH JD FP BL */
		xCofactor0 = vec_nmsub(xT2r, xT1r3r, xCofactor0);		/* xCofactor0 - KNH OJD CFP GBL */
		xResult = vec_madd( xT0, xCofactor0, xZero );
		xResult = vec_madd( vec_sld(xT0,xT0,4), vec_sld(xCofactor0,xCofactor0,4), xResult );
		xResult = vec_madd( vec_sld(xT0,xT0,8), vec_sld(xCofactor0,xCofactor0,8), xResult);
		return *(reinterpret_cast<float *>(&xResult) );
	}

	void Invert(void)
	{
		// full inverse calculation
		/* function implementation based on code from STIDC SDK:           */
		/* --------------------------------------------------------------  */
		/* PLEASE DO NOT MODIFY THIS SECTION                               */
		/* This prolog section is automatically generated.                 */
		/*                                                                 */
		/* (C)Copyright                                                    */
		/* Sony Computer Entertainment, Inc.,                              */
		/* Toshiba Corporation,                                            */
		/* International Business Machines Corporation,                    */
		/* 2001,2002.                                                      */
		/* S/T/I Confidential Information                                  */
		/* --------------------------------------------------------------  */
		vector float xIn0, xIn1, xIn2, xIn3;
		vector float xTemp0, xTemp1, xTemp2, xTemp3;
		vector float xCofactor0, xCofactor1, xCofactor2, xCofactor3;
		vector float xT0, xT1, xT2, xT3;
		vector float xT01, xT02, xT03, xT12, xT23;
		vector float xT1r, xT2r;
		vector float xT01r, xT02r, xT03r, xT12r, xT23r;
		vector float xT1r3, xT1r3r;
		vector float xDet, xDet0, xDet1, xDet2, xDet3, xInvDet;
		vector float xZero = (vector float)(0.0);
		xIn0 = m_axColumns[0].m_xVector4;
		xIn1 = m_axColumns[1].m_xVector4;
		xIn2 = m_axColumns[2].m_xVector4;
		xIn3 = m_axColumns[3].m_xVector4;
		/* Perform transform of the input matrix of the form:
		 *    A B C D
		 *    E F G H
		 *    I J K L
		 *    M N O P
		 *
		 * The pseudo transpose of the input matrix is trans:
		 *    A E I M
		 *    J N B F
		 *    C G K O
		 *    L P D H
		 */
		xTemp0 = vec_perm(xIn0, xIn1, PERM_XAZC);	/* A E C G */
		xTemp1 = vec_perm(xIn2, xIn3, PERM_XAZC);	/* I M K O */
		xTemp2 = vec_perm(xIn0, xIn1, PERM_YBWD);	/* B F D H */
		xTemp3 = vec_perm(xIn2, xIn3, PERM_YBWD);	/* J N L P */
		xT0 = vec_perm(xTemp0, xTemp1, PERM_XYAB);	/* A E I M */
		xT1 = vec_perm(xTemp3, xTemp2, PERM_XYAB);	/* J N B F */
		xT2 = vec_perm(xTemp0, xTemp1, PERM_ZWCD);	/* C G K O */
		xT3 = vec_perm(xTemp3, xTemp2, PERM_ZWCD);	/* L P D H */
		/* Generate a xCofactoractor matrix. The computed xCofactoractors reside in
		 * xCofactor0, xCofactor1, xCofactor2, xCofactor3.
		 */
		xT23 = vec_madd(xT2, xT3, xZero);		/* CL GP KD OH */
		xT23 = vec_perm(xT23, xT23, PERM_YXWZ);	/* GP CL OH KD */
		xCofactor0 = vec_nmsub(xT1, xT23, xZero);		/* -(JGP NCL FOH BKD) */
		xCofactor1 = vec_nmsub(xT0, xT23, xZero);		/* -(AGP ECL IOH MKD) */
		xT23r = vec_sld(xT23, xT23, 8);			/* OH KD GP CL */
		xCofactor0 = vec_madd(xT1, xT23r, xCofactor0);		/* JOH NKD BGP FCL + xCofactor0 */
		xCofactor1 = vec_madd(xT0, xT23r, xCofactor1);		/* AOH EKD IGP MCL + xCofactor1 */
		xCofactor1 = vec_sld(xCofactor1, xCofactor1, 8);		/* IGP MCL AOH EKD - IOH MKD AGP ECL */
		xT12 = vec_madd(xT1, xT2, xZero);		/* JC NG BK FO */
		xT12 = vec_perm(xT12, xT12, PERM_YXWZ);	/* NG JC FO BK */
		xCofactor0 = vec_madd(xT3, xT12, xCofactor0);		/* LNG PJC DFO HBK + xCofactor0 */
		xCofactor3 = vec_madd(xT0, xT12, xZero);		/* ANG EJC IFO MBK */
		xT12r = vec_sld(xT12, xT12, 8);			/* FO BK NG JC */
		xCofactor0 = vec_nmsub(xT3, xT12r, xCofactor0);		/* xCofactor0 - LFO PBK DNG HJC */
		xCofactor3 = vec_nmsub(xT0, xT12r, xCofactor3);		/* xCofactor3 - AFO EBK ING MJC */
		xCofactor3 = vec_sld(xCofactor3, xCofactor3, 8);		/* ING MJC AFO EBK - IFO MBK ANG EJC */
		xT1r = vec_sld(xT1, xT1, 8);			/* B F J N */
		xT2r = vec_sld(xT2, xT2, 8);			/* K O C G */
		xT1r3 = vec_madd(xT1r, xT3, xZero);		/* BL FP JD NH */
		xT1r3 = vec_perm(xT1r3, xT1r3, PERM_YXWZ);	/* FP BL NH JD */
		xCofactor0 = vec_madd(xT2r, xT1r3, xCofactor0);		/* KFP OBL CNH GJD + xCofactor0 */
		xCofactor2 = vec_madd(xT0, xT1r3, xZero);		/* AFP EBL INH MJD */
		xT1r3r = vec_sld(xT1r3, xT1r3, 8);		/* NH JD FP BL */
		xCofactor0 = vec_nmsub(xT2r, xT1r3r, xCofactor0);		/* xCofactor0 - KNH OJD CFP GBL */
		xCofactor2 = vec_nmsub(xT0, xT1r3r, xCofactor2);		/* xCofactor2 - ANH EJD IFP MBL */
		xCofactor2 = vec_sld(xCofactor2, xCofactor2, 8);		/* IFP MBL ANH EJD - INH MJD AFP EBL */
		xT01 = vec_madd(xT0, xT1, xZero);		/* AJ EN IB MF */
		xT01 = vec_perm(xT01, xT01, PERM_YXWZ);	/* EN AJ MF IB */
		xCofactor2 = vec_nmsub(xT3, xT01, xCofactor2);		/* xCofactor2 - LEN PAJ DMF HIB */
		xCofactor3 = vec_madd(xT2r, xT01, xCofactor3);		/* KEN OAJ CMF GIB + xCofactor3 */ 
		xT01r = vec_sld(xT01, xT01, 8);			/* MF IB EN AJ */
		xCofactor2 = vec_madd(xT3, xT01r, xCofactor2);		/* LMF PIB DEN HAJ + xCofactor2 */
		xCofactor3 = vec_nmsub(xT2r, xT01r, xCofactor3);		/* xCofactor3 - KMF OIB CEN GAJ */
		xT03 = vec_madd(xT0, xT3, xZero);		/* AL EP ID MH */
		xT03 = vec_perm(xT03, xT03, PERM_YXWZ);	/* EP AL MH ID */
		xCofactor1 = vec_nmsub(xT2r, xT03, xCofactor1);		/* xCofactor1 - KEP OAL CMH GID */
		xCofactor2 = vec_madd(xT1, xT03, xCofactor2);		/* JEP NAL BMH FID + xCofactor2 */
		xT03r = vec_sld(xT03, xT03, 8);			/* MH ID EP AL */
		xCofactor1 = vec_madd(xT2r, xT03r, xCofactor1);		/* KMH OID CEP GAL + xCofactor1 */
		xCofactor2 = vec_nmsub(xT1, xT03r, xCofactor2);		/* xCofactor2 - JMH NID BEP FAL */ 
		xT02 = vec_madd(xT0, xT2r, xZero);		/* AK EO IC MG */
		xT02 = vec_perm(xT02, xT02, PERM_YXWZ);	/* E0 AK MG IC */
		xCofactor1 = vec_madd(xT3, xT02, xCofactor1);		/* LEO PAK DMG HIC + xCofactor1 */
		xCofactor3 = vec_nmsub(xT1, xT02, xCofactor3);		/* xCofactor3 - JEO NAK BMG FIC */
		xT02r = vec_sld(xT02, xT02, 8);			/* MG IC EO AK */
		xCofactor1 = vec_nmsub(xT3, xT02r, xCofactor1);		/* xCofactor1 - LMG PIC DEO HAK */
		xCofactor3 = vec_madd(xT1, xT02r, xCofactor3);		/* JMG NIC BEO FAK + xCofactor3 */
		/* Compute the determinant of the matrix 
		 *
		 * det = sum_across(t0 * xCofactor0);
		 *
		 * We perform a sum across the entire vector so that 
		 * we don't have to splat the result when multiplying the
		 * xCofactoractors by the inverse of the determinant.
		 */
		xDet  = vec_madd(xT0, xCofactor0, xZero);
		xDet0 = vec_splat(xDet, 0);
		xDet1 = vec_splat(xDet, 1);
		xDet2 = vec_splat(xDet, 2);
		xDet3 = vec_splat(xDet, 3);
		xDet  = vec_add(xDet0, xDet1);
		xDet2 = vec_add(xDet2, xDet3);
		xDet  = vec_add(xDet, xDet2);
		/* Compute the reciprocal of the determinant.
		 */
		xInvDet = recipf4(xDet);
		/* Multiply the xCofactoractors by the reciprocal of the determinant.
		 */ 
		m_axColumns[0] = Asura_AVector4( vec_madd(xCofactor0, xInvDet, xZero) );
        m_axColumns[1] = Asura_AVector4( vec_madd(xCofactor1, xInvDet, xZero) );
        m_axColumns[2] = Asura_AVector4( vec_madd(xCofactor2, xInvDet, xZero) );
        m_axColumns[3] = Asura_AVector4( vec_madd(xCofactor3, xInvDet, xZero) );
	}
	void InvertOrthogonal(void)
	{
		// in the case that the matrix has an orthogonal top left hand corner
		// then the transpose is the inverse. We can think of this form of
		// matrix as applying a translation matrix to a rotation matrix. The 
		// inverse is then simply the inverse of the rotation matrix (the transpose)
		// applied to the inverse of the translation (replace the translation x,y,z
		// values with their negatives).

		// columns of final inverse matrix
		vector float xInv0, xInv1, xInv2, xInv3;
		
		vector float xTemp0, xTemp1;
		// vectors containing negated translation components
		vector float xX, xY, xZ;

		xTemp0 = vec_mergeh( m_axColumns[0].m_xVector4, m_axColumns[2].m_xVector4 );
		xTemp1 = vec_mergel( m_axColumns[0].m_xVector4, m_axColumns[2].m_xVector4 );

		// create inverse transform vector
		xInv3 = negatef4( m_axColumns[3].m_xVector4 );
		// create transpose of rotation matrix
		xInv0 = vec_mergeh( xTemp0, m_axColumns[1].m_xVector4 );
		xInv1 = vec_perm( xTemp0, m_axColumns[1].m_xVector4, (vector unsigned char)(0x08,0x09,0x0a,0x0b,0x14,0x15,0x16,0x17,0x0c,0x0d,0x0e,0x0f,0x00,0x01,0x02,0x03) );
		xInv2 = vec_perm( xTemp1, m_axColumns[1].m_xVector4, (vector unsigned char)(0x00,0x01,0x02,0x03,0x18,0x19,0x1a,0x1b,0x04,0x05,0x06,0x07,0x00,0x01,0x02,0x03) );

		// multiply together to get inverse
		xX = vec_splat( xInv3, 0 );
		xY = vec_splat( xInv3, 1 );
		xZ = vec_splat( xInv3, 2 );
		xInv3 = vec_madd( xInv0, xX, (vector float)(0.0f) );
		xInv3 = vec_madd( xInv1, xY, xInv3 );
		xInv3 = vec_madd( xInv2, xZ, xInv3 );

		// we also need to make sure that the bottom row consists of 0001
		vec_bint4 xSelVec = (vec_bint4)(0x0,0x0,0x0,0xFFFFFFFF);
		vector float xZeroVec = (vector float)(0.0f);
		vector float xOneVec = (vector float)(1.0f);
		xInv0 = vec_sel( xInv0, xZeroVec, xSelVec );
		xInv1 = vec_sel( xInv1, xZeroVec, xSelVec );
		xInv2 = vec_sel( xInv2, xZeroVec, xSelVec );
		xInv3 = vec_sel( xInv3, xOneVec, xSelVec );

		// construct the output matrix
		Asura_AVector4 xCol1(xInv0);
		Asura_AVector4 xCol2(xInv1);
		Asura_AVector4 xCol3(xInv2);
		Asura_AVector4 xCol4(xInv3);

		SetColumns(xCol1,xCol2,xCol3,xCol4);
	}
	void Transpose(void)
	{
		vector float xTemp0, xTemp1, xTemp2, xTemp3, xRes0, xRes1, xRes2, xRes3;
		xTemp0 = vec_mergeh( m_axColumns[0].m_xVector4, m_axColumns[2].m_xVector4 );
		xTemp1 = vec_mergeh( m_axColumns[1].m_xVector4, m_axColumns[3].m_xVector4 );
		xTemp2 = vec_mergel( m_axColumns[0].m_xVector4, m_axColumns[2].m_xVector4 );
		xTemp3 = vec_mergel( m_axColumns[1].m_xVector4, m_axColumns[3].m_xVector4 );
		xRes0 = vec_mergeh( xTemp0, xTemp1 );
		xRes1 = vec_mergel( xTemp0, xTemp1 );
		xRes2 = vec_mergeh( xTemp2, xTemp3 );
		xRes3 = vec_mergel( xTemp2, xTemp3 );
		m_axColumns[0] = Asura_AVector4(xRes0);
		m_axColumns[1] = Asura_AVector4(xRes1);
		m_axColumns[2] = Asura_AVector4(xRes2);
		m_axColumns[3] = Asura_AVector4(xRes3);	
	}

	void Invert(Asura_AMatrix4x4_OutParam xOther) const
	{
		// full inverse calculation
		/* function implementation based on code from STIDC SDK:           */
		/* --------------------------------------------------------------  */
		/* PLEASE DO NOT MODIFY THIS SECTION                               */
		/* This prolog section is automatically generated.                 */
		/*                                                                 */
		/* (C)Copyright                                                    */
		/* Sony Computer Entertainment, Inc.,                              */
		/* Toshiba Corporation,                                            */
		/* International Business Machines Corporation,                    */
		/* 2001,2002.                                                      */
		/* S/T/I Confidential Information                                  */
		/* --------------------------------------------------------------  */
		vector float xIn0, xIn1, xIn2, xIn3;
		vector float xTemp0, xTemp1, xTemp2, xTemp3;
		vector float xCofactor0, xCofactor1, xCofactor2, xCofactor3;
		vector float xT0, xT1, xT2, xT3;
		vector float xT01, xT02, xT03, xT12, xT23;
		vector float xT1r, xT2r;
		vector float xT01r, xT02r, xT03r, xT12r, xT23r;
		vector float xT1r3, xT1r3r;
		vector float xDet, xDet0, xDet1, xDet2, xDet3, xInvDet;
		vector float xZero = (vector float)(0.0f);
		xIn0 = m_axColumns[0].m_xVector4;
		xIn1 = m_axColumns[1].m_xVector4;
		xIn2 = m_axColumns[2].m_xVector4;
		xIn3 = m_axColumns[3].m_xVector4;
		/* Perform transform of the input matrix of the form:
		 *    A B C D
		 *    E F G H
		 *    I J K L
		 *    M N O P
		 *
		 * The pseudo transpose of the input matrix is trans:
		 *    A E I M
		 *    J N B F
		 *    C G K O
		 *    L P D H
		 */
		xTemp0 = vec_perm(xIn0, xIn1, PERM_XAZC);	/* A E C G */
		xTemp1 = vec_perm(xIn2, xIn3, PERM_XAZC);	/* I M K O */
		xTemp2 = vec_perm(xIn0, xIn1, PERM_YBWD);	/* B F D H */
		xTemp3 = vec_perm(xIn2, xIn3, PERM_YBWD);	/* J N L P */
		xT0 = vec_perm(xTemp0, xTemp1, PERM_XYAB);	/* A E I M */
		xT1 = vec_perm(xTemp3, xTemp2, PERM_XYAB);	/* J N B F */
		xT2 = vec_perm(xTemp0, xTemp1, PERM_ZWCD);	/* C G K O */
		xT3 = vec_perm(xTemp3, xTemp2, PERM_ZWCD);	/* L P D H */
		/* Generate a xCofactoractor matrix. The computed xCofactoractors reside in
		 * xCofactor0, xCofactor1, xCofactor2, xCofactor3.
		 */
		xT23 = vec_madd(xT2, xT3, xZero);		/* CL GP KD OH */
		xT23 = vec_perm(xT23, xT23, PERM_YXWZ);	/* GP CL OH KD */
		xCofactor0 = vec_nmsub(xT1, xT23, xZero);		/* -(JGP NCL FOH BKD) */
		xCofactor1 = vec_nmsub(xT0, xT23, xZero);		/* -(AGP ECL IOH MKD) */
		xT23r = vec_sld(xT23, xT23, 8);			/* OH KD GP CL */
		xCofactor0 = vec_madd(xT1, xT23r, xCofactor0);		/* JOH NKD BGP FCL + xCofactor0 */
		xCofactor1 = vec_madd(xT0, xT23r, xCofactor1);		/* AOH EKD IGP MCL + xCofactor1 */
		xCofactor1 = vec_sld(xCofactor1, xCofactor1, 8);		/* IGP MCL AOH EKD - IOH MKD AGP ECL */
		xT12 = vec_madd(xT1, xT2, xZero);		/* JC NG BK FO */
		xT12 = vec_perm(xT12, xT12, PERM_YXWZ);	/* NG JC FO BK */
		xCofactor0 = vec_madd(xT3, xT12, xCofactor0);		/* LNG PJC DFO HBK + xCofactor0 */
		xCofactor3 = vec_madd(xT0, xT12, xZero);		/* ANG EJC IFO MBK */
		xT12r = vec_sld(xT12, xT12, 8);			/* FO BK NG JC */
		xCofactor0 = vec_nmsub(xT3, xT12r, xCofactor0);		/* xCofactor0 - LFO PBK DNG HJC */
		xCofactor3 = vec_nmsub(xT0, xT12r, xCofactor3);		/* xCofactor3 - AFO EBK ING MJC */
		xCofactor3 = vec_sld(xCofactor3, xCofactor3, 8);		/* ING MJC AFO EBK - IFO MBK ANG EJC */
		xT1r = vec_sld(xT1, xT1, 8);			/* B F J N */
		xT2r = vec_sld(xT2, xT2, 8);			/* K O C G */
		xT1r3 = vec_madd(xT1r, xT3, xZero);		/* BL FP JD NH */
		xT1r3 = vec_perm(xT1r3, xT1r3, PERM_YXWZ);	/* FP BL NH JD */
		xCofactor0 = vec_madd(xT2r, xT1r3, xCofactor0);		/* KFP OBL CNH GJD + xCofactor0 */
		xCofactor2 = vec_madd(xT0, xT1r3, xZero);		/* AFP EBL INH MJD */
		xT1r3r = vec_sld(xT1r3, xT1r3, 8);		/* NH JD FP BL */
		xCofactor0 = vec_nmsub(xT2r, xT1r3r, xCofactor0);		/* xCofactor0 - KNH OJD CFP GBL */
		xCofactor2 = vec_nmsub(xT0, xT1r3r, xCofactor2);		/* xCofactor2 - ANH EJD IFP MBL */
		xCofactor2 = vec_sld(xCofactor2, xCofactor2, 8);		/* IFP MBL ANH EJD - INH MJD AFP EBL */
		xT01 = vec_madd(xT0, xT1, xZero);		/* AJ EN IB MF */
		xT01 = vec_perm(xT01, xT01, PERM_YXWZ);	/* EN AJ MF IB */
		xCofactor2 = vec_nmsub(xT3, xT01, xCofactor2);		/* xCofactor2 - LEN PAJ DMF HIB */
		xCofactor3 = vec_madd(xT2r, xT01, xCofactor3);		/* KEN OAJ CMF GIB + xCofactor3 */ 
		xT01r = vec_sld(xT01, xT01, 8);			/* MF IB EN AJ */
		xCofactor2 = vec_madd(xT3, xT01r, xCofactor2);		/* LMF PIB DEN HAJ + xCofactor2 */
		xCofactor3 = vec_nmsub(xT2r, xT01r, xCofactor3);		/* xCofactor3 - KMF OIB CEN GAJ */
		xT03 = vec_madd(xT0, xT3, xZero);		/* AL EP ID MH */
		xT03 = vec_perm(xT03, xT03, PERM_YXWZ);	/* EP AL MH ID */
		xCofactor1 = vec_nmsub(xT2r, xT03, xCofactor1);		/* xCofactor1 - KEP OAL CMH GID */
		xCofactor2 = vec_madd(xT1, xT03, xCofactor2);		/* JEP NAL BMH FID + xCofactor2 */
		xT03r = vec_sld(xT03, xT03, 8);			/* MH ID EP AL */
		xCofactor1 = vec_madd(xT2r, xT03r, xCofactor1);		/* KMH OID CEP GAL + xCofactor1 */
		xCofactor2 = vec_nmsub(xT1, xT03r, xCofactor2);		/* xCofactor2 - JMH NID BEP FAL */ 
		xT02 = vec_madd(xT0, xT2r, xZero);		/* AK EO IC MG */
		xT02 = vec_perm(xT02, xT02, PERM_YXWZ);	/* E0 AK MG IC */
		xCofactor1 = vec_madd(xT3, xT02, xCofactor1);		/* LEO PAK DMG HIC + xCofactor1 */
		xCofactor3 = vec_nmsub(xT1, xT02, xCofactor3);		/* xCofactor3 - JEO NAK BMG FIC */
		xT02r = vec_sld(xT02, xT02, 8);			/* MG IC EO AK */
		xCofactor1 = vec_nmsub(xT3, xT02r, xCofactor1);		/* xCofactor1 - LMG PIC DEO HAK */
		xCofactor3 = vec_madd(xT1, xT02r, xCofactor3);		/* JMG NIC BEO FAK + xCofactor3 */
		/* Compute the determinant of the matrix 
		 *
		 * det = sum_across(t0 * xCofactor0);
		 *
		 * We perform a sum across the entire vector so that 
		 * we don't have to splat the result when multiplying the
		 * xCofactoractors by the inverse of the determinant.
		 */
		xDet  = vec_madd(xT0, xCofactor0, xZero);
		xDet0 = vec_splat(xDet, 0);
		xDet1 = vec_splat(xDet, 1);
		xDet2 = vec_splat(xDet, 2);
		xDet3 = vec_splat(xDet, 3);
		xDet  = vec_add(xDet0, xDet1);
		xDet2 = vec_add(xDet2, xDet3);
		xDet  = vec_add(xDet, xDet2);
		/* Compute the reciprocal of the determinant.
		 */
		xInvDet = recipf4(xDet);
		/* Multiply the xCofactoractors by the reciprocal of the determinant.
		 */ 
		xOther.m_axColumns[0] = Asura_AVector4( vec_madd(xCofactor0, xInvDet, xZero) );
        xOther.m_axColumns[1] = Asura_AVector4( vec_madd(xCofactor1, xInvDet, xZero) );
        xOther.m_axColumns[2] = Asura_AVector4( vec_madd(xCofactor2, xInvDet, xZero) );
        xOther.m_axColumns[3] = Asura_AVector4( vec_madd(xCofactor3, xInvDet, xZero) );
	}
	void InvertOrthogonal(Asura_AMatrix4x4_OutParam xOther) const
	{
		// in the case that the matrix has an orthogonal top left hand corner
		// then the transpose of the 3x3 matrix is its inverse. We can think of this form of
		// matrix as applying a translation matrix to a rotation matrix. The 
		// inverse is then simply the inverse of the rotation matrix (the transpose)
		// applied to the inverse of the translation (replace the translation x,y,z
		// values with their negatives).

		// columns of final inverse matrix
		vector float xInv0, xInv1, xInv2, xInv3;
		
		vector float xTemp0, xTemp1;
		// vectors containing negated translation components
		vector float xX, xY, xZ;

		xTemp0 = vec_mergeh( m_axColumns[0].m_xVector4, m_axColumns[2].m_xVector4 );
		xTemp1 = vec_mergel( m_axColumns[0].m_xVector4, m_axColumns[2].m_xVector4 );

		// create inverse transform vector
		xInv3 = negatef4( m_axColumns[3].m_xVector4 );
		// create transpose of rotation matrix
		xInv0 = vec_mergeh( xTemp0, m_axColumns[1].m_xVector4 );
		xInv1 = vec_perm( xTemp0, m_axColumns[1].m_xVector4, (vector unsigned char)(0x08,0x09,0x0a,0x0b,0x14,0x15,0x16,0x17,0x0c,0x0d,0x0e,0x0f,0x00,0x01,0x02,0x03) );
		xInv2 = vec_perm( xTemp1, m_axColumns[1].m_xVector4, (vector unsigned char)(0x00,0x01,0x02,0x03,0x18,0x19,0x1a,0x1b,0x04,0x05,0x06,0x07,0x00,0x01,0x02,0x03) );

		// multiply together to get inverse
		xX = vec_splat( xInv3, 0 );
		xY = vec_splat( xInv3, 1 );
		xZ = vec_splat( xInv3, 2 );
		xInv3 = vec_madd( xInv0, xX, (vector float)(0.0f) );
		xInv3 = vec_madd( xInv1, xY, xInv3 );
		xInv3 = vec_madd( xInv2, xZ, xInv3 );

		// we also need to make sure that the bottom row consists of 0001
		vec_bint4 xSelVec = (vec_bint4)(0x0,0x0,0x0,0xFFFFFFFF);
		vector float xZeroVec = (vector float)(0.0f);
		vector float xOneVec = (vector float)(1.0f);
		xInv0 = vec_sel( xInv0, xZeroVec, xSelVec );
		xInv1 = vec_sel( xInv1, xZeroVec, xSelVec );
		xInv2 = vec_sel( xInv2, xZeroVec, xSelVec );
		xInv3 = vec_sel( xInv3, xOneVec, xSelVec );

		// construct the output matrix
		Asura_AVector4 xCol1(xInv0);
		Asura_AVector4 xCol2(xInv1);
		Asura_AVector4 xCol3(xInv2);
		Asura_AVector4 xCol4(xInv3);

		xOther.SetColumns(xCol1,xCol2,xCol3,xCol4);
	}
	void Transpose(Asura_AMatrix4x4_OutParam xOther) const
	{
		vector float xTemp0, xTemp1, xTemp2, xTemp3, xResult0, xResult1, xResult2, xResult3;
		xTemp0 = vec_mergeh( m_axColumns[0].m_xVector4, m_axColumns[2].m_xVector4 );
		xTemp1 = vec_mergeh( m_axColumns[1].m_xVector4, m_axColumns[3].m_xVector4 );
		xTemp2 = vec_mergel( m_axColumns[0].m_xVector4, m_axColumns[2].m_xVector4 );
		xTemp3 = vec_mergel( m_axColumns[1].m_xVector4, m_axColumns[3].m_xVector4 );
		xResult0 = vec_mergeh( xTemp0, xTemp1 );
		xResult1 = vec_mergel( xTemp0, xTemp1 );
		xResult2 = vec_mergeh( xTemp2, xTemp3 );
		xResult3 = vec_mergel( xTemp2, xTemp3 );
		xOther.m_axColumns[0].m_xVector4 = xResult0;
		xOther.m_axColumns[1].m_xVector4 = xResult1;
		xOther.m_axColumns[2].m_xVector4 = xResult2;
		xOther.m_axColumns[3].m_xVector4 = xResult3;
	}

	// Unary operators
	Asura_AMatrix4x4 operator-(void) const
	{
		Asura_AMatrix4x4 xTemp;
		Asura_AVector4 xCol1 = -m_axColumns[0];
		Asura_AVector4 xCol2 = -m_axColumns[1];
		Asura_AVector4 xCol3 = -m_axColumns[2];
		Asura_AVector4 xCol4 = -m_axColumns[3];
		xTemp.SetColumns( xCol1, xCol2, xCol3, xCol4 );
		return xTemp;
	}

	// Binary operators
	Asura_AMatrix4x4 operator*(const float fScale) const
	{
		Asura_AMatrix4x4 xTemp;
		Asura_AVector4 xCol1 = m_axColumns[0] * fScale;
		Asura_AVector4 xCol2 = m_axColumns[1] * fScale;
		Asura_AVector4 xCol3 = m_axColumns[2] * fScale;
		Asura_AVector4 xCol4 = m_axColumns[3] * fScale;
		xTemp.SetColumns( xCol1, xCol2, xCol3, xCol4 );
		return xTemp;
	}
	Asura_AMatrix4x4 operator/(const float fScale) const
	{
		Asura_AMatrix4x4 xTemp;
		float fInvScale = 1.0f / fScale;
		Asura_AVector4 col1 = m_axColumns[0] * fInvScale;
		Asura_AVector4 col2 = m_axColumns[1] * fInvScale;
		Asura_AVector4 col3 = m_axColumns[2] * fInvScale;
		Asura_AVector4 col4 = m_axColumns[3] * fInvScale;
		xTemp.SetColumns( col1, col2, col3, col4 );
		return xTemp;
	}
	Asura_AMatrix4x4 operator+(Asura_AMatrix4x4_InParam xOther) const
	{
		Asura_AMatrix4x4 xTemp;
		Asura_AVector4 xCol1 = m_axColumns[0] + xOther.m_axColumns[0];
		Asura_AVector4 xCol2 = m_axColumns[1] + xOther.m_axColumns[1];
		Asura_AVector4 xCol3 = m_axColumns[2] + xOther.m_axColumns[2];
		Asura_AVector4 xCol4 = m_axColumns[3] + xOther.m_axColumns[3];
		xTemp.SetColumns( xCol1, xCol2, xCol3, xCol4 );
		return xTemp;
	}

	// Assignment operators
	void operator=(Asura_AMatrix4x4_InParam xOther)
	{
		m_axColumns[0] = xOther.m_axColumns[0];
		m_axColumns[1] = xOther.m_axColumns[1];
		m_axColumns[2] = xOther.m_axColumns[2];
		m_axColumns[3] = xOther.m_axColumns[3];
	}
	void operator*=(const float fScale)
	{
		m_axColumns[0] *= fScale;
		m_axColumns[1] *= fScale;
		m_axColumns[2] *= fScale;
		m_axColumns[3] *= fScale;
	}
	
	void operator/=(const float fScale)
	{
		float fInvScale = 1.0f / fScale;
		m_axColumns[0] *= fInvScale;
		m_axColumns[1] *= fInvScale;
		m_axColumns[2] *= fInvScale;
		m_axColumns[3] *= fInvScale;
	}

	// Conditional operators
	bool operator==(Asura_AMatrix4x4_InParam xOther) const
	{
		bool bCol1Eq = vec_all_eq( m_axColumns[0].m_xVector4, xOther.m_axColumns[0].m_xVector4 );
		bool bCol2Eq = vec_all_eq( m_axColumns[1].m_xVector4, xOther.m_axColumns[1].m_xVector4 );
		bool bCol3Eq = vec_all_eq( m_axColumns[2].m_xVector4, xOther.m_axColumns[2].m_xVector4 );
		bool bCol4Eq = vec_all_eq( m_axColumns[3].m_xVector4, xOther.m_axColumns[3].m_xVector4 );
		return ( bCol1Eq && bCol2Eq && bCol3Eq && bCol4Eq );
	}
	bool ApproximatelyEquals(Asura_AMatrix4x4_InParam xOther, const float fEpsilon = 1.0e-6f) const
	{
		return (m_axColumns[0].ApproximatelyEquals(xOther.m_axColumns[0], fEpsilon))
		&& (m_axColumns[1].ApproximatelyEquals(xOther.m_axColumns[1], fEpsilon))
		&& (m_axColumns[2].ApproximatelyEquals(xOther.m_axColumns[2], fEpsilon))
		&& (m_axColumns[3].ApproximatelyEquals(xOther.m_axColumns[3], fEpsilon));
	}
	bool operator!=(Asura_AMatrix4x4_InParam xOther) const
	{
		bool bCol1Eq = vec_any_ne( m_axColumns[0].m_xVector4, xOther.m_axColumns[0].m_xVector4 );
		bool bCol2Eq = vec_any_ne( m_axColumns[1].m_xVector4, xOther.m_axColumns[1].m_xVector4 );
		bool bCol3Eq = vec_any_ne( m_axColumns[2].m_xVector4, xOther.m_axColumns[2].m_xVector4 );
		bool bCol4Eq = vec_any_ne( m_axColumns[3].m_xVector4, xOther.m_axColumns[3].m_xVector4 );
		return ( bCol1Eq || bCol2Eq || bCol3Eq || bCol4Eq );
	}

	// Elementary properties...
	float GetElement(const u_int uRow, const u_int uCol) const
	{
		Asura_Assert(uRow <= 3, "Invalid matrix row index passed to Asura_AMatrix4x4::GetElement");
		Asura_Assert(uCol <= 3, "Invalid matrix column index passed to Asura_AMatrix4x4::GetElement");
		return m_axColumns[uCol].GetElement( uRow );
	}
	void SetElement(const u_int uRow, const u_int uCol, const float fValue)
	{
		Asura_Assert(uRow <= 3, "Invalid matrix row index passed to Asura_AMatrix4x4::SetElement");
		Asura_Assert(uCol <= 3, "Invalid matrix column index passed to Asura_AMatrix4x4::SetElement");
		m_axColumns[uCol].SetElement( uRow, fValue );
	}

	Asura_AVector4 GetColumn(const u_int uColumnIndex) const
	{
		Asura_Assert(uColumnIndex <= 3, "Invalid matrix column index passed to Asura_AMatrix4x4::GetColumn");
		return m_axColumns[uColumnIndex];
	}

	Asura_AVector4 GetColumn0() const
	{
		return m_axColumns[0];
	}

	Asura_AVector4 GetColumn1() const
	{
		return m_axColumns[1];
	}

		Asura_AVector4 GetColumn2() const
	{
		return m_axColumns[2];
	}

	Asura_AVector4 GetColumn3() const
	{
		return m_axColumns[3];
	}

	void SetColumn(const u_int uColumnIndex, Asura_AVector4_InParam xColumn)
	{
		Asura_Assert(uColumnIndex <= 3, "Invalid matrix column index passed to Asura_AMatrix4x4::SetColumn");
		m_axColumns[uColumnIndex] = xColumn;
	}

	void SetColumn0(Asura_AVector4_InParam xCol0)
	{
		m_axColumns[0] = xCol0;
	}

	void SetColumn1(Asura_AVector4_InParam xCol1)
	{
		m_axColumns[1] = xCol1;
	}

	void SetColumn2(Asura_AVector4_InParam xCol2)
	{
		m_axColumns[2] = xCol2;
	}

	void SetColumn3(Asura_AVector4_InParam xCol3)
	{
		m_axColumns[3] = xCol3;
	}

    const float* GetElements( void ) const { return m_axColumns[0].GetElements(); }

	Asura_AVector4 GetRow(const u_int uRowIndex) const
	{
		Asura_Assert(uRowIndex <= 3, "Invalid matrix row index passed to Asura_AMatrix4x4::GetRow");
		return Asura_AVector4(m_axColumns[0].GetElement(uRowIndex), m_axColumns[1].GetElement(uRowIndex),
							  m_axColumns[2].GetElement(uRowIndex), m_axColumns[3].GetElement(uRowIndex));
	}

	Asura_AVector4 GetRow0() const
	{
		return Asura_AVector4(m_axColumns[0].GetElement(0), m_axColumns[1].GetElement(0),
							  m_axColumns[2].GetElement(0), m_axColumns[3].GetElement(0));
	}

	Asura_AVector4 GetRow1() const
	{
		return Asura_AVector4(m_axColumns[0].GetElement(1), m_axColumns[1].GetElement(1),
							  m_axColumns[2].GetElement(1), m_axColumns[3].GetElement(1));
	}

	Asura_AVector4 GetRow2() const
	{
		return Asura_AVector4(m_axColumns[0].GetElement(2), m_axColumns[1].GetElement(2),
							  m_axColumns[2].GetElement(2), m_axColumns[3].GetElement(2));
	}

	Asura_AVector4 GetRow3() const
	{
		return Asura_AVector4(m_axColumns[0].GetElement(3), m_axColumns[1].GetElement(3),
							  m_axColumns[2].GetElement(3), m_axColumns[3].GetElement(3));
	}

	void SetRow(const u_int uRowIndex, Asura_AVector4_InParam xRow)
	{
		Asura_Assert(uRowIndex <= 3, "Invalid matrix row index passed to Asura_AMatrix4x4::SetRow");
		m_axColumns[0].SetElement(uRowIndex, xRow.GetElement(0));
		m_axColumns[1].SetElement(uRowIndex, xRow.GetElement(1));
		m_axColumns[2].SetElement(uRowIndex, xRow.GetElement(2));
		m_axColumns[3].SetElement(uRowIndex, xRow.GetElement(3));
	}

	void SetRow0(Asura_AVector4_InParam xRow0)
	{
		m_axColumns[0].SetElement(0, xRow0.GetElement(0));
		m_axColumns[1].SetElement(0, xRow0.GetElement(1));
		m_axColumns[2].SetElement(0, xRow0.GetElement(2));
		m_axColumns[3].SetElement(0, xRow0.GetElement(3));
	}

	void SetRow1(Asura_AVector4_InParam xRow1)
	{
		m_axColumns[0].SetElement(1, xRow1.GetElement(0));
		m_axColumns[1].SetElement(1, xRow1.GetElement(1));
		m_axColumns[2].SetElement(1, xRow1.GetElement(2));
		m_axColumns[3].SetElement(1, xRow1.GetElement(3));
	}

	void SetRow2(Asura_AVector4_InParam xRow2)
	{
		m_axColumns[0].SetElement(2, xRow2.GetElement(0));
		m_axColumns[1].SetElement(2, xRow2.GetElement(1));
		m_axColumns[2].SetElement(2, xRow2.GetElement(2));
		m_axColumns[3].SetElement(2, xRow2.GetElement(3));
	}

	void SetRow3(Asura_AVector4_InParam xRow3)
	{
		m_axColumns[0].SetElement(3, xRow3.GetElement(0));
		m_axColumns[1].SetElement(3, xRow3.GetElement(1));
		m_axColumns[2].SetElement(3, xRow3.GetElement(2));
		m_axColumns[3].SetElement(3, xRow3.GetElement(3));
	}

	void SetColumns(Asura_AVector4_InParam xCol0, Asura_AVector4_InParam xCol1, Asura_AVector4_InParam xCol2, Asura_AVector4_InParam xCol3)
	{
		m_axColumns[0]	= xCol0;
		m_axColumns[1]	= xCol1;
		m_axColumns[2]	= xCol2;
		m_axColumns[3]	= xCol3;
	}
	void SetRows(Asura_AVector4_InParam xRow0, Asura_AVector4_InParam xRow1, Asura_AVector4_InParam xRow2, Asura_AVector4_InParam xRow3)
	{
		SetRow0(xRow0);
		SetRow1(xRow1);
		SetRow2(xRow2);
		SetRow3(xRow3);
	}

	// Load/store functions...
	void LoadColsAligned(const float* pfFloats)
	{
		m_axColumns[0].LoadAligned(&pfFloats[0]);
		m_axColumns[1].LoadAligned(&pfFloats[4]);
		m_axColumns[2].LoadAligned(&pfFloats[8]);
		m_axColumns[3].LoadAligned(&pfFloats[12]);
	}
	void LoadRowsAligned(const float* pfFloats)
	{
		LoadColsAligned(pfFloats);
		Transpose();
	}
	void LoadColsUnaligned(const float* pfFloats)
	{
		m_axColumns[0].LoadUnaligned(&pfFloats[0]);
		m_axColumns[1].LoadUnaligned(&pfFloats[4]);
		m_axColumns[2].LoadUnaligned(&pfFloats[8]);
		m_axColumns[3].LoadUnaligned(&pfFloats[12]);
	}
	void LoadRowsUnaligned(const float* pfFloats)
	{
		LoadColsUnaligned(pfFloats);
		Transpose();
	}

	void StoreColsAligned(float* pfFloats) const
	{
		m_axColumns[0].StoreAligned(&pfFloats[0]);
		m_axColumns[1].StoreAligned(&pfFloats[4]);
		m_axColumns[2].StoreAligned(&pfFloats[8]);
		m_axColumns[3].StoreAligned(&pfFloats[12]);
	}
	void StoreRowsAligned(float* pfFloats) const
	{
		Asura_AMatrix4x4	xTemp;

		Transpose(xTemp);
		xTemp.StoreColsAligned(pfFloats);
	}
	void StoreColsUnaligned(float* pfFloats) const
	{
		m_axColumns[0].StoreUnaligned(&pfFloats[0]);
		m_axColumns[1].StoreUnaligned(&pfFloats[4]);
		m_axColumns[2].StoreUnaligned(&pfFloats[8]);
		m_axColumns[3].StoreUnaligned(&pfFloats[12]);
	}
	void StoreRowsUnaligned(float* pfFloats) const
	{
		Asura_AMatrix4x4	xTemp;

		Transpose(xTemp);
		xTemp.StoreColsUnaligned(pfFloats);
	}

	void SetPerspectiveForLH( float fFovyRadians, float fAspect, float fZNear, float fZFar ) //match the Xbox360 documentation
	{
		float fH = tanf( PI_OVER_2 - fFovyRadians * 0.5f );
		float fW = fH / fAspect;
		float fRangeInv = 1.0f / ( fZFar - fZNear );
		float fT1 = fZFar * fRangeInv;
		float fT2 = -( fZNear * fZFar * fRangeInv );	

		m_axColumns[ 0 ].SetElements( fW,	0.0f,	0.0f,	0.0f );
		m_axColumns[ 1 ].SetElements( 0.0f,	fH,		0.0f,	0.0f );
		m_axColumns[ 2 ].SetElements( 0.0f,	0.0f,	fT1,	1.0f );
		m_axColumns[ 3 ].SetElements( 0.0f,	0.0f,	fT2,	0.0f );
	}

	void SetPerspectiveForRH( float fFovyRadians, float fAspect, float fZNear, float fZFar ) //match the Xbox360 documentation
	{
		float fH = tanf( PI_OVER_2 - fFovyRadians * 0.5f );
		float fW = fH / fAspect;
		float fRangeInv = 1.0f / ( fZFar - fZNear );
		float fT1 = fZFar * fRangeInv;
		float fT2 = fZNear * fZFar * fRangeInv;

		m_axColumns[ 0 ].SetElements( fW,	0.0f,	0.0f,	0.0f );
		m_axColumns[ 1 ].SetElements( 0.0f,	fH,		0.0f,	0.0f );
		m_axColumns[ 2 ].SetElements( 0.0f,	0.0f,	fT1,	fT2 );
		m_axColumns[ 3 ].SetElements( 0.0f,	0.0f,	-1.0f,	0.0f );
	}
	
	void SetOrthographicForLH( float fW, float fH, float fZNear, float fZFar )
	{
		float fRangeInv = 1.0f / ( fZFar - fZNear );

		m_axColumns[ 0 ].SetElements( 2.0f/fW,	0.0f,		0.0f,				0.0f ); // [-fW/2,fW/2] -> [-1,1]
		m_axColumns[ 1 ].SetElements( 0.0f,		2.0f/fH,	0.0f,				0.0f );
		m_axColumns[ 2 ].SetElements( 0.0f,		0.0f,		fRangeInv,			0.0f );
		m_axColumns[ 3 ].SetElements( 0.0f,		0.0f,		-fZNear*fRangeInv,	1.0f );
	}

private:
	Asura_AVector4 m_axColumns[4];
};


//////////////////////////////////////////////////////////////////////////
// Inline Function Definitions
//////////////////////////////////////////////////////////////////////////

inline Asura_AVector3 operator*( Asura_AVector3_InParam xVec, Asura_AMatrix3x3_InParam xMat )
{
	return xMat.Apply( xVec );
}

inline Asura_AVector3 operator*( Asura_AVector3_InParam xVec, const Asura_Matrix_3x3& xMat )
{
	const Asura_AMatrix3x3 xTemp( xMat );
	return xTemp.Apply( xVec );
}

inline Asura_AMatrix3x3 operator*( Asura_AMatrix3x3_InParam xM1, Asura_AMatrix3x3_InParam xM2 )
{
	Asura_AMatrix3x3 xResult;

	xM2.Apply( xM1, xResult );
	return xResult;
}

// scalar * 4x1 Vector - to match the member function in Asura_AVector4
inline Asura_AVector4 operator*(const float fScale, const Asura_AVector4& xOther)
{
	return(xOther * fScale);
}

inline Asura_AVector4 operator*(Asura_AVector4_InParam xVec, Asura_AMatrix4x4_InParam xMat)
{
	return xMat.Apply(xVec);
}

// To match the backwards vector-matrix multiply, we do a backwards matrix-vector multiply too.
inline Asura_AMatrix4x4 operator*(Asura_AMatrix4x4_InParam xM1, Asura_AMatrix4x4_InParam xM2)
{
	Asura_AMatrix4x4	xResult;

	xM2.Apply(xM1, xResult);
	return(xResult);
}

#endif	// defined(ASURA_USE_PLATFORM_VECTORS)
// E N D   F I L E   W R A P P E R ///////////////////////////////////////////

#endif // #ifndef ___ASURA_PS3_VECTORMATHS_H___
