PageRenderTime 129ms CodeModel.GetById 27ms app.highlight 87ms RepoModel.GetById 2ms app.codeStats 0ms

/Utilities/Datatypes/Matrix.cs

#
C# | 1820 lines | 1114 code | 118 blank | 588 comment | 86 complexity | 02a0f05e85a992c0b9fa7d8b4153d0c5 MD5 | raw file
   1using System;
   2using System.IO;
   3using System.Runtime.InteropServices;
   4using Delta.Utilities.Helpers;
   5using Delta.Utilities.Profiling;
   6using NUnit.Framework;
   7
   8
   9// Disable warnings for some tests where we don't use created values
  10#pragma warning disable 219
  11// Disable the obsolete warnings for Matrix.Multiply here, it is optimized
  12// already here and when it is not in some unit test, its ok too.
  13#pragma warning disable 618
  14
  15namespace Delta.Utilities.Datatypes
  16{
  17	/// <summary>
  18	/// Matrix struct, used mostly for 3D calculations. If possible this shares
  19	/// all the functionality of the XNA Matrix struct (if we are compiling
  20	/// only). In case we do not compile with XNA or SlimDX enabled (like on
  21	/// this class still provides all the basic functionality as a fallback. It
  22	/// might not be as optimized as some specialized xna code paths for the
  23	/// Xbox 360, but those only work on the XNA define and platform anyway.
  24	/// </summary>
  25	/// <remarks>
  26	/// Represents a row based matrix. Good wiki page about the coordinate
  27	/// system (right handed):
  28	/// http://en.wikipedia.org/wiki/Cartesian_coordinate_system
  29	/// <para />
  30	/// Note: This Matrix class was heavily profiled and optimized with help of
  31	/// the MatrixPerformance helper class included here. Many code paths
  32	/// are platform and framework dependant, this is why this class is so huge,
  33	/// but since this is very low level and highly performance critical, it is
  34	/// well worth the complexity and effort (the class is still easy to use).
  35	/// On the Xbox we will always try to use XNA implementation for performance
  36	/// critical methods, all other platforms are mostly implemented by us
  37	/// after testing performance. You can also checkout the SlimDX math classes,
  38	/// which mostly are just C#, only certain heavy methods like multiplying
  39	/// an array of matrices uses PInvoke to D3DXMatrixMultiply for example.
  40	/// <para />
  41	/// There are still performance improvements possible for this class
  42	/// like using ref for most performance critical methods like Invert,
  43	/// Multiply and most other operators, but that makes this class much harder
  44	/// to use and we can still do code-optimizations directly with the exposed
  45	/// data when needed (e.g. Bone Matrix code should be heavily optimized if
  46	/// it has to run on CPU). Currently this is fixed by marking slow methods
  47	/// as obsolete, which will give you compiler warnings, use the ref
  48	/// overload methods instead.
  49	/// <para />
  50	/// Also note that many methods are not faster by using the native XNA,
  51	/// OpenTK or SharpDX or SlimDX code paths, especially when just assigning
  52	/// values (where our code is almost always faster). But remember when any
  53	/// framework is using special tricks like low level assembly code or native
  54	/// code to calculate complex math operations, it can be faster and will be
  55	/// replaced by the build system automatically. If you notice any method
  56	/// that should be faster, contact us immediately so we can investigate!
  57	/// </remarks>
  58	[Serializable]
  59	[StructLayout(LayoutKind.Explicit)]
  60	public struct Matrix : ISaveLoadBinary, IEquatable<Matrix>
  61	{
  62		#region Constants
  63		/// <summary>
  64		/// Represents the number of values of that Matrix struct (here 4x4 = 16).
  65		/// </summary>
  66		public const int ValueCount = 16;
  67
  68		/// <summary>
  69		/// Returns a zeroed matrix.
  70		/// </summary>
  71		public static readonly Matrix Zero = new Matrix(
  72			0, 0, 0, 0,
  73			0, 0, 0, 0,
  74			0, 0, 0, 0,
  75			0, 0, 0, 0);
  76		#endregion
  77
  78		#region Invert (Static)
  79		/// <summary>
  80		/// Invert matrix, this is only slightly faster on the xna platform
  81		/// specific implementation, for performance critical code use the ref
  82		/// version!
  83		/// </summary>
  84		/// <param name="a">Matrix to invert</param>
  85		/// <returns>Inverted Matrix</returns>
  86		public static Matrix Invert(Matrix a)
  87		{
  88			//Check performance difference again:
  89
  90			//Note: Performance results from 10 million calls:
  91			//Delta: 1962ms
  92			//XNA: 1076ms
  93			//OpenTK: 8611ms
  94			//SlimDX: 981ms
  95
  96			// Ugly code that is actually lots faster than our clean implementation
  97			// below (see link on how it works basically). This is basically getting
  98			// the Determinant and then also calculating the inverse matrix in one
  99			// big ugly swoop.
 100			float m11 = a.M11;
 101			float m12 = a.M12;
 102			float m13 = a.M13;
 103			float m14 = a.M14;
 104			float m21 = a.M21;
 105			float m22 = a.M22;
 106			float m23 = a.M23;
 107			float m24 = a.M24;
 108			float m31 = a.M31;
 109			float m32 = a.M32;
 110			float m33 = a.M33;
 111			float m34 = a.M34;
 112			float m41 = a.M41;
 113			float m42 = a.M42;
 114			float m43 = a.M43;
 115			float m44 = a.M44;
 116			float m33x44m34x43 = (m33 * m44) - (m34 * m43);
 117			float m32x44m34x42 = (m32 * m44) - (m34 * m42);
 118			float m32x43m33x42 = (m32 * m43) - (m33 * m42);
 119			float m31x44m34x41 = (m31 * m44) - (m34 * m41);
 120			float m31x43m33x41 = (m31 * m43) - (m33 * m41);
 121			float m31x42m32x41 = (m31 * m42) - (m32 * m41);
 122			float m22m23m24 = ((m22 * m33x44m34x43) - (m23 * m32x44m34x42)) +
 123			                  (m24 * m32x43m33x42);
 124			float m21m23m24 = -(((m21 * m33x44m34x43) - (m23 * m31x44m34x41)) +
 125			                    (m24 * m31x43m33x41));
 126			float m21m22m24 = ((m21 * m32x44m34x42) - (m22 * m31x44m34x41)) +
 127			                  (m24 * m31x42m32x41);
 128			float m21m22m23 = -(((m21 * m32x43m33x42) - (m22 * m31x43m33x41)) +
 129			                    (m23 * m31x42m32x41));
 130
 131			// Compute the the determinant
 132			float det =
 133				(m11 * m22m23m24) + (m12 * m21m23m24) + (m13 * m21m22m24) +
 134				(m14 * m21m22m23);
 135			if (det == 0.0f)
 136			{
 137				det = 0.000001f;
 138			}
 139
 140			// Inverse(A) = 1/det(A) * B
 141			float inverseDet = 1f / det;
 142			float m23x44m24x43 = (m23 * m44) - (m24 * m43);
 143			float m22x44m24x42 = (m22 * m44) - (m24 * m42);
 144			float m22x43m23x42 = (m22 * m43) - (m23 * m42);
 145			float m21x44m24x41 = (m21 * m44) - (m24 * m41);
 146			float m21x43m23x41 = (m21 * m43) - (m23 * m41);
 147			float m21x42m22x41 = (m21 * m42) - (m22 * m41);
 148			float m23x34m24x33 = (m23 * m34) - (m24 * m33);
 149			float m22x34m24x32 = (m22 * m34) - (m24 * m32);
 150			float m22x33m23x32 = (m22 * m33) - (m23 * m32);
 151			float m21x34m23x31 = (m21 * m34) - (m24 * m31);
 152			float m21x33m23x31 = (m21 * m33) - (m23 * m31);
 153			float m21x32m22x31 = (m21 * m32) - (m22 * m31);
 154			return new Matrix(
 155				m22m23m24 * inverseDet,
 156				-(((m12 * m33x44m34x43) - (m13 * m32x44m34x42)) +
 157				  (m14 * m32x43m33x42)) * inverseDet,
 158				(((m12 * m23x44m24x43) - (m13 * m22x44m24x42)) +
 159				 (m14 * m22x43m23x42)) * inverseDet,
 160				-(((m12 * m23x34m24x33) - (m13 * m22x34m24x32)) +
 161				  (m14 * m22x33m23x32)) * inverseDet,
 162				m21m23m24 * inverseDet,
 163				(((m11 * m33x44m34x43) - (m13 * m31x44m34x41)) +
 164				 (m14 * m31x43m33x41)) * inverseDet,
 165				-(((m11 * m23x44m24x43) - (m13 * m21x44m24x41)) +
 166				  (m14 * m21x43m23x41)) * inverseDet,
 167				(((m11 * m23x34m24x33) - (m13 * m21x34m23x31)) +
 168				 (m14 * m21x33m23x31)) * inverseDet,
 169				m21m22m24 * inverseDet,
 170				-(((m11 * m32x44m34x42) - (m12 * m31x44m34x41)) +
 171				  (m14 * m31x42m32x41)) * inverseDet,
 172				(((m11 * m22x44m24x42) - (m12 * m21x44m24x41)) +
 173				 (m14 * m21x42m22x41)) * inverseDet,
 174				-(((m11 * m22x34m24x32) - (m12 * m21x34m23x31)) +
 175				  (m14 * m21x32m22x31)) * inverseDet,
 176				m21m22m23 * inverseDet,
 177				(((m11 * m32x43m33x42) - (m12 * m31x43m33x41)) +
 178				 (m13 * m31x42m32x41)) * inverseDet,
 179				-(((m11 * m22x43m23x42) - (m12 * m21x43m23x41)) +
 180				  (m13 * m21x42m22x41)) * inverseDet,
 181				(((m11 * m22x33m23x32) - (m12 * m21x33m23x31)) +
 182				 (m13 * m21x32m22x31)) * inverseDet);
 183		}
 184
 185		/// <summary>
 186		/// Invert matrix ref version, this is the fastest one, especially if
 187		/// you use matrix and result for the same value, which will copy values
 188		/// over quickly.
 189		/// </summary>
 190		/// <param name="matrix">Matrix</param>
 191		/// <param name="result">Result</param>
 192		public static void Invert(ref Matrix matrix, ref Matrix result)
 193		{
 194			//Check performance difference again:
 195
 196			// Ugly code that is actually lots faster than our clean implementation
 197			// below (see link on how it works basically). This is basically getting
 198			// the Determinant and then also calculating the inverse matrix in one
 199			// big ugly swoop.
 200			float m11 = matrix.M11;
 201			float m12 = matrix.M12;
 202			float m13 = matrix.M13;
 203			float m14 = matrix.M14;
 204			float m21 = matrix.M21;
 205			float m22 = matrix.M22;
 206			float m23 = matrix.M23;
 207			float m24 = matrix.M24;
 208			float m31 = matrix.M31;
 209			float m32 = matrix.M32;
 210			float m33 = matrix.M33;
 211			float m34 = matrix.M34;
 212			float m41 = matrix.M41;
 213			float m42 = matrix.M42;
 214			float m43 = matrix.M43;
 215			float m44 = matrix.M44;
 216			float m33x44m34x43 = (m33 * m44) - (m34 * m43);
 217			float m32x44m34x42 = (m32 * m44) - (m34 * m42);
 218			float m32x43m33x42 = (m32 * m43) - (m33 * m42);
 219			float m31x44m34x41 = (m31 * m44) - (m34 * m41);
 220			float m31x43m33x41 = (m31 * m43) - (m33 * m41);
 221			float m31x42m32x41 = (m31 * m42) - (m32 * m41);
 222			float m22m23m24 = ((m22 * m33x44m34x43) - (m23 * m32x44m34x42)) +
 223			                  (m24 * m32x43m33x42);
 224			float m21m23m24 = -(((m21 * m33x44m34x43) - (m23 * m31x44m34x41)) +
 225			                    (m24 * m31x43m33x41));
 226			float m21m22m24 = ((m21 * m32x44m34x42) - (m22 * m31x44m34x41)) +
 227			                  (m24 * m31x42m32x41);
 228			float m21m22m23 = -(((m21 * m32x43m33x42) - (m22 * m31x43m33x41)) +
 229			                    (m23 * m31x42m32x41));
 230
 231			// Computing the determinant
 232			float det = (m11 * m22m23m24) + (m12 * m21m23m24) + (m13 * m21m22m24) +
 233			            (m14 * m21m22m23);
 234			if (det == 0.0f)
 235			{
 236				det = 0.000001f;
 237			}
 238
 239			// Inverse(A) = 1/det(A) * B
 240			float inverseDet = 1f / det;
 241			float m23x44m24x43 = (m23 * m44) - (m24 * m43);
 242			float m22x44m24x42 = (m22 * m44) - (m24 * m42);
 243			float m22x43m23x42 = (m22 * m43) - (m23 * m42);
 244			float m21x44m24x41 = (m21 * m44) - (m24 * m41);
 245			float m21x43m23x41 = (m21 * m43) - (m23 * m41);
 246			float m21x42m22x41 = (m21 * m42) - (m22 * m41);
 247			float m23x34m24x33 = (m23 * m34) - (m24 * m33);
 248			float m22x34m24x32 = (m22 * m34) - (m24 * m32);
 249			float m22x33m23x32 = (m22 * m33) - (m23 * m32);
 250			float m21x34m23x31 = (m21 * m34) - (m24 * m31);
 251			float m21x33m23x31 = (m21 * m33) - (m23 * m31);
 252			float m21x32m22x31 = (m21 * m32) - (m22 * m31);
 253			// Now just set everything to the output result matrix, which even can
 254			// be the same as the input matrix (for best performance).
 255			result.M11 = m22m23m24 * inverseDet;
 256			result.M12 = -(((m12 * m33x44m34x43) - (m13 * m32x44m34x42)) +
 257			               (m14 * m32x43m33x42)) * inverseDet;
 258			result.M13 = (((m12 * m23x44m24x43) - (m13 * m22x44m24x42)) +
 259			              (m14 * m22x43m23x42)) * inverseDet;
 260			result.M14 = -(((m12 * m23x34m24x33) - (m13 * m22x34m24x32)) +
 261			               (m14 * m22x33m23x32)) * inverseDet;
 262			result.M21 = m21m23m24 * inverseDet;
 263			result.M22 = (((m11 * m33x44m34x43) - (m13 * m31x44m34x41)) +
 264			              (m14 * m31x43m33x41)) * inverseDet;
 265			result.M23 = -(((m11 * m23x44m24x43) - (m13 * m21x44m24x41)) +
 266			               (m14 * m21x43m23x41)) * inverseDet;
 267			result.M24 = (((m11 * m23x34m24x33) - (m13 * m21x34m23x31)) +
 268			              (m14 * m21x33m23x31)) * inverseDet;
 269			result.M31 = m21m22m24 * inverseDet;
 270			result.M32 = -(((m11 * m32x44m34x42) - (m12 * m31x44m34x41)) +
 271			               (m14 * m31x42m32x41)) * inverseDet;
 272			result.M33 = (((m11 * m22x44m24x42) - (m12 * m21x44m24x41)) +
 273			              (m14 * m21x42m22x41)) * inverseDet;
 274			result.M34 = -(((m11 * m22x34m24x32) - (m12 * m21x34m23x31)) +
 275			               (m14 * m21x32m22x31)) * inverseDet;
 276			result.M41 = m21m22m23 * inverseDet;
 277			result.M42 = (((m11 * m32x43m33x42) - (m12 * m31x43m33x41)) +
 278			              (m13 * m31x42m32x41)) * inverseDet;
 279			result.M43 = -(((m11 * m22x43m23x42) - (m12 * m21x43m23x41)) +
 280			               (m13 * m21x42m22x41)) * inverseDet;
 281			result.M44 = (((m11 * m22x33m23x32) - (m12 * m21x33m23x31)) +
 282			              (m13 * m21x32m22x31)) * inverseDet;
 283		}
 284		#endregion
 285
 286		#region Transpose (Static)
 287		/// <summary>
 288		/// Transpose, which basically just means we switch from a row to a column
 289		/// based matrix.
 290		/// </summary>
 291		/// <param name="matrix">Matrix</param>
 292		/// <returns>Transposed matrix (row/columns switched)</returns>
 293		public static Matrix Transpose(Matrix matrix)
 294		{
 295			//Check performance difference again:
 296
 297			// Note: In all cases our implementation is slighly faster/or a lot faster
 298			//Delta: 569ms (10 mio calls)
 299			//XNA: 760ms (10 mio calls)
 300			//OpenTK: 904ms (10 mio calls)
 301			//SlimDX: 658ms (10 mio calls)
 302
 303			return new Matrix(
 304				// First row
 305				matrix.M11, matrix.M21, matrix.M31, matrix.M41,
 306				// Second row
 307				matrix.M12, matrix.M22, matrix.M32, matrix.M42,
 308				// Third row
 309				matrix.M13, matrix.M23, matrix.M33, matrix.M43,
 310				// Forth row
 311				matrix.M14, matrix.M24, matrix.M34, matrix.M44);
 312		}
 313
 314		/// <summary>
 315		/// Transpose, which basically just means we switch from a row to a column
 316		/// based matrix. Faster version using ref input and output values.
 317		/// </summary>
 318		/// <param name="matrix">Matrix</param>
 319		/// <param name="result">result</param>
 320		public static void Transpose(ref Matrix matrix, ref Matrix result)
 321		{
 322			result.M11 = matrix.M11;
 323			result.M12 = matrix.M21;
 324			result.M13 = matrix.M31;
 325			result.M14 = matrix.M41;
 326			result.M21 = matrix.M12;
 327			result.M22 = matrix.M22;
 328			result.M23 = matrix.M32;
 329			result.M24 = matrix.M42;
 330			result.M31 = matrix.M13;
 331			result.M32 = matrix.M23;
 332			result.M33 = matrix.M33;
 333			result.M34 = matrix.M43;
 334			result.M41 = matrix.M14;
 335			result.M42 = matrix.M24;
 336			result.M43 = matrix.M34;
 337			result.M44 = matrix.M44;
 338		}
 339		#endregion
 340
 341		#region CreateColumnBased (Static)
 342		/// <summary>
 343		/// Create column based matrix, also used for creating an exported collada
 344		/// matrix, because the normal XNA/OpenTK/whatever matrix is row based.
 345		/// </summary>
 346		/// <param name="floatValues">Float values</param>
 347		/// <returns>Matrix</returns>
 348		public static Matrix CreateColumnBased(float[] floatValues)
 349		{
 350			// If we get no valid or an incorrect number of float values, then we
 351			// return only the identity matrix
 352			if (floatValues == null ||
 353			    floatValues.Length != 16)
 354			{
 355				return Identity;
 356			}
 357
 358			return new Matrix(
 359				floatValues[0], floatValues[4], floatValues[8], floatValues[12],
 360				floatValues[1], floatValues[5], floatValues[9], floatValues[13],
 361				floatValues[2], floatValues[6], floatValues[10], floatValues[14],
 362				floatValues[3], floatValues[7], floatValues[11], floatValues[15]);
 363		}
 364		#endregion
 365
 366		#region CreateScale (Static)
 367		/// <summary>
 368		/// Create a scale matrix with the specified value for X, Y and Z.
 369		/// </summary>
 370		/// <param name="scale">Scale</param>
 371		/// <returns>New Matrix with the given scale</returns>
 372		public static Matrix CreateScale(float scale)
 373		{
 374			return new Matrix(
 375				scale, 0f, 0f, 0f,
 376				0f, scale, 0f, 0f,
 377				0f, 0f, scale, 0f,
 378				0f, 0f, 0f, 1.0f);
 379		}
 380
 381		/// <summary>
 382		/// Create scale
 383		/// </summary>
 384		/// <param name="scale">Scale</param>
 385		/// <returns>New Matrix with the given scale</returns>
 386		public static Matrix CreateScale(Vector scale)
 387		{
 388			return CreateScale(scale.X, scale.Y, scale.Z);
 389		}
 390
 391		/// <summary>
 392		/// Create a scale matrix with the specified X, Y and Z scaling values.
 393		/// </summary>
 394		/// <param name="scaleX">The amount of scaling in X dimension.</param>
 395		/// <param name="scaleY">The amount of scaling in Y dimension.</param>
 396		/// <param name="scaleZ">The amount of scaling in Z dimension.</param>
 397		/// <returns>New Matrix with the given scale</returns>
 398		public static Matrix CreateScale(float scaleX, float scaleY, float scaleZ)
 399		{
 400			return new Matrix(
 401				scaleX, 0f, 0f, 0f,
 402				0f, scaleY, 0f, 0f,
 403				0f, 0f, scaleZ, 0f,
 404				0f, 0f, 0f, 1f);
 405		}
 406		#endregion
 407
 408		#region CreateRotationX (Static)
 409		/// <summary>
 410		/// Create rotation around x-axis (pitch)
 411		/// <para />
 412		/// It is important to know which axis is meant with Yaw, Pitch and Roll:
 413		/// http://www.spotimage.com/dimap/spec/dictionary/Spot_Scene/Illustrations/YAW.gif
 414		/// </summary>
 415		public static Matrix CreateRotationX(float degrees)
 416		{
 417			//Check performance difference again:
 418
 419			// Note: All CreateRotation implementations are slow, there are slightly
 420			// faster on XNA, but overall not much of a difference between all
 421			// platforms and frameworks (see below in MatrixPerformance class).
 422
 423			float cosValue = MathHelper.Cos(degrees);
 424			float sinValue = MathHelper.Sin(degrees);
 425			return new Matrix(
 426				1f, 0f, 0f, 0f,
 427				0f, cosValue, sinValue, 0f,
 428				0f, -sinValue, cosValue, 0f,
 429				0f, 0f, 0f, 1f);
 430		}
 431		#endregion
 432
 433		#region CreateRotationY (Static)
 434		/// <summary>
 435		/// Create rotation around y-axis (yaw)
 436		/// <para />
 437		/// It is important to know which axis is meant with Yaw, Pitch and Roll:
 438		/// http://www.spotimage.com/dimap/spec/dictionary/Spot_Scene/Illustrations/YAW.gif
 439		/// </summary>
 440		/// <param name="degrees">Degrees to rotate around the Y axis</param>
 441		/// <returns>New matrix with the given rotation</returns>
 442		public static Matrix CreateRotationY(float degrees)
 443		{
 444			//Check performance difference again:
 445			// Note: All CreateRotation implementations are slow, there are slightly
 446			// faster on XNA, but overall not much of a difference between all
 447			// platforms and frameworks (see below in MatrixPerformance class).
 448
 449			float cosValue = MathHelper.Cos(degrees);
 450			float sinValue = MathHelper.Sin(degrees);
 451			return new Matrix(
 452				cosValue, 0f, -sinValue, 0f,
 453				0f, 1f, 0f, 0f,
 454				sinValue, 0f, cosValue, 0f,
 455				0f, 0f, 0f, 1f);
 456		}
 457		#endregion
 458
 459		#region CreateRotationZ (Static)
 460		/// <summary>
 461		/// Create rotation around z-axis (roll)
 462		/// <para />
 463		/// It is important to know which axis is meant with Yaw, Pitch and Roll:
 464		/// http://www.spotimage.com/dimap/spec/dictionary/Spot_Scene/Illustrations/YAW.gif
 465		/// </summary>
 466		/// <param name="degrees">Degrees to rotate around the Z axis</param>
 467		/// <returns>New matrix with the given rotation</returns>
 468		public static Matrix CreateRotationZ(float degrees)
 469		{
 470			//Check performance difference again:
 471
 472			// Note: All CreateRotation implementations are slow, there are slightly
 473			// faster on XNA, but overall not much of a difference between all
 474			// platforms and frameworks (see below in MatrixPerformance class).
 475
 476			float cosValue = MathHelper.Cos(degrees);
 477			float sinValue = MathHelper.Sin(degrees);
 478			return new Matrix(
 479				cosValue, sinValue, 0f, 0f,
 480				-sinValue, cosValue, 0f, 0f,
 481				0f, 0f, 1f, 0f,
 482				0f, 0f, 0f, 1f);
 483		}
 484		#endregion
 485
 486		#region CreateRotationZYX (Static)
 487		/// <summary>
 488		/// Rotate around the Z axis, then the Y axis, and finally the X axis
 489		/// (rotX * rotY * rotZ). Please note that this is not the same as
 490		/// FromYawPitchRoll, which uses Z (yaw), then X (pitch), then Y (roll).
 491		/// <para />
 492		/// Note: This is the same as calling the CreateRotationX, Y, Z methods,
 493		/// just combined all together without having to multiply matrices).
 494		/// </summary>
 495		/// <param name="x">Degrees to rotate around the X axis</param>
 496		/// <param name="y">Degrees to rotate around the Y axis</param>
 497		/// <param name="z">Degrees to rotate around the Z axis</param>
 498		/// <returns>New matrix with the given rotation</returns>
 499		public static Matrix CreateRotationZYX(float x, float y, float z)
 500		{
 501			float cx = 1;
 502			float sx = 0;
 503			float cy = 1;
 504			float sy = 0;
 505			float cz = 1;
 506			float sz = 0;
 507			if (x != 0.0f)
 508			{
 509				cx = MathHelper.Cos(x);
 510				sx = MathHelper.Sin(x);
 511			}
 512			if (y != 0.0f)
 513			{
 514				cy = MathHelper.Cos(y);
 515				sy = MathHelper.Sin(y);
 516			}
 517			if (z != 0.0f)
 518			{
 519				cz = MathHelper.Cos(z);
 520				sz = MathHelper.Sin(z);
 521			}
 522
 523			return new Matrix(
 524				cy * cz, cy * sz, -sy, 0.0f,
 525				(sx * sy * cz) + (cx * -sz), (sx * sy * sz) + (cx * cz), sx * cy, 0.0f,
 526				(cx * sy * cz) + (sx * sz), (cx * sy * sz) + (-sx * cz), cx * cy, 0.0f,
 527				0.0f, 0.0f, 0.0f, 1.0f);
 528		}
 529		#endregion
 530
 531		#region CreateRotationZY (Static)
 532		/// <summary>
 533		/// Rotate around Z axis, then around Y axis.
 534		/// (equals rotY * rotZ)
 535		/// </summary>
 536		/// <param name="y">Degrees to rotate around the Y axis</param>
 537		/// <param name="z">Degrees to rotate around the Z axis</param>
 538		/// <returns>New matrix with the given rotation</returns>
 539		public static Matrix CreateRotationZY(float y, float z)
 540		{
 541			float cy = 1;
 542			float sy = 0;
 543			float cz = 1;
 544			float sz = 0;
 545			// We only calculate sin/cos if given input is not zero
 546			if (y != 0.0f)
 547			{
 548				cy = MathHelper.Cos(y);
 549				sy = MathHelper.Sin(y);
 550			}
 551			if (z != 0.0f)
 552			{
 553				cz = MathHelper.Cos(z);
 554				sz = MathHelper.Sin(z);
 555			}
 556			return new Matrix(
 557				cy * cz, cy * sz, -sy, 0.0f,
 558				-sz, cz, 0.0f, 0.0f,
 559				sy * cz, sy * sz, cy, 0.0f,
 560				0.0f, 0.0f, 0.0f, 1.0f);
 561		}
 562		#endregion
 563
 564		#region FromYawPitchRoll (Static)
 565		/// <summary>
 566		/// From yaw pitch roll (yaw = Y, pitch = X, roll = Z)
 567		/// http://www.euclideanspace.com/maths/geometry/rotations/conversions/eulerToMatrix/index.htm
 568		/// It is important to know which axis is meant with Yaw, Pitch and Roll:
 569		/// http://www.spotimage.com/dimap/spec/dictionary/Spot_Scene/Illustrations/YAW.gif
 570		/// <para />
 571		/// Note: This method calculates Y * X * Z where
 572		/// Y = value.X = yaw
 573		/// X = value.Y = pitch
 574		/// Z = value.Z = roll
 575		/// </summary>
 576		public static Matrix FromYawPitchRoll(Vector value)
 577		{
 578			return FromYawPitchRoll(value.X, value.Y, value.Z);
 579		}
 580
 581		/// <summary>
 582		/// From yaw pitch roll (yaw = Y, pitch = X, roll = Z)
 583		/// It is important to know which axis is meant with Yaw, Pitch and Roll:
 584		/// http://www.spotimage.com/dimap/spec/dictionary/Spot_Scene/Illustrations/YAW.gif
 585		/// <para />
 586		/// Note: This method calculates Y * X * Z where
 587		/// Y = yaw
 588		/// X = pitch
 589		/// Z = roll
 590		/// </summary>
 591		public static Matrix FromYawPitchRoll(float yaw, float pitch, float roll)
 592		{
 593			//Check performance difference again:
 594
 595			float cy = MathHelper.Cos(yaw);
 596			float sy = MathHelper.Sin(yaw);
 597			float cp = MathHelper.Cos(pitch);
 598			float sp = MathHelper.Sin(pitch);
 599			float cr = MathHelper.Cos(roll);
 600			float sr = MathHelper.Sin(roll);
 601
 602			return new Matrix(
 603				cy * cr - sr * sy * sp, cy * sr + sy * sp * cr, -sy * cp, 0.0f,
 604				-sr * cp, cp * cr, sp, 0.0f,
 605				sy * cr + sp * cy * sr, sy * sr - sp * cy * cr, cy * cp, 0.0f,
 606				0.0f, 0.0f, 0.0f, 1.0f);
 607		}
 608		#endregion
 609
 610		#region CreateTranslation (Static)
 611		/// <summary>
 612		/// Create translation matrix
 613		/// </summary>
 614		/// <param name="position">Position</param>
 615		/// <returns>Identity matrix with the translation (M41, M42, M43)</returns>
 616		public static Matrix CreateTranslation(Vector position)
 617		{
 618			return new Matrix(
 619				1f, 0f, 0f, 0f,
 620				0f, 1f, 0f, 0f,
 621				0f, 0f, 1f, 0f,
 622				position.X, position.Y, position.Z, 1.0f);
 623		}
 624
 625		/// <summary>
 626		/// Create translation matrix
 627		/// </summary>
 628		/// <param name="x">x</param>
 629		/// <param name="y">y</param>
 630		/// <param name="z">z</param>
 631		/// <returns>Identity matrix with the translation (M41, M42, M43)</returns>
 632		public static Matrix CreateTranslation(float x, float y, float z)
 633		{
 634			return new Matrix(
 635				1f, 0f, 0f, 0f,
 636				0f, 1f, 0f, 0f,
 637				0f, 0f, 1f, 0f,
 638				x, y, z, 1.0f);
 639		}
 640		#endregion
 641
 642		#region CreateScaleAndTranslation (Static)
 643		/// <summary>
 644		/// Create scale and translation, basically just does CreateScale and
 645		/// CreateTranslation, but much more efficient.
 646		/// </summary>
 647		/// <param name="scale">Scale for the new matrix</param>
 648		/// <param name="position">Position for the matrix translation.</param>
 649		/// <returns>New matrix with the given scale and translation.</returns>
 650		public static Matrix CreateScaleAndTranslation(float scale,
 651			Vector position)
 652		{
 653			return new Matrix(
 654				scale, 0f, 0f, 0f,
 655				0f, scale, 0f, 0f,
 656				0f, 0f, scale, 0f,
 657				position.X, position.Y, position.Z, 1f);
 658		}
 659		#endregion
 660
 661		#region CreateLookAt (Static)
 662		/// <summary>
 663		/// Creates a right-handed look at matrix for a camera.
 664		/// </summary>
 665		/// <param name="cameraPosition">camera position</param>
 666		/// <param name="cameraTarget">camera target</param>
 667		/// <param name="cameraUpVector">camera up vector</param>
 668		/// <returns>New look at matrix</returns>
 669		public static Matrix CreateLookAt(Vector cameraPosition,
 670			Vector cameraTarget, Vector cameraUpVector)
 671		{
 672			//Check performance difference again:
 673
 674			// Slightly faster in XNA (5%), but overall all these implementations
 675			// are roughly the same.
 676			// Note: All Create view implementations are slow, there are slighly
 677			// faster on XNA, but overall not much of a difference between all
 678			// platforms and frameworks (see below in MatrixPerformance class).
 679
 680
 681			// Simplified and optimized formula to create look matrix
 682			Vector dir = Vector.Normalize(cameraPosition - cameraTarget);
 683			Vector up = Vector.Normalize(Vector.Cross(cameraUpVector, dir));
 684			// Switch up around if cameraUpVector is the same as dir!
 685			if (up.LengthSquared == 0)
 686			{
 687				up = Vector.UnitY;
 688			}
 689			Vector right = Vector.Cross(dir, up);
 690			return new Matrix(
 691				up.X, right.X, dir.X, 0f,
 692				up.Y, right.Y, dir.Y, 0f,
 693				up.Z, right.Z, dir.Z, 0f,
 694				-Vector.Dot(up, cameraPosition),
 695				-Vector.Dot(right, cameraPosition),
 696				-Vector.Dot(dir, cameraPosition), 1f);
 697		}
 698		#endregion
 699
 700		#region CreateFromAxisAngle (Static)
 701		/// <summary>
 702		/// Create from axis angle (in degrees)
 703		/// </summary>
 704		/// <param name="angle">Angle in degrees</param>
 705		/// <param name="axis">Axis to rotate around</param>
 706		/// <returns>Rotated matrix around axis for rotations</returns>
 707		public static Matrix CreateFromAxisAngle(Vector axis, float angle)
 708		{
 709			// Note: Variables in this method need to be beautifulized!
 710			//Check performance difference again:
 711
 712			Matrix matrix = new Matrix();
 713
 714			float x = axis.X;
 715			float y = axis.Y;
 716			float z = axis.Z;
 717			float num2 = MathHelper.Sin(angle);
 718			float num = MathHelper.Cos(angle);
 719			float num11 = x * x;
 720			float num10 = y * y;
 721			float num9 = z * z;
 722			float num8 = x * y;
 723			float num7 = x * z;
 724			float num6 = y * z;
 725
 726			matrix.M11 = num11 + (num * (1f - num11));
 727			matrix.M12 = (num8 - (num * num8)) + (num2 * z);
 728			matrix.M13 = (num7 - (num * num7)) - (num2 * y);
 729			matrix.M14 = 0f;
 730			matrix.M21 = (num8 - (num * num8)) - (num2 * z);
 731			matrix.M22 = num10 + (num * (1f - num10));
 732			matrix.M23 = (num6 - (num * num6)) + (num2 * x);
 733			matrix.M24 = 0f;
 734			matrix.M31 = (num7 - (num * num7)) + (num2 * y);
 735			matrix.M32 = (num6 - (num * num6)) - (num2 * x);
 736			matrix.M33 = num9 + (num * (1f - num9));
 737			matrix.M34 = 0f;
 738			matrix.M41 = 0f;
 739			matrix.M42 = 0f;
 740			matrix.M43 = 0f;
 741			matrix.M44 = 1f;
 742
 743			return matrix;
 744		}
 745
 746		/// <summary>
 747		/// Create from axis angle (in degrees)
 748		/// </summary>
 749		/// <param name="angle">Angle in degrees</param>
 750		/// <param name="axis">Axis to rotate around</param>
 751		/// <param name="result">Rotated matrix around axis for rotations</param>
 752		public static void CreateFromAxisAngle(ref Vector axis, float angle,
 753			ref Matrix result)
 754		{
 755			float x = axis.X;
 756			float y = axis.Y;
 757			float z = axis.Z;
 758			float num2 = MathHelper.Sin(angle);
 759			float num = MathHelper.Cos(angle);
 760			float num11 = x * x;
 761			float num10 = y * y;
 762			float num9 = z * z;
 763			float num8 = x * y;
 764			float num7 = x * z;
 765			float num6 = y * z;
 766
 767			result.M11 = num11 + (num * (1f - num11));
 768			result.M12 = (num8 - (num * num8)) + (num2 * z);
 769			result.M13 = (num7 - (num * num7)) - (num2 * y);
 770			result.M14 = 0f;
 771			result.M21 = (num8 - (num * num8)) - (num2 * z);
 772			result.M22 = num10 + (num * (1f - num10));
 773			result.M23 = (num6 - (num * num6)) + (num2 * x);
 774			result.M24 = 0f;
 775			result.M31 = (num7 - (num * num7)) + (num2 * y);
 776			result.M32 = (num6 - (num * num6)) - (num2 * x);
 777			result.M33 = num9 + (num * (1f - num9));
 778			result.M34 = 0f;
 779			result.M41 = 0f;
 780			result.M42 = 0f;
 781			result.M43 = 0f;
 782			result.M44 = 1f;
 783		}
 784		#endregion
 785
 786		#region CreateOrthographic (Static)
 787		/// <summary>
 788		/// Create orthographic
 789		/// </summary>
 790		/// <param name="bottom">bottom</param>
 791		/// <param name="farPlane">farPlane</param>
 792		/// <param name="left">left</param>
 793		/// <param name="nearPlane">near plane</param>
 794		/// <param name="right">right</param>
 795		/// <param name="top">top</param>
 796		/// <returns>New Orthographic matrix</returns>
 797		public static Matrix CreateOrthographic(float left, float right,
 798			float bottom, float top, float nearPlane, float farPlane)
 799		{
 800			//Check performance difference again:
 801
 802			Matrix matrix = new Matrix();
 803
 804			matrix.M11 = 2f / (right - left);
 805			matrix.M12 = matrix.M13 = matrix.M14 = 0f;
 806			matrix.M22 = 2f / (top - bottom);
 807			matrix.M21 = matrix.M23 = matrix.M24 = 0f;
 808			matrix.M33 = 1f / (nearPlane - farPlane);
 809			matrix.M31 = matrix.M32 = matrix.M34 = 0f;
 810			matrix.M41 = (left + right) / (left - right);
 811			matrix.M42 = (top + bottom) / (bottom - top);
 812			matrix.M43 = nearPlane / (nearPlane - farPlane);
 813			matrix.M44 = 1f;
 814
 815			return matrix;
 816		}
 817		#endregion
 818
 819		#region CreatePerspective (Static)
 820		/// <summary>
 821		/// Creates a right-handed perspective field of view matrix for a camera.
 822		/// </summary>
 823		/// <param name="fieldOfView">Field of view (Warning: In radians)</param>
 824		/// <param name="aspectRatio">Aspect ratio</param>
 825		/// <param name="nearPlaneDistance">Near plane distance</param>
 826		/// <param name="farPlaneDistance">Far plane distance</param>
 827		/// <returns>Matrix</returns>
 828		public static Matrix CreatePerspective(float fieldOfView,
 829			float aspectRatio, float nearPlaneDistance, float farPlaneDistance)
 830		{
 831			//Check performance difference again:
 832			// Slightly faster in Delta (10%), but overall all this implementations
 833			// are ruffly the same. Same as the following, just faster: 
 834			//XnaMatrix.CreatePerspectiveFieldOfView(...)
 835			//TKMatrix.CreatePerspectiveFieldOfView(...)
 836			//SlimDX.Matrix.PerspectiveFovRH(...)
 837
 838			float m22 = 1f / ((float)Math.Tan(fieldOfView * 0.5f));
 839			return new Matrix(
 840				m22 / aspectRatio, 0f, 0f, 0f,
 841				0f, m22, 0f, 0f,
 842				0f, 0f, farPlaneDistance / (nearPlaneDistance - farPlaneDistance), -1f,
 843				0f, 0f, (nearPlaneDistance * farPlaneDistance) /
 844				        (nearPlaneDistance - farPlaneDistance), 0f);
 845		}
 846		#endregion
 847
 848		#region Multiply (Static)
 849		/// <summary>
 850		/// Ugly Matrix Multiply method (normally you would just use
 851		/// matrix1*matrix2 to multiply two matrices), but this one is a little
 852		/// faster because we don't need to copy over the matrix values and we can
 853		/// even provide a faster fallback on the xbox platform (using xna).
 854		/// Note: matrix1 or matrix2 can be the same as result, which results in
 855		/// even better performance for performance critical matrix multiply code.
 856		/// </summary>
 857		/// <param name="matrix1">matrix1</param>
 858		/// <param name="matrix2">matrix2</param>
 859		/// <param name="result">result</param>
 860		public static void Multiply(ref Matrix matrix1, ref Matrix matrix2,
 861			ref Matrix result)
 862		{
 863			//Check performance difference again:
 864			// Note: Values are stored in local variables to allow matrix1 or matrix2
 865			// and result be the same matrices (we won't overwrite M11-M44 while
 866			// calculating the result), which happens quite often in optimized code!
 867
 868			float m11 = matrix1.M11 * matrix2.M11 + matrix1.M12 * matrix2.M21 +
 869			            matrix1.M13 * matrix2.M31 + matrix1.M14 * matrix2.M41;
 870			float m12 = matrix1.M11 * matrix2.M12 + matrix1.M12 * matrix2.M22 +
 871			            matrix1.M13 * matrix2.M32 + matrix1.M14 * matrix2.M42;
 872			float m13 = matrix1.M11 * matrix2.M13 + matrix1.M12 * matrix2.M23 +
 873			            matrix1.M13 * matrix2.M33 + matrix1.M14 * matrix2.M43;
 874			float m14 = matrix1.M11 * matrix2.M14 + matrix1.M12 * matrix2.M24 +
 875			            matrix1.M13 * matrix2.M34 + matrix1.M14 * matrix2.M44;
 876			float m21 = matrix1.M21 * matrix2.M11 + matrix1.M22 * matrix2.M21 +
 877			            matrix1.M23 * matrix2.M31 + matrix1.M24 * matrix2.M41;
 878			float m22 = matrix1.M21 * matrix2.M12 + matrix1.M22 * matrix2.M22 +
 879			            matrix1.M23 * matrix2.M32 + matrix1.M24 * matrix2.M42;
 880			float m23 = matrix1.M21 * matrix2.M13 + matrix1.M22 * matrix2.M23 +
 881			            matrix1.M23 * matrix2.M33 + matrix1.M24 * matrix2.M43;
 882			float m24 = matrix1.M21 * matrix2.M14 + matrix1.M22 * matrix2.M24 +
 883			            matrix1.M23 * matrix2.M34 + matrix1.M24 * matrix2.M44;
 884			float m31 = matrix1.M31 * matrix2.M11 + matrix1.M32 * matrix2.M21 +
 885			            matrix1.M33 * matrix2.M31 + matrix1.M34 * matrix2.M41;
 886			float m32 = matrix1.M31 * matrix2.M12 + matrix1.M32 * matrix2.M22 +
 887			            matrix1.M33 * matrix2.M32 + matrix1.M34 * matrix2.M42;
 888			float m33 = matrix1.M31 * matrix2.M13 + matrix1.M32 * matrix2.M23 +
 889			            matrix1.M33 * matrix2.M33 + matrix1.M34 * matrix2.M43;
 890			float m34 = matrix1.M31 * matrix2.M14 + matrix1.M32 * matrix2.M24 +
 891			            matrix1.M33 * matrix2.M34 + matrix1.M34 * matrix2.M44;
 892			float m41 = matrix1.M41 * matrix2.M11 + matrix1.M42 * matrix2.M21 +
 893			            matrix1.M43 * matrix2.M31 + matrix1.M44 * matrix2.M41;
 894			float m42 = matrix1.M41 * matrix2.M12 + matrix1.M42 * matrix2.M22 +
 895			            matrix1.M43 * matrix2.M32 + matrix1.M44 * matrix2.M42;
 896			float m43 = matrix1.M41 * matrix2.M13 + matrix1.M42 * matrix2.M23 +
 897			            matrix1.M43 * matrix2.M33 + matrix1.M44 * matrix2.M43;
 898			float m44 = matrix1.M41 * matrix2.M14 + matrix1.M42 * matrix2.M24 +
 899			            matrix1.M43 * matrix2.M34 + matrix1.M44 * matrix2.M44;
 900			result.M11 = m11;
 901			result.M12 = m12;
 902			result.M13 = m13;
 903			result.M14 = m14;
 904			result.M21 = m21;
 905			result.M22 = m22;
 906			result.M23 = m23;
 907			result.M24 = m24;
 908			result.M31 = m31;
 909			result.M32 = m32;
 910			result.M33 = m33;
 911			result.M34 = m34;
 912			result.M41 = m41;
 913			result.M42 = m42;
 914			result.M43 = m43;
 915			result.M44 = m44;
 916		}
 917
 918		/// <summary>
 919		/// Ugly Matrix*Vector multiply method (basically transforming a vector
 920		/// with a matrix). Normally you would just use matrix*position to
 921		/// multiply, but this one is a little faster because we don't need to
 922		/// copy over the matrix values and we can even provide a faster fallback
 923		/// on the xbox platform (using xna).
 924		/// Note: position and result can be the same as result, which results in
 925		/// even better performance for performance critical matrix multiply code.
 926		/// </summary>
 927		/// <param name="matrix">matrix</param>
 928		/// <param name="position">position</param>
 929		/// <param name="result">result</param>
 930		public static void Multiply(ref Matrix matrix, ref Vector position,
 931			ref Vector result)
 932		{
 933			//Check performance difference again:
 934
 935			float x = position.X * matrix.M11 + position.Y * matrix.M21 +
 936			          position.Z * matrix.M31 + matrix.M41;
 937			float y = position.X * matrix.M12 + position.Y * matrix.M22 +
 938			          position.Z * matrix.M32 + matrix.M42;
 939			float z = position.X * matrix.M13 + position.Y * matrix.M23 +
 940			          position.Z * matrix.M33 + matrix.M43;
 941			result.X = x;
 942			result.Y = y;
 943			result.Z = z;
 944		}
 945
 946		/// <summary>
 947		/// Multiply matrix with a scale factor. Note: This is the slow version.
 948		/// This method needs to copy 128 bytes around (2 matrices, each 64 bytes).
 949		/// Try to use the ref version instead.
 950		/// </summary>
 951		/// <param name="matrix">Matrix</param>
 952		/// <param name="scaleFactor">Scale factor</param>
 953		/// <returns>
 954		/// New matrix with each value multiplied with scaleFactor
 955		/// </returns>
 956		public static Matrix Multiply(Matrix matrix, float scaleFactor)
 957		{
 958			return new Matrix(
 959				matrix.M11 * scaleFactor, matrix.M12 * scaleFactor,
 960				matrix.M13 * scaleFactor, matrix.M14 * scaleFactor,
 961				matrix.M21 * scaleFactor, matrix.M22 * scaleFactor,
 962				matrix.M23 * scaleFactor, matrix.M24 * scaleFactor,
 963				matrix.M31 * scaleFactor, matrix.M32 * scaleFactor,
 964				matrix.M33 * scaleFactor, matrix.M34 * scaleFactor,
 965				matrix.M41 * scaleFactor, matrix.M42 * scaleFactor,
 966				matrix.M43 * scaleFactor, matrix.M44 * scaleFactor);
 967		}
 968		#endregion
 969
 970		#region Equal (Static)
 971		/// <summary>
 972		/// Equality check. Optimized for speed via ref parameters, which
 973		/// are much faster than copying 128 bytes (2 matrices each 64 bytes).
 974		/// </summary>
 975		public static bool Equal(ref Matrix matrix1, ref Matrix matrix2)
 976		{
 977			// Note: The checks are sorted by the diagonal, which is most important
 978			// and usually the most different for 2 matrices (faster this way!)
 979			// Note: Calling XNA (and OpenTK or SlimDX or SharpDX too) here is
 980			// slower! So it is disabled.
 981			return
 982				matrix1.M11 == matrix2.M11 &&
 983				matrix1.M22 == matrix2.M22 &&
 984				matrix1.M33 == matrix2.M33 &&
 985				matrix1.M44 == matrix2.M44 &&
 986				// Now check the rest if the diagonal matrices were the same.
 987				matrix1.M12 == matrix2.M12 &&
 988				matrix1.M13 == matrix2.M13 &&
 989				matrix1.M14 == matrix2.M14 &&
 990				matrix1.M21 == matrix2.M21 &&
 991				matrix1.M23 == matrix2.M23 &&
 992				matrix1.M24 == matrix2.M24 &&
 993				matrix1.M31 == matrix2.M31 &&
 994				matrix1.M32 == matrix2.M32 &&
 995				matrix1.M34 == matrix2.M34 &&
 996				matrix1.M41 == matrix2.M41 &&
 997				matrix1.M42 == matrix2.M42 &&
 998				matrix1.M43 == matrix2.M43;
 999		}
1000		#endregion
1001
1002		#region Unequal (Static)
1003		/// <summary>
1004		/// Inequality check. Optimized for speed via ref parameters, which
1005		/// are much faster than copying 128 bytes (2 matrices each 64 bytes).
1006		/// </summary>
1007		public static bool Unequal(ref Matrix matrix1, ref Matrix matrix2)
1008		{
1009			// Check if everything is equal, if not return true (this is the unequal
1010			// check), might look a bit strange, but it is fast because of early out.
1011			if (matrix1.M11 == matrix2.M11 &&
1012			    matrix1.M12 == matrix2.M12 &&
1013			    matrix1.M13 == matrix2.M13 &&
1014			    matrix1.M14 == matrix2.M14 &&
1015			    matrix1.M21 == matrix2.M21 &&
1016			    matrix1.M22 == matrix2.M22 &&
1017			    matrix1.M23 == matrix2.M23 &&
1018			    matrix1.M24 == matrix2.M24 &&
1019			    matrix1.M31 == matrix2.M31 &&
1020			    matrix1.M32 == matrix2.M32 &&
1021			    matrix1.M33 == matrix2.M33 &&
1022			    matrix1.M34 == matrix2.M34 &&
1023			    matrix1.M41 == matrix2.M41 &&
1024			    matrix1.M42 == matrix2.M42 &&
1025			    matrix1.M43 == matrix2.M43)
1026			{
1027				return (matrix1.M44 != matrix2.M44);
1028			}
1029
1030			return true;
1031		}
1032		#endregion
1033
1034		#region TranslationNearlyEqual (Static)
1035		/// <summary>
1036		/// Nearly equals (only checks the Translation of the matrices), two
1037		/// matrices are considered nearly equal if their translation distance
1038		/// is below MathHelper.Epsilon. Use the MatrixNearlyEqual method to
1039		/// compare two complete matrices (much slower than this check).
1040		/// </summary>
1041		/// <param name="matrix1">Matrix 1 to compare</param>
1042		/// <param name="matrix2">Matrix 2 to compare</param>
1043		/// <returns>True if the translations of the matrices are almost the
1044		/// same, false otherwise.</returns>
1045		public static bool TranslationNearlyEqual(ref Matrix matrix1,
1046			ref Matrix matrix2)
1047		{
1048			// Check to see if the absolute difference of each value is smaller than
1049			// the Nearly Equal value.
1050			return (matrix1.Translation - matrix2.Translation).LengthSquared <
1051			       MathHelper.Epsilon * MathHelper.Epsilon;
1052		}
1053		#endregion
1054
1055		#region MatrixNearlyEqual (Static)
1056		/// <summary>
1057		/// Matrix nearly equal helper method to compare two matrices that might
1058		/// not be 100% the same because of rounding errors, but they produce
1059		/// pretty much the same results.
1060		/// </summary>
1061		/// <param name="matrix1">Matrix 1 to compare</param>
1062		/// <param name="matrix2">Matrix 2 to compare</param>
1063		/// <returns>True if the matrices are almost the same</returns>
1064		public static bool MatrixNearlyEqual(ref Matrix mat1, ref Matrix mat2)
1065		{
1066			return
1067				mat1.M11.NearlyEqual(mat2.M11) &&
1068				mat1.M12.NearlyEqual(mat2.M12) &&
1069				mat1.M13.NearlyEqual(mat2.M13) &&
1070				mat1.M14.NearlyEqual(mat2.M14) &&
1071				mat1.M21.NearlyEqual(mat2.M21) &&
1072				mat1.M22.NearlyEqual(mat2.M22) &&
1073				mat1.M23.NearlyEqual(mat2.M23) &&
1074				mat1.M24.NearlyEqual(mat2.M24) &&
1075				mat1.M31.NearlyEqual(mat2.M31) &&
1076				mat1.M32.NearlyEqual(mat2.M32) &&
1077				mat1.M33.NearlyEqual(mat2.M33) &&
1078				mat1.M34.NearlyEqual(mat2.M34) &&
1079				mat1.M41.NearlyEqual(mat2.M41) &&
1080				mat1.M42.NearlyEqual(mat2.M42) &&
1081				mat1.M43.NearlyEqual(mat2.M43) &&
1082				mat1.M44.NearlyEqual(mat2.M44);
1083		}
1084		#endregion
1085
1086		#region FromQuaternion (Static)
1087		/// <summary>
1088		/// From quaternion
1089		/// </summary>
1090		public static Matrix FromQuaternion(Quaternion quaternion)
1091		{
1092			//Check performance difference again:
1093			// Note: We can only be sure that this is faster on the Xbox, not used
1094			// often anyway.
1095
1096			//Note: This could need some more comments, help, etc.
1097			float qxx = quaternion.X * quaternion.X;
1098			float qyy = quaternion.Y * quaternion.Y;
1099			float qzz = quaternion.Z * quaternion.Z;
1100			float qxy = quaternion.X * quaternion.Y;
1101			float qzw = quaternion.Z * quaternion.W;
1102			float qwx = quaternion.Z * quaternion.X;
1103			float qyw = quaternion.Y * quaternion.W;
1104			float qyz = quaternion.Y * quaternion.Z;
1105			float qxw = quaternion.X * quaternion.W;
1106
1107			return new Matrix(
1108				1f - (2f * (qyy + qzz)), 2f * (qxy + qzw), 2f * (qwx - qyw), 0f,
1109				2f * (qxy - qzw), 1f - (2f * (qzz + qxx)), 2f * (qyz + qxw), 0f,
1110				2f * (qwx + qyw), 2f * (qyz - qxw), 1f - (2f * (qyy + qxx)), 0f,
1111				0f, 0f, 0f, 1f);
1112		}
1113		#endregion
1114
1115		#region FromString (Static)
1116		/// <summary>
1117		/// Convert a string to a Matrix. The expected format is
1118		/// (M11, M12, M13, M14)
1119		/// (M21, M22, M23, M24)
1120		/// ...
1121		/// </summary>
1122		/// <param name="matrixString">The string containing the values in the
1123		/// correct format.</param>
1124		public static Matrix FromString(string matrixString)
1125		{
1126			// First remove the brackets
1127			matrixString = matrixString.
1128				Replace("(", "").
1129				Replace(")", "");
1130			// And split the string up into seperate values.
1131			string[] matrixStrings = matrixString.SplitAndTrim(',');
1132			//StringSplitOptions.RemoveEmptyEntries);
1133
1134			// Then check if the length is 16 for 16 values and return the new matrix.
1135			// If the length is not 16 than return Matrix.Identity.
1136			if (matrixStrings.Length == 16)
1137			{
1138				return new Matrix(
1139					// First row
1140					matrixStrings[0].FromInvariantString(1.0f),
1141					matrixStrings[1].FromInvariantString(0.0f),
1142					matrixStrings[2].FromInvariantString(0.0f),
1143					matrixStrings[3].FromInvariantString(0.0f),
1144					// Second row
1145					matrixStrings[4].FromInvariantString(0.0f),
1146					matrixStrings[5].FromInvariantString(1.0f),
1147					matrixStrings[6].FromInvariantString(0.0f),
1148					matrixStrings[7].FromInvariantString(0.0f),
1149					// Third row
1150					matrixStrings[8].FromInvariantString(0.0f),
1151					matrixStrings[9].FromInvariantString(0.0f),
1152					matrixStrings[10].FromInvariantString(1.0f),
1153					matrixStrings[11].FromInvariantString(0.0f),
1154					// Forth row
1155					matrixStrings[12].FromInvariantString(0.0f),
1156					matrixStrings[13].FromInvariantString(0.0f),
1157					matrixStrings[14].FromInvariantString(0.0f),
1158					matrixStrings[15].FromInvariantString(1.0f));
1159			}
1160
1161			Log.Warning("Unable to convert matrixString=" + matrixString +
1162			            " because it has not exactly 16 float values!");
1163			return Identity;
1164		}
1165		#endregion
1166
1167		#region Identity (Static)
1168		/// <summary>
1169		/// Returns a new identity matrix (no rotation, translation or scaling)!
1170		/// Note: Not readonly to allow passing it by ref, but this should never
1171		/// be set!
1172		/// </summary>
1173		public static Matrix Identity = new Matrix(
1174			1, 0, 0, 0,
1175			0, 1, 0, 0,
1176			0, 0, 1, 0,
1177			0, 0, 0, 1);
1178		#endregion
1179
1180		#region Framework Union Defines (Public)
1181		#endregion
1182
1183		#region M11 (Public)
1184		/// <summary>
1185		/// 1st row - 1st column value of the Matrix.
1186		/// </summary>
1187		[FieldOffset(0)]
1188		public float M11;
1189		#endregion
1190
1191		#region M12 (Public)
1192		/// <summary>
1193		/// 1st row - 2nd column value of the Matrix.
1194		/// </summary>
1195		[FieldOffset(4)]
1196		public float M12;
1197		#endregion
1198
1199		#region M13 (Public)
1200		/// <summary>
1201		/// 1st row - 3rd column value of the Matrix.
1202		/// </summary>
1203		[FieldOffset(8)]
1204		public float M13;
1205		#endregion
1206
1207		#region M14 (Public)
1208		/// <summary>
1209		/// 1st row - 4th column value of the Matrix.
1210		/// </summary>
1211		[FieldOffset(12)]
1212		public float M14;
1213		#endregion
1214
1215		#region M21 (Public)
1216		/// <summary>
1217		/// 2nd row - 1st column value of the Matrix.
1218		/// </summary>
1219		[FieldOffset(16)]
1220		public float M21;
1221		#endregion
1222
1223		#region M22 (Public)
1224		/// <summary>
1225		/// 2nd row - 2nd column value of the Matrix.
1226		/// </summary>
1227		[FieldOffset(20)]
1228		public float M22;
1229		#endregion
1230
1231		#region M23 (Public)
1232		/// <summary>
1233		/// 2nd row - 3nd column value of the Matrix.
1234		/// </summary>
1235		[FieldOffset(24)]
1236		public float M23;
1237		#endregion
1238
1239		#region M24 (Public)
1240		/// <summary>
1241		/// 2nd row - 4nd column value of the Matrix.
1242		/// </summary>
1243		[FieldOffset(28)]
1244		public float M24;
1245		#endregion
1246
1247		#region M31 (Public)
1248		/// <summary>
1249		/// 3rd row - 1st column value of the Matrix.
1250		/// </summary>
1251		[FieldOffset(32)]
1252		public float M31;
1253		#endregion
1254
1255		#region M32 (Public)
1256		/// <summary>
1257		/// 3rd row - 2nd column value of the Matrix.
1258		/// </summary>
1259		[FieldOffset(36)]
1260		public float M32;
1261		#endregion
1262
1263		#region M33 (Public)
1264		/// <summary>
1265		/// 3rd row - 3rd column value of the Matrix.
1266		/// </summary>
1267		[FieldOffset(40)]
1268		public float M33;
1269		#endregion
1270
1271		#region M34 (Public)
1272		/// <summary>
1273		/// 3rd row - 4th column value of the Matrix.
1274		/// </summary>
1275		[FieldOffset(44)]
1276		public float M34;
1277		#endregion
1278
1279		#region M41 (Public)
1280		/// <summary>
1281		/// 4th row - 1st column value of the Matrix.
1282		/// </summary>
1283		[FieldOffset(48)]
1284		public float M41;
1285		#endregion
1286
1287		#region M42 (Public)
1288		/// <summary>
1289		/// 4th row - 2nd column value of the Matrix.
1290		/// </summary>
1291		[FieldOffset(52)]
1292		public float M42;
1293		#endregion
1294
1295		#region M43 (Public)
1296		/// <summary>
1297		/// 4th row - 3rd column value of the Matrix.
1298		/// </summary>
1299		[FieldOffset(56)]
1300		public float M43;
1301		#endregion
1302
1303		#region M44 (Public)
1304		/// <summary>
1305		/// 4th row - 4th column value of the Matrix.
1306		/// </summary>
1307		[FieldOffset(60)]
1308		public float M44;
1309		#endregion
1310
1311		#region Translation (Public)
1312		/// <summary>
1313		/// The translation amount representing by the matrix.
1314		/// This vector shares the same data as M41, M42, M43, see above!
1315		/// </summary>
1316		[FieldOffset(48)]
1317		public Vector Translation;
1318		#endregion
1319
1320		#region Scaling (Public)
1321		/// <summary>
1322		/// The scaling amount representing by the matrix.
1323		/// </summary>
1324		public Vector Scaling
1325		{
1326			get
1327			{
1328				return new Vector(M11, M22, M33);
1329			}
1330			set
1331			{
1332				M11 = value.X;
1333				M22 = value.Y;
1334				M33 = value.Z;
1335			}
1336		}
1337		#endregion
1338
1339		#region Right (Public)
1340		/// <summary>
1341		/// Right vector of the matrix (M11, M12, M13)
1342		/// </summary>
1343		[FieldOffset(0)]
1344		public Vector Right;
1345		#endregion
1346
1347		#region Up (Public)
1348		/// <summary>
1349		/// Up vector of the matrix (M21, M22, M23)
1350		/// </summary>
1351		[FieldOffset(16)]
1352		public Vector Up;
1353		#endregion
1354
1355		#region Front (Public)
1356		/// <summary>
1357		/// Front vector of the matrix (M31, M32, M33)
1358		/// </summary>
1359		[FieldOffset(32)]
1360		public Vector Front;
1361		#endregion
1362
1363		#region Determinant (Public)
1364		/// <summary>
1365		/// Returns the determinant of the matrix. Uses optimized code paths
1366		/// because XNAs implementation is twice as fast as ours, OpenTK is
1367		/// slightly faster (see MatrixPerformance class).
1368		/// </summary>
1369		public float Determinant
1370		{
1371			get
1372			{
1373				//Check performance difference again:
1374
1375				// Ugly hacky code that is actually lots faster than our clean
1376				// implementation below (see link on how it works basically).
1377				float m44x33m43x34 = (M44 * M33) - (M43 * M34);
1378				float m32x44m42x34 = (M32 * M44) - (M42 * M34);
1379				float m32x43m42x33 = (M32 * M43) - (M42 * M33);
1380				float m31x44m41x34 = (M31 * M44) - (M41 * M34);
1381				float m31x43m41x33 = (M31 * M43) - (M41 * M33);
1382				float m31x42m41x32 = (M31 * M42) - (M41 * M32);
1383				return (((((((M22 * m44x33m43x34) - (M23 * m32x44m42x34)) +
1384				            (M24 * m32x43m42x33)) * M11) - ((((M21 * m44x33m43x34) -
1385				                                              (M23 * m31x44m41x34)) +
1386				                                             (M24 * m31x43m41x33)) * M12)) +
1387				         ((((M21 * m32x44m42x34) - (M22 * m31x44m41x34)) +
1388				           (M24 * m31x42m41x32)) * M13)) - ((((M21 * m32x43m42x33) -
1389				                                              (M22 * m31x43m41x33)) +
1390				                                             (M23 * m31x42m41x32)) * M14));
1391			}
1392		}
1393		#endregion
1394
1395		#region Inverse (Public)
1396		/// <summary>
1397		/// Returns the inverse of the current matrix, which can be slightly
1398		/// faster than using the static Invert method below depending on the
1399		/// platform (e.g. with XNA), with Delta matrix code Invert() is a little
1400		/// faster however, so just use whatever method fits best for you.
1401		/// Both methods are the slowest matrix methods however, use with care!
1402		/// </summary>
1403		public Matrix Inverse
1404		{
1405			get
1406			{
1407				//Check performance difference again:
1408				// Ugly hacky code that is actually lots faster than our clean
1409				// implementation below (see link on how it works basically).
1410				// This is basically getting the Determinant and then also calculating
1411				// the inverse matrix in one big ugly swoop.
1412				float m33x44m34x43 = (M33 * M44) - (M34 * M43);
1413				float m32x44m34x42 = (M32 * M44) - (M34 * M42);
1414				float m32x43m33x42 = (M32 * M43) - (M33 * M42);
1415				float m31x44m34x41 = (M31 * M44) - (M34 * M41);
1416				float m31x43m33x41 = (M31 * M43) - (M33 * M41);
1417				float m31x42m32x41 = (M31 * M42) - (M32 * M41);
1418				float m22m23m24 = ((M22 * m33x44m34x43) - (M23 * m32x44m34x42)) +
1419				                  (M24 * m32x43m33x42);
1420				float m21m23m24 = -(((M21 * m33x44m34x43) - (M23 * m31x44m34x41)) +
1421				                    (M24 * m31x43m33x41));
1422				float m21m22m24 = ((M21 * m32x44m34x42) - (M22 * m31x44m34x41)) +
1423				                  (M24 * m31x42m32x41);
1424				float m21m22m23 = -(((M21 * m32x43m33x42) - (M22 * m31x43m33x41)) +
1425				                    (M23 * m31x42m32x41));
1426				// Inverse(A) = 1/det(A) * B
1427				float inverseDet = 1f / ((((M11 * m22m23m24) + (M12 * m21m23m24)) +
1428				                          (M13 * m21m22m24)) + (M14 * m21m22m23));
1429				float m23x44m24x43 = (M23 * M44) - (M24 * M43);
1430				float m22x44m24x42 = (M22 * M44) - (M24 * M42);
1431				float m22x43m23x42 = (M22 * M43) - (M23 * M42);
1432				float m21x44m24x41 = (M21 * M44) - (M24 * M41);
1433				float m21x43m23x41 = (M21 * M43) - (M23 * M41);
1434				float m21x42m22x41 = (M21 * M42) - (M22 * M41);
1435				float m23x34m24x33 = (M23 * M34) - (M24 * M33);
1436				float m22x34m24x32 = (M22 * M34) - (M24 * M32);
1437				float m22x33m23x32 = (M22 * M33) - (M23 * M32);
1438				float m21x34m23x31 = (M21 * M34) - (M24 * M31);
1439				float m21x33m23x31 = (M21 * M33) - (M23 * M31);
1440				float m21x32m22x31 = (M21 * M32) - (M22 * M31);
1441				return new Matrix(
1442					m22m23m24 * inverseDet,
1443					-(((M12 * m33x44m34x43) - (M13 * m32x44m34x42)) +
1444					  (M14 * m32x43m33x42)) * inverseDet,
1445					(((M12 * m23x44m24x43) - (M13 * m22x44m24x42)) +
1446					 (M14 * m22x43m23x42)) * inverseDet,
1447					-(((M12 * m23x34m24x33) - (M13 * m22x34m24x32)) +
1448					  (M14 * m22x33m23x32)) * inverseDet,
1449					m21m23m24 * inverseDet,
1450					(((M11 * m33x44m34x43) - (M13 * m31x44m34x41)) +
1451					 (M14 * m31x43m33x41)) * inverseDet,
1452					-(((M11 * m23x44m24x43) - (M13 * m21x44m24x41)) +
1453					  (M14 * m21x43m23x41)) * inverseDet,
1454					(((M11 * m23x34m24x33) - (M13 * m21x34m23x31)) +
1455					 (M14 * m21x33m23x31)) * inverseDet,
1456					m21m22m24 * inverseDet,
1457					-(((M11 * m32x44m34x42) - (M12 * m31x44m34x41)) +
1458					  (M14 * m31x42m32x41)) * inverseDet,
1459					(((M11 * m22x44m24x42) - (M12 * m21x44m24x41)) +
1460					 (M14 * m21x42m22x41)) * inverseDet,
1461					-(((M11 * m22x34m24x32) - (M12 * m21x34m23x31)) +
1462					  (M14 * m21x32m22x31)) * inverseDet,
1463					m21m22m23 * inverseDet,
1464					(((M11 * m32x43m33x42) - (M12 * m31x43m33x41)) +
1465					 (M13 * m31x42m32x41)) * inverseDet,
1466					-(((M11 * m22x43m23x42) - (M12 * m21x43m23x41)) +
1467					  (M13 * m21x42m22x41)) * inverseDet,
1468					(((M11 * m22x33m23x32) - (M12 * m21x33m23x31)) +
1469					 (M13 * m21x32m22x31)) * inverseDet);
1470			}
1471		}
1472		#endregion
1473
1474		#region Constructors
1475		/// <summary>
1476		/// Create matrix by passing all the necessary values
1477		/// </summary>
1478		/// <param name="setM11">Set M11 value (first row)</param>
1479		/// <param name="setM12">Set M12 value (first row)</param>
1480		/// <param name="setM13">Set M13 value (first row)</param>
1481		/// <param name="setM14">Set M14 value (first row)</param>
1482		/// <param name="setM21">Set M21 value (second row)</param>
1483		/// <param name="setM22">Set M22 value (second row)</param>
1484		/// <param name="setM23">Set M23 value (second row)</param>
1485		/// <param name="setM24">Set M24 value (second row)</param>
1486		/// <param name="setM31">Set M31 value (third row)</param>
1487		/// <param name="setM32">Set M32 value (third row)</param>
1488		/// <param name="setM33">Set M33 value (third row)</param>
1489		/// <param name="setM34">Set M34 value (third row)</param>
1490		/// <param name="setM41">Set M41 value (forth row)</param>
1491		/// <param name="setM42">Set M42 value (forth row)</param>
1492		/// <param name="setM43">Set M43 value (forth row)</param>
1493		/// <param name="setM44">Set M44 value (forth row)</param>
1494		public Matrix(float setM11, float setM12, float setM13, float setM14,
1495			float setM21, float setM22, float setM23, float setM24,
1496			float setM31, float setM32, float setM33, float setM34,
1497			float setM41, float setM42, float setM43, float setM44)
1498			: this()
1499		{
1500			M11 = setM11;
1501			M12 = setM12;
1502			M13 = setM13;
1503			M14 = setM14;
1504			M21 = setM21;
1505			M22 = setM22;
1506			M23 = setM23;
1507			M24 = setM24;
1508			M31 = setM31;
1509			M32 = setM32;
1510			M33 = setM33;
1511			M34 = setM34;
1512			M41 = setM41;
1513			M42 = setM42;
1514			M43 = setM43;
1515			M44 = setM44;
1516		}
1517
1518		/// <summary>
1519		/// Create matrix by passing all 16 matrix values as array.
1520		/// m11, m12, m13, m14, 
1521		/// m21, m22, m23, m24,
1522		/// m31, m32, m33, m34, 
1523		/// m41, m42, m43, m44, 
1524		/// </summary>
1525		/// <param name="setMatrix">Set Matrix</param>
1526		public Matrix(float[] setMatrix)
1527			: this()
1528		{
1529			int length = setMatrix.Length;
1530			if (length < 16)
1531			{
1532				Log.Warning("The matrix has too little values.\n" +
1533				            "Length: " + length);
1534			}
1535			if (length > 16)
1536			{
1537				Log.Warning("The matrix has too many values.\n" +
1538				            "Length: " + length);
1539			}
1540
1541			// Copy all values.
1542			M11 = setMatrix[0];
1543			M12 = setMatrix[1];
1544			M13 = setMatrix[2];
1545			M14 = setMatrix[3];
1546			M21 = setMatrix[4];
1547			M22 = setMatrix[5];
1548			M23 = setMatrix[6];
1549			M24 = setMatrix[7];
1550			M31 = setMatrix[8];
1551			M32 = setMatrix[9];
1552			M33 = setMatrix[10];
1553			M34 = setMatrix[11];
1554			M41 = setMatrix[12];
1555			M42 = setMatrix[13];
1556			M43 = setMatrix[14];
1557			M44 = setMatrix[15];
1558		}
1559
1560		/// <summary>
1561		/// Create matrix by passing in 3 vectors for M11-M33, rest stays default.
1562		/// </summary>
1563		/// <param name="firstRowVector">Vector for the first row (M11-M13)
1564		/// </param>
1565		/// <param name="secondRowVector">Vector for the second row (M21-M23)
1566		/// </param>
1567		/// <param name="thirdRowVector">Vector for the third row (M31-M33)
1568		/// </param>
1569		public Matrix(Vector firstRowVector, Vector secondRowVector,
1570			Vector thirdRowVector)
1571			: this()
1572		{
1573			// Copy all values and fill in the blanks
1574			M11 = firstRowVector.X;
1575			M12 = firstRowVector.Y;
1576			M13 = firstRowVector.Z;
1577			M14 = 0.0f;
1578			M21 = secondRowVector.X;
1579			M22 = secondRowVector.Y;
1580			M23 = secondRowVector.Z;
1581			M24 = 0.0f;
1582			M31 = thirdRowVector.X;
1583			M32 = thirdRowVector.Y;
1584			M33 = thirdRowVector.Z;
1585			M34 = 0.0f;
1586			M41 = 0.0f;
1587			M42 = 0.0f;
1588			M43 = 0.0f;
1589			M44 = 1.0f;
1590		}
1591
1592		/// <summary>
1593		/// Create matrix with help of BinaryReader, we will read 16 floats from
1594		/// the stream that is linked up with this BinaryReader.
1595		/// </summary>
1596		/// <param name="reader">Reader</param>
1597		public Matrix(BinaryReader reader)
1598			: this()
1599		{
1600			Load(reader);
1601		}
1602		#endregion
1603
1604		#region IEquatable<Matrix> Members
1605		/// <summary>
1606		/// Equals, quickly checks if another matrix has the exact same values.
1607		/// This methods needs to copy 64 bytes around (1 matrix of 64 bytes).
1608		/// Try to use the ref version instead.
1609		/// </summary>
1610		/// <param name="other">Other Matrix to compare to</param>
1611		/// <returns>True if both matrices are equal, false otherwise.</returns>
1612		public bool Equals(Matrix other)
1613		{
1614			// First check the diagonal values, those mostly differ.
1615			return
1616				M11 == other.M11 &&
1617				M22 == other.M22 &&
1618				M33 == other.M33 &&
1619				M44 == other.M44 &&
1620				// Now check the rest if the diagonal values were the same.
1621				M12 == other.M12 &&
1622				M13 == other.M13 &&
1623				M14 == other.M14 &&
1624				M21 == other.M21 &&
1625				M23 == other.M23 &&
1626				M24 == other.M24 &&
1627				M31 == other.M31 &&
1628				M32 == other.M32 &&
1629				M34 == other.M34 &&
1630				M41 == other.M41 &&
1631				M42 == other.M42 &&
1632				M43 == other.M43;
1633		}
1634		#endregion
1635
1636		#region ISaveLoadBinary Members
1637		/// <summary>
1638		/// Load the matrix values from a stream (64 bytes, 16 floats).
1639		/// </summary>
1640		/// <param name="reader">The stream that will be used.</param>
1641		public void Load(BinaryReader reader)
1642		{
1643			// Note XNA and OpenTK matrices share the same data, no need to set them
1644			M11 = reader.ReadSingle();
1645			M12 = reader.ReadSingle();
1646			M13 = reader.ReadSingle();
1647			M14 = reader.ReadSingle();
1648
1649			M21 = reader.ReadSingle();
1650			M22 = reader.ReadSingle();
1651			M23 = reader.ReadSingle();
1652			M24 = reader.ReadSingle();
1653
1654			M31 = reader.ReadSingle();
1655			M32 = reader.ReadSingle();
1656			M33 = reader.ReadSingle();
1657			M34 = reader.ReadSingle();
1658
1659			M41 = reader.ReadSingle();
1660			M42 = reader.ReadSingle();
1661			M43 = reader.ReadSingle();
1662			M44 = reader.ReadSingle();
1663		}
1664
1665		/// <summary>
1666		/// Saves the matrix to a stream (64 bytes, 16 floats).
1667		/// </summary>
1668		/// <param name="writer">The stream that will be used.</param>
1669		public void Save(BinaryWriter writer)
1670		{
1671			writer.Write(M11);
1672			writer.Write(M12);
1673			writer.Write(M13);
1674			writer.Write(M14);
1675
1676			writer.Write(M21);
1677			writer.Write(M22);
1678			writer.Write(M23);
1679			writer.Write(M24);
1680
1681			writer.Write(M31);
1682			writer.Write(M32);
1683			writer.Write(M33);
1684			writer.Write(M34);
1685
1686			writer.Write(M41);
1687			writer.Write(M42);
1688			writer.Write(M43);
1689			writer.Write(M44);
1690		}
1691		#endregion
1692
1693		#region op_Multiply (Operator)
1694		/// <summary>
1695		/// Multiply operator. Warning: This operator is slower than using ref
1696		/// version, it needs to copy 192 bytes around (3 matrices, each 64 bytes).
1697		/// </summary>
1698		/// <param name="matrix1">Matrix 1 to multiply</param>
1699		/// <param name="matrix2">Matrix 2 to multiply</param>
1700		public static Matrix operator *(Matrix matrix1, Matrix matrix2)
1701		{
1702			// Other discussions about XNA Multiply performance
1703			// http://forums.xna.com/forums/t/7061.aspx?PageIndex=2
1704			// Machine: Core 2 Duo E6600 clocked at 3.0GHz
1705			// Managed XNA - Matrix.Multiply(ref,ref,out): 557ms
1706			// Managed w/ P/Invoke to SSE2 function:       341ms
1707			// Native C++ - Full Compiler Optimizations:   300ms
1708			// Native C++ - hand-optimized SSE2:           138ms
1709			// Note: We got like 97ms with this C# code, but a faster PC ^^
1710
1711			// Try to do some early out code in case matrix2 is an identity matrix
1712			if (matrix2.M11 == 1.0f &&
1713			    matrix2.M22 == 1.0f &&
1714			    matrix2.M33 == 1.0f &&
1715			    matrix2.M44 == 1.0f &&
1716			    // Ok, just to be sure, we need to check if all the other values
1717			    // are really 0, only then this is an identity matrix
1718			    matrix2.M21 == 0.0f &&
1719			    matrix2.M23 == 0.0f &&
1720			    matrix2.M24 == 0.0f &&
1721			    matrix2.M31 == 0.0f &&
1722			    matrix2.M32 == 0.0f &&
1723			    matrix2.M34 == 0.0f &&
1724			    matrix2.M41 == 0.0f &&
1725			    matrix2.M42 == 0.0f &&
1726			    matrix2.M43 == 0.0f)
1727			{
1728				return matrix1;
1729			}
1730
1731			return new Matrix(
1732				// M11
1733				(matrix1.M11 * matrix2.M11) + (matrix1.M12 * matrix2.M21) +
1734				(matrix1.M13 * matrix2.M31) + (matrix1.M14 * matrix2.M41),
1735				// M12
1736				(matrix1.M11 * matrix2.M12) + (matrix1.M12 * matrix2.M22) +
1737				(matrix1.M13 * matrix2.M32) + (matrix1.M14 * matrix2.M42),
1738				// M13
1739				(matrix1.M11 * matrix2.M13) + (matrix1.M12 * matrix2.M23) +
1740				(matrix1.M13 * matrix2.M33) + (matrix1.M14 * matrix2.M43),
1741				// M14
1742				(matrix1.M11 * matrix2.M14) + (matrix1.M12 * matrix2.M24) +
1743				(matrix1.M13 * matrix2.M34) + (matrix1.M14 * matrix2.M44),
1744				// M21
1745				(matrix1.M21 * matrix2.M11) + (matrix1.M22 * matrix2.M21) +
1746				(matrix1.M23 * matrix2.M31) + (matrix1.M24 * matrix2.M41),
1747				// M22
1748				(matrix1.M21 * matrix2.M12) + (matrix1.M22 * matrix2.M22) +
1749				(matrix1.M23 * matrix2.M32) + (matrix1.M24 * matrix2.M42),
1750				// M23
1751				(matrix1.M21 * matrix2.M13) + (matrix1.M22 * matrix2.M23) +
1752				(matrix1.M23 * matrix2.M33) + (matrix1.M24 * matrix2.M43),
1753				// M24
1754				(matrix1.M21 * matrix2.M14) + (matrix1.M22 * matrix2.M24) +
1755				(matrix1.M23 * matrix2.M34) + (matrix1.M24 * matrix2.M44),
1756				// M31
1757				(matrix1.M31 * matrix2.M11) + (matrix1.M32 * matrix2.M21) +
1758				(matrix1.M33 * matrix2.M31) + (matrix1.M34 * matrix2.M41),
1759				// M32
1760				(matrix1.M31 * matrix2.M12) + (matrix1.M32 * matrix2.M22) +
1761				(matrix1.M33 * matrix2.M32) + (matrix1.M34 * matrix2.M42),
1762				// M33
1763				(matrix1.M31 * matrix2.M13) + (matrix1.M32 * matrix2.M23) +
1764				(matrix1.M33 * matrix2.M33) + (matrix1.M34 * matrix2.M43),
1765				// M34
1766				(matrix1.M31 * matrix2.M14) + (matrix1.M32 * matrix2.M24) +
1767				(matrix1.M33 * matrix2.M34) + (matrix1.M34 * matrix2.M44),
1768				// M41
1769				(matrix1.M41 * matrix2.M11) + (matrix1.M42 * matrix2.M21) +
1770				(matrix1.M43 * matrix2.M31) + (matrix1.M44 * matrix2.M41),
1771				// M42
1772				(matrix1.M41 * matrix2.M12) + (matrix1.M42 * matrix2.M22) +
1773				(matrix1.M43 * matrix2.M32) + (matrix1.M44 * matrix2.M42),
1774				// M43
1775				(matrix1.M41 * matrix2.M13) + (matrix1.M42 * matrix2.M23) +
1776				(matrix1.M43 * matrix2.M33) + (matrix1.M44 * matrix2.M43),
1777				// M44
1778				(matrix1.M41 * matrix2.M14) + (matrix1.M42 * matrix2.M24) +
1779				(matrix1.M43 * matrix2.M34) + (matrix1.M44 * matrix2.M44));
1780		}
1781
1782		/// <summary>
1783		/// Operator to multiply matrix with a vector. This is very useful to
1784		/// transform vectors with the given matrix. Warning: This is slower than
1785		/// the ref version, this operator needs to copy 88 bytes around (1 matrix,
1786		/// 64 bytes and 2 vectors, each 12 bytes)
1787		/// </summary>
1788		/// <param name="matrix">Matrix to multiply</param>
1789		/// <param name="vector">Vector to multiply with</param>
1790		/// <returns>Resulting vector from the multiplication.</returns>
1791		public static Vector operator *(Matrix matrix, Vector vector)
1792		{
1793			// Note: Seems like calling XNA (OpenTK or SlimDX too?) is here slower!
1794			//Delta: 209ms (10 million calls)
1795			//XNA: 353ms (10 million calls)
1796			//OpenTK: 245ms (10 million calls)
1797			//SlimDX: 265ms (10 million calls)
1798			// Just use the normal vector * matrix formula (see Vector class)
1799			return new Vector(
1800				(vector.X * matrix.M11) + (vector.Y * matrix.M21) +
1801				(vector.Z * matrix.M31) + matrix.M41,
1802				(vector.X * matrix.M12) + (vector.Y * matrix.M22) +
1803				(vector.Z * matrix.M32) + matrix.M42,
1804				(vector.X * matrix.M13) + (vector.Y * matrix.M23) +
1805				(vector.Z * matrix.M33) + matrix.M43);
1806		}
1807
1808		/// <summary>
1809		/// Operator to multiply matrix with a scale factor. Try to use the ref
1810		/// version instead, this needs to copy 128 bytes around (2 matrices, each
1811		/// 64 bytes).
1812		/// </summary>
1813		/// <param name="matrix">Matrix to multiply</param>
1814		/// <param name="scaleFactor">Scale factor to multiply with</param>
1815		/// <returns>
1816		/// New matrix with all values multiplied by scaleFactor.
1817		/// </returns>
1818		public static Matrix operator *(Matrix matrix, float scaleFactor)
1819		{
1820			//Note: XNA, OpenTK, and SlimDX path