Mercurial > Math

view include/Matrix.h @ 15:ea32c94fc495

find changesets by author, revision, files, or words in the commit message
General minor improvements. - Spelling corrections. - Made the top level constants like PI static. - Changed instances of literal 1 and 0 to T(1) and T(0)
author Eris Caffee <discordia@eldalin.com>
date Sun, 19 May 2013 20:17:16 -0500
parents e00321d11fe1
children 2ac83ae9129b
line source
1 #ifndef MATRIX_H_
2 #define MATRIX_H_
3
4 #include "Vector.h"
5
6 #include <cmath>
7 #include <cstring>
8 #include <string>
9 #include <sstream>
10
11 namespace arda
12 {
13 namespace Math
14 {
15 //////////////////////////////////////////////////////////////////////////
16 // Maxtrices
17 //
18 // Supported operations on matrices
19 //
20 // []
21 // If you have a Matrix named m you can access it's member elements as
22 // m[0], m[1], m[2], etc. The indexing is COLUMN MAJOR. So for a 2x2
23 // matrix, m[0] and m[1] are the elements of the first column, not the
24 // first row.
25 // == != + += - -=
26 // Defined for operations on two matrices of the same type.
27 // TODO Is there a better way to do this? With integer types this is fine,
28 // but with floating point types this == and != are really worthless. Is it possible to
29 // create variant methods in a template class that vary by the type? Or do I just have
30 // to suck it up and create type specific external functions for comparing Matrices?
31 // * *= / /=
32 // Defined for scalar multiplication/division with int, float, and double.
33 // * is defined as a standalone template operator for the scalar * Matrix
34 // form.
35 // * *=
36 // Also defined for Matrix * Matrix
37 //
38 // assign() Assign the value of one matrix to another.
39 // getcol() Returns the vector that is column n of the matrix.
40 // to_string() Returns a printable string representation of the matrix.
41 // setcol() Set the vector that is column n of the matrix.
42 // Overloaded so that you can set from a vector or from
43 // specific values.
44 // setidentity() Sets a matrix to the identity matrix.
45 //
46 //
47 //
48 // The following are functions, not class methods. This is because I think
49 // that the notation det(m) and transpose(m), for example, is more natural
50 // to write and read than m.det() and m.transpose(). Using functional
51 // notation just more naturally mimics the standard mathematical notation.
52 //
53 // det() Calculates the determinant of a matrix.
54 // transpose() Returns the transpose of a matrix.
55
56 //////////////////////////////////////////////////////////////////////////
57 template <typename T>
58 class Matrix22
59 {
60 public:
61 T m[4];
62
63 // Constructors
64 Matrix22() {}
65 explicit Matrix22(T a) { m[0] = m[1] = m[2] = m[3] = a; }
66 Matrix22(T a0, T a1, T a2, T a3)
67 {
68 m[0] = a0; m[2] = a2;
69 m[1] = a1; m[3] = a3;
70 }
71
72 // Array indexing
73 // Remember: column major order is used.
74 inline T& operator[](unsigned int const i)
75 { assert (i<4); return m[i]; }
76 inline T operator[](unsigned int const i) const
77 { assert (i<4); return m[i]; }
78
79 // Assignment
80 inline Matrix22<T>& operator=(Matrix22<T> const & m2)
81 {
82 m[0] = m2[0]; m[2] = m2[2];
83 m[1] = m2[1]; m[3] = m2[3];
84 return *this;
85 }
86 inline Matrix22<T>& assign(T const a0, T const a1,
87 T const a2, T const a3)
88 {
89 m[0] = a0; m[2] = a2;
90 m[1] = a1; m[3] = a3;
91 return *this;
92 }
93
94 // Comparison
95 inline bool operator==(Matrix22<T> const & m2) const
96 {
97 return
98 m[0] == m2[0] && m[2] == m2[2] &&
99 m[1] == m2[1] && m[3] == m2[3];
100 }
101 inline bool operator!=(Matrix22<T> const & m2) const
102 { return ! (*this == m2); }
103
104 // Matrix addition
105 inline Matrix22<T>& operator+=(Matrix22<T> const & m2)
106 {
107 m[0] += m2[0]; m[2] += m2[2];
108 m[1] += m2[1]; m[3] += m2[3];
109 return *this;
110 }
111 inline Matrix22<T> operator+(Matrix22<T> const & m2) const
112 { return Matrix22<T>(*this) += m2; }
113
114 // Matrix subtraction
115 inline Matrix22<T>& operator-=(Matrix22<T> const & m2)
116 {
117 m[0] -= m2[0]; m[2] -= m2[2];
118 m[1] -= m2[1]; m[3] -= m2[3];
119 return *this;
120 }
121 inline Matrix22<T> operator-(Matrix22<T> const & m2) const
122 { return Matrix22<T>(*this) -= m2; }
123
124 // Scalar multiplication
125 inline Matrix22<T>& operator*=(int const a)
126 {
127 m[0] *= a; m[2] *= a;
128 m[1] *= a; m[3] *= a;
129 return *this;
130 }
131 inline Matrix22<T>& operator*=(float const a)
132 {
133 m[0] *= a; m[2] *= a;
134 m[1] *= a; m[3] *= a;
135 return *this;
136 }
137 inline Matrix22<T>& operator*=(double const a)
138 {
139 m[0] *= a; m[2] *= a;
140 m[1] *= a; m[3] *= a;
141 return *this;
142 }
143
144 // Scalar division
145 inline Matrix22<T>& operator/=(int const a)
146 {
147 assert(a!=0);
148 m[0] /= a; m[2] /= a;
149 m[1] /= a; m[3] /= a;
150 return *this;
151 }
152 inline Matrix22<T>& operator/=(float const a)
153 {
154 assert(a!=0);
155 m[0] /= a; m[2] /= a;
156 m[1] /= a; m[3] /= a;
157 return *this;
158 }
159 inline Matrix22<T>& operator/=(double const a)
160 {
161 assert(a!=0);
162 m[0] /= a; m[2] /= a;
163 m[1] /= a; m[3] /= a;
164 return *this;
165 }
166 inline Matrix22<T> operator/(int const a)
167 { return Matrix22<T>(*this) /= a; }
168 inline Matrix22<T> operator/(float const a)
169 { return Matrix22<T>(*this) /= a; }
170 inline Matrix22<T> operator/(double const a)
171 { return Matrix22<T>(*this) /= a; }
172
173 // Matrix multiplication
174 // Not sure if this should be inlined at all. Sure the compiler will
175 // probably ignore the inline request here, but maybe it won't and maybe
176 // that would be bad. Needs real testing.
177 inline Matrix22<T>& operator*=(Matrix22<T> const & m2)
178 {
179 const int size=2;
180 Matrix22<T> mres(0);
181 int i, j, k;
182 for (i=0; i<size; ++i)
183 for (j=0; j<size; ++j)
184 for (k=0; k<size; ++k)
185 mres[size*i+j] += (*this)[size*k+j] * m2[size*i+k];
186 *this = mres;
187 return *this;
188 }
189 inline Matrix22<T> operator*(Matrix22<T> const & m2) const
190 { return Matrix22<T>(*this) *= m2; }
191
192
193 inline Vector2<T> getcol(unsigned int const i) const
194 {
195 assert(i<2);
196 int const size=2;
197 return Vector2<T>(m[i*size], m[i*size+1]);
198 }
199
200 inline Vector2<T> getrow(unsigned int const i) const
201 {
202 assert(i<2);
203 int const size=2;
204 return Vector2<T>(m[i], m[i+size]);
205 }
206
207 inline Matrix22<T>& setcol(unsigned int const i, Vector2<T> const & v)
208 {
209 assert(i<2);
210 int const size=2;
211 m[i*size] = v[0];
212 m[i*size+1] = v[1];
213 return *this;
214 }
215
216 inline Matrix22<T>& setcol(unsigned int const i, T const a, T const b)
217 {
218 assert(i<2);
219 const int size=2;
220 m[i*size] = a;
221 m[i*size+1] = b;
222 return *this;
223 }
224
225 inline Matrix22<T>& setrow(unsigned int const i, Vector2<T> const & v)
226 {
227 assert(i<2);
228 int const size=2;
229 m[i] = v[0];
230 m[i+size] = v[1];
231 return *this;
232 }
233
234 inline Matrix22<T>& setrow(unsigned int const i, T const a, T const b)
235 {
236 assert(i<2);
237 int const size=2;
238 m[i] = a;
239 m[i+size] = b;
240 return *this;
241 }
242
243 std::string to_string(void) const;
244
245 inline Matrix22<T>& setidentity()
246 {
247 memset(m, 0, sizeof(m));
248 m[0] = m[3] = T(1);
249 return *this;
250 }
251
252 };
253
254 // Scalar multiplication continued
255 template <typename T>
256 inline Matrix22<T> operator*(Matrix22<T> const & m, int const a)
257 { return Matrix22<T>(m) *= a; }
258
259 template <typename T>
260 inline Matrix22<T> operator*(int const a, Matrix22<T> const & m)
261 { return Matrix22<T>(m) *= a; }
262
263 template <typename T>
264 inline Matrix22<T> operator*(Matrix22<T> const & m, float const a)
265 { return Matrix22<T>(m) *= a; }
266
267 template <typename T>
268 inline Matrix22<T> operator*(float const a, Matrix22<T> const & m)
269 { return Matrix22<T>(m) *= a; }
270
271 template <typename T>
272 inline Matrix22<T> operator*(Matrix22<T> const & m, double const a)
273 { return Matrix22<T>(m) *= a; }
274
275 template <typename T>
276 inline Matrix22<T> operator*(double const a, Matrix22<T> const & m)
277 { return Matrix22<T>(m) *= a; }
278
279
280
281
282 //////////////////////////////////////////////////////////////////////////
283 template <typename T>
284 class Matrix33
285 {
286 public:
287 T m[9];
288
289 // Constructors
290 Matrix33() {}
291 explicit Matrix33(T a)
292 {
293 // TODO: Is this really more efficient than memset?
294 m[0] = m[1] = m[2] = m[3] = m[4] = m[5] = m[6] = m[7] = m[8] = a;
295 }
296 Matrix33(T a0, T a1, T a2,
297 T a3, T a4, T a5,
298 T a6, T a7, T a8)
299 {
300 m[0] = a0; m[3] = a3; m[6] = a6;
301 m[1] = a1; m[4] = a4; m[7] = a7;
302 m[2] = a2; m[5] = a5; m[8] = a8;
303 }
304
305 // Array indexing
306 // Remember: column major order is used.
307 inline T& operator[](unsigned int const i)
308 { assert (i<9); return m[i]; }
309 inline T operator[](unsigned int const i) const
310 { assert (i<9); return m[i]; }
311
312 // Assignment
313 inline Matrix33<T>& operator=(Matrix33<T> const & m2)
314 {
315 m[0] = m2[0]; m[3] = m2[3]; m[6] = m2[6];
316 m[1] = m2[1]; m[4] = m2[4]; m[7] = m2[7];
317 m[2] = m2[2]; m[5] = m2[5]; m[8] = m2[8];
318 return *this;
319 }
320 inline Matrix33<T>& assign(T const a0, T const a1, T const a2,
321 T const a3, T const a4, T const a5,
322 T const a6, T const a7, T const a8)
323 {
324 m[0] = a0; m[3] = a3; m[6] = a6;
325 m[1] = a1; m[4] = a4; m[7] = a7;
326 m[2] = a2; m[5] = a5; m[8] = a8;
327 return *this;
328 }
329
330 // Comparison
331 inline bool operator==(Matrix33<T> const & m2) const
332 {
333 return
334 m[0] == m2[0] && m[3] == m2[3] && m[6] == m2[6] &&
335 m[1] == m2[1] && m[4] == m2[4] && m[7] == m2[7] &&
336 m[2] == m2[2] && m[5] == m2[5] && m[8] == m2[8];
337 }
338 inline bool operator!=(Matrix33<T> const & m2) const
339 { return ! (*this == m2); }
340
341 // Matrix addition
342 inline Matrix33<T>& operator+=(Matrix33<T> const & m2)
343 {
344 m[0] += m2[0]; m[3] += m2[3]; m[6] += m2[6];
345 m[1] += m2[1]; m[4] += m2[4]; m[7] += m2[7];
346 m[2] += m2[2]; m[5] += m2[5]; m[8] += m2[8];
347 return *this;
348 }
349 inline Matrix33<T> operator+(Matrix33<T> const & m2) const
350 { return Matrix33<T>(*this) += m2; }
351
352 // Matrix subtraction
353 inline Matrix33<T>& operator-=(Matrix33<T> const & m2)
354 {
355 m[0] -= m2[0]; m[3] -= m2[3]; m[6] -= m2[6];
356 m[1] -= m2[1]; m[4] -= m2[4]; m[7] -= m2[7];
357 m[2] -= m2[2]; m[5] -= m2[5]; m[8] -= m2[8];
358 return *this;
359 }
360 inline Matrix33<T> operator-(Matrix33<T> const & m2) const
361 { return Matrix33<T>(*this) -= m2; }
362
363 // Scalar multiplication
364 inline Matrix33<T>& operator*=(int const a)
365 {
366 m[0] *= a; m[3] *= a; m[6] *= a;
367 m[1] *= a; m[4] *= a; m[7] *= a;
368 m[2] *= a; m[5] *= a; m[8] *= a;
369 return *this;
370 }
371 inline Matrix33<T>& operator*=(float const a)
372 {
373 m[0] *= a; m[3] *= a; m[6] *= a;
374 m[1] *= a; m[4] *= a; m[7] *= a;
375 m[2] *= a; m[5] *= a; m[8] *= a;
376 return *this;
377 }
378 inline Matrix33<T>& operator*=(double const a)
379 {
380 m[0] *= a; m[3] *= a; m[6] *= a;
381 m[1] *= a; m[4] *= a; m[7] *= a;
382 m[2] *= a; m[5] *= a; m[8] *= a;
383 return *this;
384 }
385
386 // Scalar division
387 inline Matrix33<T>& operator/=(int const a)
388 {
389 assert(a!=0);
390 m[0] /= a; m[3] /= a; m[6] /= a;
391 m[1] /= a; m[4] /= a; m[7] /= a;
392 m[2] /= a; m[5] /= a; m[8] /= a;
393 return *this;
394 }
395 inline Matrix33<T>& operator/=(float const a)
396 {
397 assert(a!=0);
398 m[0] /= a; m[3] /= a; m[6] /= a;
399 m[1] /= a; m[4] /= a; m[7] /= a;
400 m[2] /= a; m[5] /= a; m[8] /= a;
401 return *this;
402 }
403 inline Matrix33<T>& operator/=(double const a)
404 {
405 assert(a!=0);
406 m[0] /= a; m[3] /= a; m[6] /= a;
407 m[1] /= a; m[4] /= a; m[7] /= a;
408 m[2] /= a; m[5] /= a; m[8] /= a;
409 return *this;
410 }
411 inline Matrix33<T> operator/(int const a)
412 { return Matrix33<T>(*this) /= a; }
413 inline Matrix33<T> operator/(float const a)
414 { return Matrix33<T>(*this) /= a; }
415 inline Matrix33<T> operator/(double const a)
416 { return Matrix33<T>(*this) /= a; }
417
418 // Matrix multiplication
419 // Not sure if this should be inlined at all. Sure the compiler will
420 // probably ignore the inline request here, but maybe it won't and maybe
421 // that would be bad. Needs real testing.
422 inline Matrix33<T>& operator*=(Matrix33<T> const & m2)
423 {
424 const int size=3;
425 Matrix33<T> mres(0);
426 int i, j, k;
427 for (i=0; i<size; ++i)
428 for (j=0; j<size; ++j)
429 for (k=0; k<size; ++k)
430 mres[size*i+j] += (*this)[size*k+j] * m2[size*i+k];
431 *this = mres;
432 return *this;
433 }
434 inline Matrix33<T> operator*(Matrix33<T> const & m2) const
435 { return Matrix33<T>(*this) *= m2; }
436
437
438 inline Vector3<T> getcol(unsigned int const i) const
439 {
440 assert(i<3);
441 const int size=3;
442 return Vector3<T>(m[i*size], m[i*size+1], m[i*size+2]);
443 }
444
445 inline Vector3<T> getrow(unsigned int const i) const
446 {
447 assert(i<3);
448 const int size=3;
449 return Vector3<T>(m[i], m[i+size], m[i+2*size]);
450 }
451
452 inline Matrix33<T>& setcol(unsigned int const i, Vector3<T> const & v)
453 {
454 assert(i<3);
455 const int size=3;
456 m[i*size] = v[0];
457 m[i*size+1] = v[1];
458 m[i*size+2] = v[2];
459 return *this;
460 }
461
462 inline Matrix33<T>& setcol(unsigned int const i,
463 T const a0, T const a1, T const a2)
464 {
465 assert(i<3);
466 const int size=3;
467 m[i*size] = a0;
468 m[i*size+1] = a1;
469 m[i*size+2] = a2;
470 return *this;
471 }
472
473 inline Matrix33<T>& setrow(unsigned int const i, Vector3<T> const & v)
474 {
475 assert(i<3);
476 const int size=3;
477 m[i] = v[0];
478 m[i+size] = v[1];
479 m[i+2*size] = v[2];
480 return *this;
481 }
482
483 inline Matrix33<T>& setrow(unsigned int const i,
484 T const a0, T const a1, T const a2)
485 {
486 assert(i<3);
487 const int size=3;
488 m[i] = a0;
489 m[i+size] = a1;
490 m[i+2*size] = a2;
491 return *this;
492 }
493
494 std::string to_string(void) const;
495
496 inline Matrix33<T>& setidentity()
497 {
498 memset(m, 0, sizeof(m));
499 m[0] = m[4] = m[8] = T(1);
500 return *this;
501 }
502
503 };
504
505 // Scalar multiplication continued
506 template <typename T>
507 inline Matrix33<T> operator*(Matrix33<T> const & m, int const a)
508 { return Matrix33<T>(m) *= a; }
509
510 template <typename T>
511 inline Matrix33<T> operator*(int const a, Matrix33<T> const & m)
512 { return Matrix33<T>(m) *= a; }
513
514 template <typename T>
515 inline Matrix33<T> operator*(Matrix33<T> const & m, float const a)
516 { return Matrix33<T>(m) *= a; }
517
518 template <typename T>
519 inline Matrix33<T> operator*(float const a, Matrix33<T> const & m)
520 { return Matrix33<T>(m) *= a; }
521
522 template <typename T>
523 inline Matrix33<T> operator*(Matrix33<T> const & m, double const a)
524 { return Matrix33<T>(m) *= a; }
525
526 template <typename T>
527 inline Matrix33<T> operator*(double const a, Matrix33<T> const & m)
528 { return Matrix33<T>(m) *= a; }
529
530
531
532
533 //////////////////////////////////////////////////////////////////////////
534 template <typename T>
535 class Matrix44
536 {
537 public:
538 T m[16];
539
540 // Constructors
541 Matrix44() {}
542 explicit Matrix44(T a)
543 {
544 // TODO: Is this really more efficient than memset?
545 m[0] = m[1] = m[2] = m[3] = m[4] = m[5] = m[6] = m[7] = m[8] =
546 m[9] = m[10] = m[11] = m[12] = m[13] = m[14] = m[15] = a;
547 }
548 Matrix44(T a0, T a1, T a2, T a3,
549 T a4, T a5, T a6, T a7,
550 T a8, T a9, T a10, T a11,
551 T a12, T a13, T a14, T a15)
552 {
553 m[0] = a0; m[4] = a4; m[8] = a8; m[12] = a12;
554 m[1] = a1; m[5] = a5; m[9] = a9; m[13] = a13;
555 m[2] = a2; m[6] = a6; m[10] = a10; m[14] = a14;
556 m[3] = a3; m[7] = a7; m[11] = a11; m[15] = a15;
557 }
558
559 // Array indexing
560 // Remember: column major order is used.
561 inline T& operator[](unsigned int const i)
562 { assert (i<16); return m[i]; }
563 inline T operator[](unsigned int const i) const
564 { assert (i<16); return m[i]; }
565
566 // Assignment
567 inline Matrix44<T>& operator=(Matrix44<T> const & m2)
568 {
569 m[0] = m2[0]; m[4] = m2[4]; m[8] = m2[8]; m[12] = m2[12];
570 m[1] = m2[1]; m[5] = m2[5]; m[9] = m2[9]; m[13] = m2[13];
571 m[2] = m2[2]; m[6] = m2[6]; m[10] = m2[10]; m[14] = m2[14];
572 m[3] = m2[3]; m[7] = m2[7]; m[11] = m2[11]; m[15] = m2[15];
573 return *this;
574 }
575 inline Matrix44<T>& assign(T const a0, T const a1, T const a2, T const a3,
576 T const a4, T const a5, T const a6, T const a7,
577 T const a8, T const a9, T const a10, T const a11,
578 T const a12, T const a13, T const a14, T const a15)
579 {
580 m[0] = a0; m[4] = a4; m[8] = a8; m[12] = a12;
581 m[1] = a1; m[5] = a5; m[9] = a9; m[13] = a13;
582 m[2] = a2; m[6] = a6; m[10] = a10; m[14] = a14;
583 m[3] = a3; m[7] = a7; m[11] = a11; m[15] = a15;
584 return *this;
585 }
586
587 // Comparison
588 inline bool operator==(Matrix44<T> const & m2) const
589 {
590 // TODO: Should this use memcmp?
591 return
592 m[0] == m2[0] && m[4] == m2[4] && m[8] == m2[8] && m[12] == m2[12] &&
593 m[1] == m2[1] && m[5] == m2[5] && m[9] == m2[9] && m[13] == m2[13] &&
594 m[2] == m2[2] && m[6] == m2[6] && m[10] == m2[10] && m[14] == m2[14] &&
595 m[3] == m2[3] && m[7] == m2[7] && m[11] == m2[11] && m[15] == m2[15];
596 }
597 inline bool operator!=(Matrix44<T> const & m2) const
598 { return ! (*this == m2); }
599
600 // Matrix addition
601 inline Matrix44<T>& operator+=(Matrix44<T> const & m2)
602 {
603 m[0] += m2[0]; m[4] += m2[4]; m[8] += m2[8]; m[12] += m2[12];
604 m[1] += m2[1]; m[5] += m2[5]; m[9] += m2[9]; m[13] += m2[13];
605 m[2] += m2[2]; m[6] += m2[6]; m[10] += m2[10]; m[14] += m2[14];
606 m[3] += m2[3]; m[7] += m2[7]; m[11] += m2[11]; m[15] += m2[15];
607 return *this;
608 }
609 inline Matrix44<T> operator+(Matrix44<T> const & m2) const
610 { return Matrix44<T>(*this) += m2; }
611
612 // Matrix subtraction
613 inline Matrix44<T>& operator-=(Matrix44<T> const & m2)
614 {
615 m[0] -= m2[0]; m[4] -= m2[4]; m[8] -= m2[8]; m[12] -= m2[12];
616 m[1] -= m2[1]; m[5] -= m2[5]; m[9] -= m2[9]; m[13] -= m2[13];
617 m[2] -= m2[2]; m[6] -= m2[6]; m[10] -= m2[10]; m[14] -= m2[14];
618 m[3] -= m2[3]; m[7] -= m2[7]; m[11] -= m2[11]; m[15] -= m2[15];
619 return *this;
620 }
621 inline Matrix44<T> operator-(Matrix44<T> const & m2) const
622 { return Matrix44<T>(*this) -= m2; }
623
624 // Scalar multiplication
625 inline Matrix44<T>& operator*=(int const a)
626 {
627 m[0] *= a; m[4] *= a; m[8] *= a; m[12] *= a;
628 m[1] *= a; m[5] *= a; m[9] *= a; m[13] *= a;
629 m[2] *= a; m[6] *= a; m[10] *= a; m[14] *= a;
630 m[3] *= a; m[7] *= a; m[11] *= a; m[15] *= a;
631 return *this;
632 }
633 inline Matrix44<T>& operator*=(float const a)
634 {
635 m[0] *= a; m[4] *= a; m[8] *= a; m[12] *= a;
636 m[1] *= a; m[5] *= a; m[9] *= a; m[13] *= a;
637 m[2] *= a; m[6] *= a; m[10] *= a; m[14] *= a;
638 m[3] *= a; m[7] *= a; m[11] *= a; m[15] *= a;
639 return *this;
640 }
641 inline Matrix44<T>& operator*=(double const a)
642 {
643 m[0] *= a; m[4] *= a; m[8] *= a; m[12] *= a;
644 m[1] *= a; m[5] *= a; m[9] *= a; m[13] *= a;
645 m[2] *= a; m[6] *= a; m[10] *= a; m[14] *= a;
646 m[3] *= a; m[7] *= a; m[11] *= a; m[15] *= a;
647 return *this;
648 }
649
650 // Scalar division
651 inline Matrix44<T>& operator/=(int const a)
652 {
653 assert(a!=0);
654 m[0] /= a; m[4] /= a; m[8] /= a; m[12] /= a;
655 m[1] /= a; m[5] /= a; m[9] /= a; m[13] /= a;
656 m[2] /= a; m[6] /= a; m[10] /= a; m[14] /= a;
657 m[3] /= a; m[7] /= a; m[11] /= a; m[15] /= a;
658 return *this;
659 }
660 inline Matrix44<T>& operator/=(float const a)
661 {
662 assert(a!=0);
663 m[0] /= a; m[4] /= a; m[8] /= a; m[12] /= a;
664 m[1] /= a; m[5] /= a; m[9] /= a; m[13] /= a;
665 m[2] /= a; m[6] /= a; m[10] /= a; m[14] /= a;
666 m[3] /= a; m[7] /= a; m[11] /= a; m[15] /= a;
667 return *this;
668 }
669 inline Matrix44<T>& operator/=(double const a)
670 {
671 assert(a!=0);
672 m[0] /= a; m[4] /= a; m[8] /= a; m[12] /= a;
673 m[1] /= a; m[5] /= a; m[9] /= a; m[13] /= a;
674 m[2] /= a; m[6] /= a; m[10] /= a; m[14] /= a;
675 m[3] /= a; m[7] /= a; m[11] /= a; m[15] /= a;
676 return *this;
677 }
678 inline Matrix44<T> operator/(int const a)
679 { return Matrix44<T>(*this) /= a; }
680 inline Matrix44<T> operator/(float const a)
681 { return Matrix44<T>(*this) /= a; }
682 inline Matrix44<T> operator/(double const a)
683 { return Matrix44<T>(*this) /= a; }
684
685 // Matrix multiplication
686 // Not sure if this should be inlined at all. Sure the compiler will
687 // probably ignore the inline request here, but maybe it won't and maybe
688 // that would be bad. Needs real testing.
689 inline Matrix44<T>& operator*=(Matrix44<T> const & m2)
690 {
691 const int size=4;
692 Matrix44<T> mres(0);
693 int i, j, k;
694 for (i=0; i<size; ++i)
695 for (j=0; j<size; ++j)
696 for (k=0; k<size; ++k)
697 mres[size*i+j] += (*this)[size*k+j] * m2[size*i+k];
698 *this = mres;
699 return *this;
700 }
701 inline Matrix44<T> operator*(Matrix44<T> const & m2) const
702 { return Matrix44<T>(*this) *= m2; }
703
704
705 inline Vector4<T> getcol(unsigned int const i) const
706 {
707 assert(i<4);
708 const int size=4;
709 return Vector4<T>(m[i*size], m[i*size+1], m[i*size+2], m[i*size+3]);
710 }
711
712 inline Vector4<T> getrow(unsigned int const i) const
713 {
714 assert(i<4);
715 const int size=4;
716 return Vector4<T>(m[i], m[i+size], m[i+2*size], m[i+3*size]);
717 }
718
719 inline Matrix44<T>& setcol(unsigned int const i, Vector4<T> const & v)
720 {
721 assert(i<4);
722 const int size=4;
723 m[i*size] = v[0];
724 m[i*size+1] = v[1];
725 m[i*size+2] = v[2];
726 m[i*size+3] = v[3];
727 return *this;
728 }
729
730 inline Matrix44<T>& setcol(unsigned int const i,
731 T const a0, T const a1, T const a2, T const a3)
732 {
733 assert(i<4);
734 const int size=4;
735 m[i*size] = a0;
736 m[i*size+1] = a1;
737 m[i*size+2] = a2;
738 m[i*size+3] = a3;
739 return *this;
740 }
741
742 inline Matrix44<T>& setrow(unsigned int const i, Vector4<T> const & v)
743 {
744 assert(i<4);
745 const int size=4;
746 m[i] = v[0];
747 m[i+size] = v[1];
748 m[i+2*size] = v[2];
749 m[i+3*size] = v[3];
750 return *this;
751 }
752
753 inline Matrix44<T>& setrow(unsigned int const i,
754 T const a0, T const a1, T const a2, T const a3)
755 {
756 assert(i<4);
757 const int size=4;
758 m[i] = a0;
759 m[i+size] = a1;
760 m[i+2*size] = a2;
761 m[i+3*size] = a3;
762 return *this;
763 }
764
765 std::string to_string(void) const;
766
767 inline Matrix44<T>& setidentity()
768 {
769 memset(m, 0, sizeof(m));
770 m[0] = m[5] = m[10] = m[15] = T(1);
771 return *this;
772 }
773
774 };
775
776 // Scalar multiplication continued
777 template <typename T>
778 inline Matrix44<T> operator*(Matrix44<T> const & m, int const a)
779 { return Matrix44<T>(m) *= a; }
780
781 template <typename T>
782 inline Matrix44<T> operator*(int const a, Matrix44<T> const & m)
783 { return Matrix44<T>(m) *= a; }
784
785 template <typename T>
786 inline Matrix44<T> operator*(Matrix44<T> const & m, float const a)
787 { return Matrix44<T>(m) *= a; }
788
789 template <typename T>
790 inline Matrix44<T> operator*(float const a, Matrix44<T> const & m)
791 { return Matrix44<T>(m) *= a; }
792
793 template <typename T>
794 inline Matrix44<T> operator*(Matrix44<T> const & m, double const a)
795 { return Matrix44<T>(m) *= a; }
796
797 template <typename T>
798 inline Matrix44<T> operator*(double const a, Matrix44<T> const & m)
799 { return Matrix44<T>(m) *= a; }
800
801
802
803
804 //////////////////////////////////////////////////////////////////////////
805 typedef Matrix22<int> Matrix22i;
806 typedef Matrix22<float> Matrix22f;
807 typedef Matrix22<double> Matrix22d;
808
809 typedef Matrix33<int> Matrix33i;
810 typedef Matrix33<float> Matrix33f;
811 typedef Matrix33<double> Matrix33d;
812
813 typedef Matrix44<int> Matrix44i;
814 typedef Matrix44<float> Matrix44f;
815 typedef Matrix44<double> Matrix44d;
816
817
818 //////////////////////////////////////////////////////////////////////////
819 // Matrix functions
820 // TODO: The determinant functions suffer from some bad round off error.
821 // If it significant? I don't know. None of the other game engines
822 // I've checked bother to do anything to calculate det more accurately.
823 // So maybe it just doesn't matter enough for game purposes.
824
825 template <typename T>
826 inline double det(Matrix22<T> const & m)
827 {
828 return (double) (m[0] * m[3]) - (double) (m[1] * m[2]);
829 }
830
831 template <typename T>
832 inline Matrix22<T> transpose(Matrix22<T> const & m)
833 {
834 const int size=2;
835 Matrix22<T> mres(0);
836 int i, j;
837 for (i=0; i<size; ++i)
838 for (j=0; j<size; j++)
839 mres[size*i+j] = m[size*j+i];
840 return mres;
841 }
842
843 template <typename T>
844 inline double det(Matrix33<T> const & m)
845 {
846 return m[0] * ( (double) m[4] * m[8] - (double) m[5] * m[7])
847 - m[3] * ( (double) m[1] * m[8] - (double) m[2] * m[7])
848 + m[6] * ( (double) m[1] * m[5] - (double) m[2] * m[4]);
849 }
850 template <typename T>
851 inline Matrix33<T> transpose(Matrix33<T> const & m)
852 {
853 const int size=3;
854 Matrix33<T> mres(0);
855 int i, j;
856 for (i=0; i<size; ++i)
857 for (j=0; j<size; j++)
858 mres[size*i+j] = m[size*j+i];
859 return mres;
860 }
861
862 template <typename T>
863 inline double det(Matrix44<T> const & m)
864 {
865 return m[0] * ( m[5] * (m[10] * m[15] - m[11] * m[14])
866 - m[9] * (m[6] * m[15] - m[7] * m[14])
867 + m[13] * (m[6] * m[11] - m[7] * m[10]))
868 - m[4] * ( m[1] * (m[10] * m[15] - m[11] * m[14])
869 - m[9] * (m[2] * m[15] - m[3] * m[14])
870 + m[13] * (m[2] * m[11] - m[3] * m[10]))
871 + m[8] * ( m[1] * (m[6] * m[15] - m[7] * m[14])
872 - m[5] * (m[2] * m[15] - m[3] * m[14])
873 + m[13] * (m[2] * m[7] - m[3] * m[6]))
874 + m[12] * ( m[1] * (m[6] * m[11] - m[7] * m[10])
875 - m[5] * (m[2] * m[11] - m[3] * m[10])
876 + m[9] * (m[2] * m[7] - m[3] * m[6]));
877 }
878 template <typename T>
879 inline Matrix44<T> transpose(Matrix44<T> const & m)
880 {
881 const int size=4;
882 Matrix44<T> mres(0);
883 int i, j;
884 for (i=0; i<size; ++i)
885 for (j=0; j<size; j++)
886 mres[size*i+j] = m[size*j+i];
887 return mres;
888 }
889
890 } // namespace Math
891 } // namespace arda
892
893
894 ////////////////////////////////////////////////////////////////////////////////
895 // Matrix methods
896
897 template <typename T>
898 std::string arda::Math::Matrix22<T>::to_string(void) const
899 {
900 std::stringstream ss;
901 ss << "[ "
902 "[ " << m[0] << ", " << m[1] << " ], "
903 "[ " << m[2] << ", " << m[3] << " ] ]";
904 return ss.str();
905 }
906
907 template <typename T>
908 std::string arda::Math::Matrix33<T>::to_string(void) const
909 {
910 std::stringstream ss;
911 ss << "[ "
912 "[ " << m[0] << ", " << m[1] << ", " << m[2] << " ], "
913 "[ " << m[3] << ", " << m[4] << ", " << m[5] << " ], "
914 "[ " << m[6] << ", " << m[7] << ", " << m[8] << " ], ""]";
915 return ss.str();
916 }
917
918
919 template <typename T>
920 std::string arda::Math::Matrix44<T>::to_string(void) const
921 {
922 std::stringstream ss;
923 ss << "[ "
924 "[ " << m[0] << ", " << m[1] << ", " << m[2] << ", " << m[3] << " ], "
925 "[ " << m[4] << ", " << m[5] << ", " << m[6] << ", " << m[7] << " ], "
926 "[ " << m[8] << ", " << m[9] << ", " << m[10] << ", " << m[11] << " ], "
927 "[ " << m[12] << ", " << m[13] << ", " << m[14] << ", " << m[15] << " ], ""]";
928 return ss.str();
929 }
930
931 #endif // MATRIX_H_