OpenVDB  7.0.0
Mat4.h
Go to the documentation of this file.
1 // Copyright Contributors to the OpenVDB Project
2 // SPDX-License-Identifier: MPL-2.0
3 
4 #ifndef OPENVDB_MATH_MAT4_H_HAS_BEEN_INCLUDED
5 #define OPENVDB_MATH_MAT4_H_HAS_BEEN_INCLUDED
6 
7 #include <openvdb/Exceptions.h>
8 #include <openvdb/Platform.h>
9 #include "Math.h"
10 #include "Mat3.h"
11 #include "Vec3.h"
12 #include "Vec4.h"
13 #include <algorithm> // for std::copy(), std::swap()
14 #include <cassert>
15 #include <iomanip>
16 #include <cmath>
17 
18 
19 namespace openvdb {
21 namespace OPENVDB_VERSION_NAME {
22 namespace math {
23 
24 template<typename T> class Vec4;
25 
26 
29 template<typename T>
30 class Mat4: public Mat<4, T>
31 {
32 public:
34  using value_type = T;
35  using ValueType = T;
36  using MyBase = Mat<4, T>;
37 
39  Mat4() {}
40 
42 
48  template<typename Source>
49  Mat4(Source *a)
50  {
51  for (int i = 0; i < 16; i++) {
52  MyBase::mm[i] = static_cast<T>(a[i]);
53  }
54  }
55 
57 
63  template<typename Source>
64  Mat4(Source a, Source b, Source c, Source d,
65  Source e, Source f, Source g, Source h,
66  Source i, Source j, Source k, Source l,
67  Source m, Source n, Source o, Source p)
68  {
69  MyBase::mm[ 0] = static_cast<T>(a);
70  MyBase::mm[ 1] = static_cast<T>(b);
71  MyBase::mm[ 2] = static_cast<T>(c);
72  MyBase::mm[ 3] = static_cast<T>(d);
73 
74  MyBase::mm[ 4] = static_cast<T>(e);
75  MyBase::mm[ 5] = static_cast<T>(f);
76  MyBase::mm[ 6] = static_cast<T>(g);
77  MyBase::mm[ 7] = static_cast<T>(h);
78 
79  MyBase::mm[ 8] = static_cast<T>(i);
80  MyBase::mm[ 9] = static_cast<T>(j);
81  MyBase::mm[10] = static_cast<T>(k);
82  MyBase::mm[11] = static_cast<T>(l);
83 
84  MyBase::mm[12] = static_cast<T>(m);
85  MyBase::mm[13] = static_cast<T>(n);
86  MyBase::mm[14] = static_cast<T>(o);
87  MyBase::mm[15] = static_cast<T>(p);
88  }
89 
92  template<typename Source>
93  Mat4(const Vec4<Source> &v1, const Vec4<Source> &v2,
94  const Vec4<Source> &v3, const Vec4<Source> &v4, bool rows = true)
95  {
96  if (rows) {
97  this->setRows(v1, v2, v3, v4);
98  } else {
99  this->setColumns(v1, v2, v3, v4);
100  }
101  }
102 
104  Mat4(const Mat<4, T> &m)
105  {
106  for (int i = 0; i < 4; ++i) {
107  for (int j = 0; j < 4; ++j) {
108  MyBase::mm[i*4 + j] = m[i][j];
109  }
110  }
111  }
112 
114  template<typename Source>
115  explicit Mat4(const Mat4<Source> &m)
116  {
117  const Source *src = m.asPointer();
118 
119  for (int i=0; i<16; ++i) {
120  MyBase::mm[i] = static_cast<T>(src[i]);
121  }
122  }
123 
125  static const Mat4<T>& identity() {
126  static const Mat4<T> sIdentity = Mat4<T>(
127  1, 0, 0, 0,
128  0, 1, 0, 0,
129  0, 0, 1, 0,
130  0, 0, 0, 1
131  );
132  return sIdentity;
133  }
134 
136  static const Mat4<T>& zero() {
137  static const Mat4<T> sZero = Mat4<T>(
138  0, 0, 0, 0,
139  0, 0, 0, 0,
140  0, 0, 0, 0,
141  0, 0, 0, 0
142  );
143  return sZero;
144  }
145 
147  void setRow(int i, const Vec4<T> &v)
148  {
149  // assert(i>=0 && i<4);
150  int i4 = i * 4;
151  MyBase::mm[i4+0] = v[0];
152  MyBase::mm[i4+1] = v[1];
153  MyBase::mm[i4+2] = v[2];
154  MyBase::mm[i4+3] = v[3];
155  }
156 
158  Vec4<T> row(int i) const
159  {
160  // assert(i>=0 && i<3);
161  return Vec4<T>((*this)(i,0), (*this)(i,1), (*this)(i,2), (*this)(i,3));
162  }
163 
165  void setCol(int j, const Vec4<T>& v)
166  {
167  // assert(j>=0 && j<4);
168  MyBase::mm[ 0+j] = v[0];
169  MyBase::mm[ 4+j] = v[1];
170  MyBase::mm[ 8+j] = v[2];
171  MyBase::mm[12+j] = v[3];
172  }
173 
175  Vec4<T> col(int j) const
176  {
177  // assert(j>=0 && j<4);
178  return Vec4<T>((*this)(0,j), (*this)(1,j), (*this)(2,j), (*this)(3,j));
179  }
180 
182  T* operator[](int i) { return &(MyBase::mm[i<<2]); }
185  const T* operator[](int i) const { return &(MyBase::mm[i<<2]); }
187 
189  T* asPointer() {return MyBase::mm;}
190  const T* asPointer() const {return MyBase::mm;}
191 
195  T& operator()(int i, int j)
196  {
197  // assert(i>=0 && i<4);
198  // assert(j>=0 && j<4);
199  return MyBase::mm[4*i+j];
200  }
201 
205  T operator()(int i, int j) const
206  {
207  // assert(i>=0 && i<4);
208  // assert(j>=0 && j<4);
209  return MyBase::mm[4*i+j];
210  }
211 
213  void setRows(const Vec4<T> &v1, const Vec4<T> &v2,
214  const Vec4<T> &v3, const Vec4<T> &v4)
215  {
216  MyBase::mm[ 0] = v1[0];
217  MyBase::mm[ 1] = v1[1];
218  MyBase::mm[ 2] = v1[2];
219  MyBase::mm[ 3] = v1[3];
220 
221  MyBase::mm[ 4] = v2[0];
222  MyBase::mm[ 5] = v2[1];
223  MyBase::mm[ 6] = v2[2];
224  MyBase::mm[ 7] = v2[3];
225 
226  MyBase::mm[ 8] = v3[0];
227  MyBase::mm[ 9] = v3[1];
228  MyBase::mm[10] = v3[2];
229  MyBase::mm[11] = v3[3];
230 
231  MyBase::mm[12] = v4[0];
232  MyBase::mm[13] = v4[1];
233  MyBase::mm[14] = v4[2];
234  MyBase::mm[15] = v4[3];
235  }
236 
238  void setColumns(const Vec4<T> &v1, const Vec4<T> &v2,
239  const Vec4<T> &v3, const Vec4<T> &v4)
240  {
241  MyBase::mm[ 0] = v1[0];
242  MyBase::mm[ 1] = v2[0];
243  MyBase::mm[ 2] = v3[0];
244  MyBase::mm[ 3] = v4[0];
245 
246  MyBase::mm[ 4] = v1[1];
247  MyBase::mm[ 5] = v2[1];
248  MyBase::mm[ 6] = v3[1];
249  MyBase::mm[ 7] = v4[1];
250 
251  MyBase::mm[ 8] = v1[2];
252  MyBase::mm[ 9] = v2[2];
253  MyBase::mm[10] = v3[2];
254  MyBase::mm[11] = v4[2];
255 
256  MyBase::mm[12] = v1[3];
257  MyBase::mm[13] = v2[3];
258  MyBase::mm[14] = v3[3];
259  MyBase::mm[15] = v4[3];
260  }
261 
262  // Set this matrix to zero
263  void setZero()
264  {
265  MyBase::mm[ 0] = 0;
266  MyBase::mm[ 1] = 0;
267  MyBase::mm[ 2] = 0;
268  MyBase::mm[ 3] = 0;
269  MyBase::mm[ 4] = 0;
270  MyBase::mm[ 5] = 0;
271  MyBase::mm[ 6] = 0;
272  MyBase::mm[ 7] = 0;
273  MyBase::mm[ 8] = 0;
274  MyBase::mm[ 9] = 0;
275  MyBase::mm[10] = 0;
276  MyBase::mm[11] = 0;
277  MyBase::mm[12] = 0;
278  MyBase::mm[13] = 0;
279  MyBase::mm[14] = 0;
280  MyBase::mm[15] = 0;
281  }
282 
284  void setIdentity()
285  {
286  MyBase::mm[ 0] = 1;
287  MyBase::mm[ 1] = 0;
288  MyBase::mm[ 2] = 0;
289  MyBase::mm[ 3] = 0;
290 
291  MyBase::mm[ 4] = 0;
292  MyBase::mm[ 5] = 1;
293  MyBase::mm[ 6] = 0;
294  MyBase::mm[ 7] = 0;
295 
296  MyBase::mm[ 8] = 0;
297  MyBase::mm[ 9] = 0;
298  MyBase::mm[10] = 1;
299  MyBase::mm[11] = 0;
300 
301  MyBase::mm[12] = 0;
302  MyBase::mm[13] = 0;
303  MyBase::mm[14] = 0;
304  MyBase::mm[15] = 1;
305  }
306 
307 
309  void setMat3(const Mat3<T> &m)
310  {
311  for (int i = 0; i < 3; i++)
312  for (int j=0; j < 3; j++)
313  MyBase::mm[i*4+j] = m[i][j];
314  }
315 
316  Mat3<T> getMat3() const
317  {
318  Mat3<T> m;
319 
320  for (int i = 0; i < 3; i++)
321  for (int j = 0; j < 3; j++)
322  m[i][j] = MyBase::mm[i*4+j];
323 
324  return m;
325  }
326 
329  {
330  return Vec3<T>(MyBase::mm[12], MyBase::mm[13], MyBase::mm[14]);
331  }
332 
333  void setTranslation(const Vec3<T> &t)
334  {
335  MyBase::mm[12] = t[0];
336  MyBase::mm[13] = t[1];
337  MyBase::mm[14] = t[2];
338  }
339 
341  template<typename Source>
342  const Mat4& operator=(const Mat4<Source> &m)
343  {
344  const Source *src = m.asPointer();
345 
346  // don't suppress warnings when assigning from different numerical types
347  std::copy(src, (src + this->numElements()), MyBase::mm);
348  return *this;
349  }
350 
352  bool eq(const Mat4 &m, T eps=1.0e-8) const
353  {
354  for (int i = 0; i < 16; i++) {
355  if (!isApproxEqual(MyBase::mm[i], m.mm[i], eps))
356  return false;
357  }
358  return true;
359  }
360 
363  {
364  return Mat4<T>(
365  -MyBase::mm[ 0], -MyBase::mm[ 1], -MyBase::mm[ 2], -MyBase::mm[ 3],
366  -MyBase::mm[ 4], -MyBase::mm[ 5], -MyBase::mm[ 6], -MyBase::mm[ 7],
367  -MyBase::mm[ 8], -MyBase::mm[ 9], -MyBase::mm[10], -MyBase::mm[11],
368  -MyBase::mm[12], -MyBase::mm[13], -MyBase::mm[14], -MyBase::mm[15]
369  );
370  } // trivial
371 
373  template <typename S>
374  const Mat4<T>& operator*=(S scalar)
375  {
376  MyBase::mm[ 0] *= scalar;
377  MyBase::mm[ 1] *= scalar;
378  MyBase::mm[ 2] *= scalar;
379  MyBase::mm[ 3] *= scalar;
380 
381  MyBase::mm[ 4] *= scalar;
382  MyBase::mm[ 5] *= scalar;
383  MyBase::mm[ 6] *= scalar;
384  MyBase::mm[ 7] *= scalar;
385 
386  MyBase::mm[ 8] *= scalar;
387  MyBase::mm[ 9] *= scalar;
388  MyBase::mm[10] *= scalar;
389  MyBase::mm[11] *= scalar;
390 
391  MyBase::mm[12] *= scalar;
392  MyBase::mm[13] *= scalar;
393  MyBase::mm[14] *= scalar;
394  MyBase::mm[15] *= scalar;
395  return *this;
396  }
397 
399  template <typename S>
400  const Mat4<T> &operator+=(const Mat4<S> &m1)
401  {
402  const S* s = m1.asPointer();
403 
404  MyBase::mm[ 0] += s[ 0];
405  MyBase::mm[ 1] += s[ 1];
406  MyBase::mm[ 2] += s[ 2];
407  MyBase::mm[ 3] += s[ 3];
408 
409  MyBase::mm[ 4] += s[ 4];
410  MyBase::mm[ 5] += s[ 5];
411  MyBase::mm[ 6] += s[ 6];
412  MyBase::mm[ 7] += s[ 7];
413 
414  MyBase::mm[ 8] += s[ 8];
415  MyBase::mm[ 9] += s[ 9];
416  MyBase::mm[10] += s[10];
417  MyBase::mm[11] += s[11];
418 
419  MyBase::mm[12] += s[12];
420  MyBase::mm[13] += s[13];
421  MyBase::mm[14] += s[14];
422  MyBase::mm[15] += s[15];
423 
424  return *this;
425  }
426 
428  template <typename S>
429  const Mat4<T> &operator-=(const Mat4<S> &m1)
430  {
431  const S* s = m1.asPointer();
432 
433  MyBase::mm[ 0] -= s[ 0];
434  MyBase::mm[ 1] -= s[ 1];
435  MyBase::mm[ 2] -= s[ 2];
436  MyBase::mm[ 3] -= s[ 3];
437 
438  MyBase::mm[ 4] -= s[ 4];
439  MyBase::mm[ 5] -= s[ 5];
440  MyBase::mm[ 6] -= s[ 6];
441  MyBase::mm[ 7] -= s[ 7];
442 
443  MyBase::mm[ 8] -= s[ 8];
444  MyBase::mm[ 9] -= s[ 9];
445  MyBase::mm[10] -= s[10];
446  MyBase::mm[11] -= s[11];
447 
448  MyBase::mm[12] -= s[12];
449  MyBase::mm[13] -= s[13];
450  MyBase::mm[14] -= s[14];
451  MyBase::mm[15] -= s[15];
452 
453  return *this;
454  }
455 
457  template <typename S>
458  const Mat4<T> &operator*=(const Mat4<S> &m1)
459  {
460  Mat4<T> m0(*this);
461 
462  const T* s0 = m0.asPointer();
463  const S* s1 = m1.asPointer();
464 
465  for (int i = 0; i < 4; i++) {
466  int i4 = 4 * i;
467  MyBase::mm[i4+0] = static_cast<T>(s0[i4+0] * s1[ 0] +
468  s0[i4+1] * s1[ 4] +
469  s0[i4+2] * s1[ 8] +
470  s0[i4+3] * s1[12]);
471 
472  MyBase::mm[i4+1] = static_cast<T>(s0[i4+0] * s1[ 1] +
473  s0[i4+1] * s1[ 5] +
474  s0[i4+2] * s1[ 9] +
475  s0[i4+3] * s1[13]);
476 
477  MyBase::mm[i4+2] = static_cast<T>(s0[i4+0] * s1[ 2] +
478  s0[i4+1] * s1[ 6] +
479  s0[i4+2] * s1[10] +
480  s0[i4+3] * s1[14]);
481 
482  MyBase::mm[i4+3] = static_cast<T>(s0[i4+0] * s1[ 3] +
483  s0[i4+1] * s1[ 7] +
484  s0[i4+2] * s1[11] +
485  s0[i4+3] * s1[15]);
486  }
487  return *this;
488  }
489 
491  Mat4 transpose() const
492  {
493  return Mat4<T>(
494  MyBase::mm[ 0], MyBase::mm[ 4], MyBase::mm[ 8], MyBase::mm[12],
495  MyBase::mm[ 1], MyBase::mm[ 5], MyBase::mm[ 9], MyBase::mm[13],
496  MyBase::mm[ 2], MyBase::mm[ 6], MyBase::mm[10], MyBase::mm[14],
497  MyBase::mm[ 3], MyBase::mm[ 7], MyBase::mm[11], MyBase::mm[15]
498  );
499  }
500 
501 
504  Mat4 inverse(T tolerance = 0) const
505  {
506  //
507  // inv [ A | b ] = [ E | f ] A: 3x3, b: 3x1, c': 1x3 d: 1x1
508  // [ c' | d ] [ g' | h ]
509  //
510  // If A is invertible use
511  //
512  // E = A^-1 + p*h*r
513  // p = A^-1 * b
514  // f = -p * h
515  // g' = -h * c'
516  // h = 1 / (d - c'*p)
517  // r' = c'*A^-1
518  //
519  // Otherwise use gauss-jordan elimination
520  //
521 
522  //
523  // We create this alias to ourself so we can easily use own subscript
524  // operator.
525  const Mat4<T>& m(*this);
526 
527  T m0011 = m[0][0] * m[1][1];
528  T m0012 = m[0][0] * m[1][2];
529  T m0110 = m[0][1] * m[1][0];
530  T m0210 = m[0][2] * m[1][0];
531  T m0120 = m[0][1] * m[2][0];
532  T m0220 = m[0][2] * m[2][0];
533 
534  T detA = m0011 * m[2][2] - m0012 * m[2][1] - m0110 * m[2][2]
535  + m0210 * m[2][1] + m0120 * m[1][2] - m0220 * m[1][1];
536 
537  bool hasPerspective =
538  (!isExactlyEqual(m[0][3], T(0.0)) ||
539  !isExactlyEqual(m[1][3], T(0.0)) ||
540  !isExactlyEqual(m[2][3], T(0.0)) ||
541  !isExactlyEqual(m[3][3], T(1.0)));
542 
543  T det;
544  if (hasPerspective) {
545  det = m[0][3] * det3(m, 1,2,3, 0,2,1)
546  + m[1][3] * det3(m, 2,0,3, 0,2,1)
547  + m[2][3] * det3(m, 3,0,1, 0,2,1)
548  + m[3][3] * detA;
549  } else {
550  det = detA * m[3][3];
551  }
552 
553  Mat4<T> inv;
554  bool invertible;
555 
556  if (isApproxEqual(det,T(0.0),tolerance)) {
557  invertible = false;
558 
559  } else if (isApproxEqual(detA,T(0.0),T(1e-8))) {
560  // det is too small to rely on inversion by subblocks
561  invertible = m.invert(inv, tolerance);
562 
563  } else {
564  invertible = true;
565  detA = 1.0 / detA;
566 
567  //
568  // Calculate A^-1
569  //
570  inv[0][0] = detA * ( m[1][1] * m[2][2] - m[1][2] * m[2][1]);
571  inv[0][1] = detA * (-m[0][1] * m[2][2] + m[0][2] * m[2][1]);
572  inv[0][2] = detA * ( m[0][1] * m[1][2] - m[0][2] * m[1][1]);
573 
574  inv[1][0] = detA * (-m[1][0] * m[2][2] + m[1][2] * m[2][0]);
575  inv[1][1] = detA * ( m[0][0] * m[2][2] - m0220);
576  inv[1][2] = detA * ( m0210 - m0012);
577 
578  inv[2][0] = detA * ( m[1][0] * m[2][1] - m[1][1] * m[2][0]);
579  inv[2][1] = detA * ( m0120 - m[0][0] * m[2][1]);
580  inv[2][2] = detA * ( m0011 - m0110);
581 
582  if (hasPerspective) {
583  //
584  // Calculate r, p, and h
585  //
586  Vec3<T> r;
587  r[0] = m[3][0] * inv[0][0] + m[3][1] * inv[1][0]
588  + m[3][2] * inv[2][0];
589  r[1] = m[3][0] * inv[0][1] + m[3][1] * inv[1][1]
590  + m[3][2] * inv[2][1];
591  r[2] = m[3][0] * inv[0][2] + m[3][1] * inv[1][2]
592  + m[3][2] * inv[2][2];
593 
594  Vec3<T> p;
595  p[0] = inv[0][0] * m[0][3] + inv[0][1] * m[1][3]
596  + inv[0][2] * m[2][3];
597  p[1] = inv[1][0] * m[0][3] + inv[1][1] * m[1][3]
598  + inv[1][2] * m[2][3];
599  p[2] = inv[2][0] * m[0][3] + inv[2][1] * m[1][3]
600  + inv[2][2] * m[2][3];
601 
602  T h = m[3][3] - p.dot(Vec3<T>(m[3][0],m[3][1],m[3][2]));
603  if (isApproxEqual(h,T(0.0),tolerance)) {
604  invertible = false;
605 
606  } else {
607  h = 1.0 / h;
608 
609  //
610  // Calculate h, g, and f
611  //
612  inv[3][3] = h;
613  inv[3][0] = -h * r[0];
614  inv[3][1] = -h * r[1];
615  inv[3][2] = -h * r[2];
616 
617  inv[0][3] = -h * p[0];
618  inv[1][3] = -h * p[1];
619  inv[2][3] = -h * p[2];
620 
621  //
622  // Calculate E
623  //
624  p *= h;
625  inv[0][0] += p[0] * r[0];
626  inv[0][1] += p[0] * r[1];
627  inv[0][2] += p[0] * r[2];
628  inv[1][0] += p[1] * r[0];
629  inv[1][1] += p[1] * r[1];
630  inv[1][2] += p[1] * r[2];
631  inv[2][0] += p[2] * r[0];
632  inv[2][1] += p[2] * r[1];
633  inv[2][2] += p[2] * r[2];
634  }
635  } else {
636  // Equations are much simpler in the non-perspective case
637  inv[3][0] = - (m[3][0] * inv[0][0] + m[3][1] * inv[1][0]
638  + m[3][2] * inv[2][0]);
639  inv[3][1] = - (m[3][0] * inv[0][1] + m[3][1] * inv[1][1]
640  + m[3][2] * inv[2][1]);
641  inv[3][2] = - (m[3][0] * inv[0][2] + m[3][1] * inv[1][2]
642  + m[3][2] * inv[2][2]);
643  inv[0][3] = 0.0;
644  inv[1][3] = 0.0;
645  inv[2][3] = 0.0;
646  inv[3][3] = 1.0;
647  }
648  }
649 
650  if (!invertible) OPENVDB_THROW(ArithmeticError, "Inversion of singular 4x4 matrix");
651  return inv;
652  }
653 
654 
656  T det() const
657  {
658  const T *ap;
659  Mat3<T> submat;
660  T det;
661  T *sp;
662  int i, j, k, sign;
663 
664  det = 0;
665  sign = 1;
666  for (i = 0; i < 4; i++) {
667  ap = &MyBase::mm[ 0];
668  sp = submat.asPointer();
669  for (j = 0; j < 4; j++) {
670  for (k = 0; k < 4; k++) {
671  if ((k != i) && (j != 0)) {
672  *sp++ = *ap;
673  }
674  ap++;
675  }
676  }
677 
678  det += sign * MyBase::mm[i] * submat.det();
679  sign = -sign;
680  }
681 
682  return det;
683  }
684 
686  static Mat4 translation(const Vec3d& v)
687  {
688  return Mat4(
689  T(1), T(0), T(0), T(0),
690  T(0), T(1), T(0), T(0),
691  T(0), T(0), T(1), T(0),
692  T(v.x()), T(v.y()),T(v.z()), T(1));
693  }
694 
696  template <typename T0>
697  void setToTranslation(const Vec3<T0>& v)
698  {
699  MyBase::mm[ 0] = 1;
700  MyBase::mm[ 1] = 0;
701  MyBase::mm[ 2] = 0;
702  MyBase::mm[ 3] = 0;
703 
704  MyBase::mm[ 4] = 0;
705  MyBase::mm[ 5] = 1;
706  MyBase::mm[ 6] = 0;
707  MyBase::mm[ 7] = 0;
708 
709  MyBase::mm[ 8] = 0;
710  MyBase::mm[ 9] = 0;
711  MyBase::mm[10] = 1;
712  MyBase::mm[11] = 0;
713 
714  MyBase::mm[12] = v.x();
715  MyBase::mm[13] = v.y();
716  MyBase::mm[14] = v.z();
717  MyBase::mm[15] = 1;
718  }
719 
721  template <typename T0>
722  void preTranslate(const Vec3<T0>& tr)
723  {
724  Vec3<T> tmp(tr.x(), tr.y(), tr.z());
725  Mat4<T> Tr = Mat4<T>::translation(tmp);
726 
727  *this = Tr * (*this);
728 
729  }
730 
732  template <typename T0>
733  void postTranslate(const Vec3<T0>& tr)
734  {
735  Vec3<T> tmp(tr.x(), tr.y(), tr.z());
736  Mat4<T> Tr = Mat4<T>::translation(tmp);
737 
738  *this = (*this) * Tr;
739 
740  }
741 
742 
744  template <typename T0>
745  void setToScale(const Vec3<T0>& v)
746  {
747  this->setIdentity();
748  MyBase::mm[ 0] = v.x();
749  MyBase::mm[ 5] = v.y();
750  MyBase::mm[10] = v.z();
751  }
752 
753  // Left multiples by the specified scale matrix, i.e. Sc * (*this)
754  template <typename T0>
755  void preScale(const Vec3<T0>& v)
756  {
757  MyBase::mm[ 0] *= v.x();
758  MyBase::mm[ 1] *= v.x();
759  MyBase::mm[ 2] *= v.x();
760  MyBase::mm[ 3] *= v.x();
761 
762  MyBase::mm[ 4] *= v.y();
763  MyBase::mm[ 5] *= v.y();
764  MyBase::mm[ 6] *= v.y();
765  MyBase::mm[ 7] *= v.y();
766 
767  MyBase::mm[ 8] *= v.z();
768  MyBase::mm[ 9] *= v.z();
769  MyBase::mm[10] *= v.z();
770  MyBase::mm[11] *= v.z();
771  }
772 
773 
774 
775  // Right multiples by the specified scale matrix, i.e. (*this) * Sc
776  template <typename T0>
777  void postScale(const Vec3<T0>& v)
778  {
779 
780  MyBase::mm[ 0] *= v.x();
781  MyBase::mm[ 1] *= v.y();
782  MyBase::mm[ 2] *= v.z();
783 
784  MyBase::mm[ 4] *= v.x();
785  MyBase::mm[ 5] *= v.y();
786  MyBase::mm[ 6] *= v.z();
787 
788  MyBase::mm[ 8] *= v.x();
789  MyBase::mm[ 9] *= v.y();
790  MyBase::mm[10] *= v.z();
791 
792  MyBase::mm[12] *= v.x();
793  MyBase::mm[13] *= v.y();
794  MyBase::mm[14] *= v.z();
795 
796  }
797 
798 
802  void setToRotation(Axis axis, T angle) {*this = rotation<Mat4<T> >(axis, angle);}
803 
807  void setToRotation(const Vec3<T>& axis, T angle) {*this = rotation<Mat4<T> >(axis, angle);}
808 
811  void setToRotation(const Vec3<T>& v1, const Vec3<T>& v2) {*this = rotation<Mat4<T> >(v1, v2);}
812 
813 
817  void preRotate(Axis axis, T angle)
818  {
819  T c = static_cast<T>(cos(angle));
820  T s = -static_cast<T>(sin(angle)); // the "-" makes it clockwise
821 
822  switch (axis) {
823  case X_AXIS:
824  {
825  T a4, a5, a6, a7;
826 
827  a4 = c * MyBase::mm[ 4] - s * MyBase::mm[ 8];
828  a5 = c * MyBase::mm[ 5] - s * MyBase::mm[ 9];
829  a6 = c * MyBase::mm[ 6] - s * MyBase::mm[10];
830  a7 = c * MyBase::mm[ 7] - s * MyBase::mm[11];
831 
832 
833  MyBase::mm[ 8] = s * MyBase::mm[ 4] + c * MyBase::mm[ 8];
834  MyBase::mm[ 9] = s * MyBase::mm[ 5] + c * MyBase::mm[ 9];
835  MyBase::mm[10] = s * MyBase::mm[ 6] + c * MyBase::mm[10];
836  MyBase::mm[11] = s * MyBase::mm[ 7] + c * MyBase::mm[11];
837 
838  MyBase::mm[ 4] = a4;
839  MyBase::mm[ 5] = a5;
840  MyBase::mm[ 6] = a6;
841  MyBase::mm[ 7] = a7;
842  }
843  break;
844 
845  case Y_AXIS:
846  {
847  T a0, a1, a2, a3;
848 
849  a0 = c * MyBase::mm[ 0] + s * MyBase::mm[ 8];
850  a1 = c * MyBase::mm[ 1] + s * MyBase::mm[ 9];
851  a2 = c * MyBase::mm[ 2] + s * MyBase::mm[10];
852  a3 = c * MyBase::mm[ 3] + s * MyBase::mm[11];
853 
854  MyBase::mm[ 8] = -s * MyBase::mm[ 0] + c * MyBase::mm[ 8];
855  MyBase::mm[ 9] = -s * MyBase::mm[ 1] + c * MyBase::mm[ 9];
856  MyBase::mm[10] = -s * MyBase::mm[ 2] + c * MyBase::mm[10];
857  MyBase::mm[11] = -s * MyBase::mm[ 3] + c * MyBase::mm[11];
858 
859 
860  MyBase::mm[ 0] = a0;
861  MyBase::mm[ 1] = a1;
862  MyBase::mm[ 2] = a2;
863  MyBase::mm[ 3] = a3;
864  }
865  break;
866 
867  case Z_AXIS:
868  {
869  T a0, a1, a2, a3;
870 
871  a0 = c * MyBase::mm[ 0] - s * MyBase::mm[ 4];
872  a1 = c * MyBase::mm[ 1] - s * MyBase::mm[ 5];
873  a2 = c * MyBase::mm[ 2] - s * MyBase::mm[ 6];
874  a3 = c * MyBase::mm[ 3] - s * MyBase::mm[ 7];
875 
876  MyBase::mm[ 4] = s * MyBase::mm[ 0] + c * MyBase::mm[ 4];
877  MyBase::mm[ 5] = s * MyBase::mm[ 1] + c * MyBase::mm[ 5];
878  MyBase::mm[ 6] = s * MyBase::mm[ 2] + c * MyBase::mm[ 6];
879  MyBase::mm[ 7] = s * MyBase::mm[ 3] + c * MyBase::mm[ 7];
880 
881  MyBase::mm[ 0] = a0;
882  MyBase::mm[ 1] = a1;
883  MyBase::mm[ 2] = a2;
884  MyBase::mm[ 3] = a3;
885  }
886  break;
887 
888  default:
889  assert(axis==X_AXIS || axis==Y_AXIS || axis==Z_AXIS);
890  }
891  }
892 
893 
897  void postRotate(Axis axis, T angle)
898  {
899  T c = static_cast<T>(cos(angle));
900  T s = -static_cast<T>(sin(angle)); // the "-" makes it clockwise
901 
902 
903 
904  switch (axis) {
905  case X_AXIS:
906  {
907  T a2, a6, a10, a14;
908 
909  a2 = c * MyBase::mm[ 2] - s * MyBase::mm[ 1];
910  a6 = c * MyBase::mm[ 6] - s * MyBase::mm[ 5];
911  a10 = c * MyBase::mm[10] - s * MyBase::mm[ 9];
912  a14 = c * MyBase::mm[14] - s * MyBase::mm[13];
913 
914 
915  MyBase::mm[ 1] = c * MyBase::mm[ 1] + s * MyBase::mm[ 2];
916  MyBase::mm[ 5] = c * MyBase::mm[ 5] + s * MyBase::mm[ 6];
917  MyBase::mm[ 9] = c * MyBase::mm[ 9] + s * MyBase::mm[10];
918  MyBase::mm[13] = c * MyBase::mm[13] + s * MyBase::mm[14];
919 
920  MyBase::mm[ 2] = a2;
921  MyBase::mm[ 6] = a6;
922  MyBase::mm[10] = a10;
923  MyBase::mm[14] = a14;
924  }
925  break;
926 
927  case Y_AXIS:
928  {
929  T a2, a6, a10, a14;
930 
931  a2 = c * MyBase::mm[ 2] + s * MyBase::mm[ 0];
932  a6 = c * MyBase::mm[ 6] + s * MyBase::mm[ 4];
933  a10 = c * MyBase::mm[10] + s * MyBase::mm[ 8];
934  a14 = c * MyBase::mm[14] + s * MyBase::mm[12];
935 
936  MyBase::mm[ 0] = c * MyBase::mm[ 0] - s * MyBase::mm[ 2];
937  MyBase::mm[ 4] = c * MyBase::mm[ 4] - s * MyBase::mm[ 6];
938  MyBase::mm[ 8] = c * MyBase::mm[ 8] - s * MyBase::mm[10];
939  MyBase::mm[12] = c * MyBase::mm[12] - s * MyBase::mm[14];
940 
941  MyBase::mm[ 2] = a2;
942  MyBase::mm[ 6] = a6;
943  MyBase::mm[10] = a10;
944  MyBase::mm[14] = a14;
945  }
946  break;
947 
948  case Z_AXIS:
949  {
950  T a1, a5, a9, a13;
951 
952  a1 = c * MyBase::mm[ 1] - s * MyBase::mm[ 0];
953  a5 = c * MyBase::mm[ 5] - s * MyBase::mm[ 4];
954  a9 = c * MyBase::mm[ 9] - s * MyBase::mm[ 8];
955  a13 = c * MyBase::mm[13] - s * MyBase::mm[12];
956 
957  MyBase::mm[ 0] = c * MyBase::mm[ 0] + s * MyBase::mm[ 1];
958  MyBase::mm[ 4] = c * MyBase::mm[ 4] + s * MyBase::mm[ 5];
959  MyBase::mm[ 8] = c * MyBase::mm[ 8] + s * MyBase::mm[ 9];
960  MyBase::mm[12] = c * MyBase::mm[12] + s * MyBase::mm[13];
961 
962  MyBase::mm[ 1] = a1;
963  MyBase::mm[ 5] = a5;
964  MyBase::mm[ 9] = a9;
965  MyBase::mm[13] = a13;
966 
967  }
968  break;
969 
970  default:
971  assert(axis==X_AXIS || axis==Y_AXIS || axis==Z_AXIS);
972  }
973  }
974 
979  void setToShear(Axis axis0, Axis axis1, T shearby)
980  {
981  *this = shear<Mat4<T> >(axis0, axis1, shearby);
982  }
983 
984 
987  void preShear(Axis axis0, Axis axis1, T shear)
988  {
989  int index0 = static_cast<int>(axis0);
990  int index1 = static_cast<int>(axis1);
991 
992  // to row "index1" add a multiple of the index0 row
993  MyBase::mm[index1 * 4 + 0] += shear * MyBase::mm[index0 * 4 + 0];
994  MyBase::mm[index1 * 4 + 1] += shear * MyBase::mm[index0 * 4 + 1];
995  MyBase::mm[index1 * 4 + 2] += shear * MyBase::mm[index0 * 4 + 2];
996  MyBase::mm[index1 * 4 + 3] += shear * MyBase::mm[index0 * 4 + 3];
997  }
998 
999 
1002  void postShear(Axis axis0, Axis axis1, T shear)
1003  {
1004  int index0 = static_cast<int>(axis0);
1005  int index1 = static_cast<int>(axis1);
1006 
1007  // to collumn "index0" add a multiple of the index1 row
1008  MyBase::mm[index0 + 0] += shear * MyBase::mm[index1 + 0];
1009  MyBase::mm[index0 + 4] += shear * MyBase::mm[index1 + 4];
1010  MyBase::mm[index0 + 8] += shear * MyBase::mm[index1 + 8];
1011  MyBase::mm[index0 + 12] += shear * MyBase::mm[index1 + 12];
1012 
1013  }
1014 
1016  template<typename T0>
1017  Vec4<T0> transform(const Vec4<T0> &v) const
1018  {
1019  return static_cast< Vec4<T0> >(v * *this);
1020  }
1021 
1023  template<typename T0>
1024  Vec3<T0> transform(const Vec3<T0> &v) const
1025  {
1026  return static_cast< Vec3<T0> >(v * *this);
1027  }
1028 
1030  template<typename T0>
1032  {
1033  return static_cast< Vec4<T0> >(*this * v);
1034  }
1035 
1037  template<typename T0>
1039  {
1040  return static_cast< Vec3<T0> >(*this * v);
1041  }
1042 
1044  template<typename T0>
1045  Vec3<T0> transformH(const Vec3<T0> &p) const
1046  {
1047  T0 w;
1048 
1049  // w = p * (*this).col(3);
1050  w = static_cast<T0>(p[0] * MyBase::mm[ 3] + p[1] * MyBase::mm[ 7]
1051  + p[2] * MyBase::mm[11] + MyBase::mm[15]);
1052 
1053  if ( !isExactlyEqual(w , 0.0) ) {
1054  return Vec3<T0>(static_cast<T0>((p[0] * MyBase::mm[ 0] + p[1] * MyBase::mm[ 4] +
1055  p[2] * MyBase::mm[ 8] + MyBase::mm[12]) / w),
1056  static_cast<T0>((p[0] * MyBase::mm[ 1] + p[1] * MyBase::mm[ 5] +
1057  p[2] * MyBase::mm[ 9] + MyBase::mm[13]) / w),
1058  static_cast<T0>((p[0] * MyBase::mm[ 2] + p[1] * MyBase::mm[ 6] +
1059  p[2] * MyBase::mm[10] + MyBase::mm[14]) / w));
1060  }
1061 
1062  return Vec3<T0>(0, 0, 0);
1063  }
1064 
1066  template<typename T0>
1068  {
1069  T0 w;
1070 
1071  // w = p * (*this).col(3);
1072  w = p[0] * MyBase::mm[12] + p[1] * MyBase::mm[13] + p[2] * MyBase::mm[14] + MyBase::mm[15];
1073 
1074  if ( !isExactlyEqual(w , 0.0) ) {
1075  return Vec3<T0>(static_cast<T0>((p[0] * MyBase::mm[ 0] + p[1] * MyBase::mm[ 1] +
1076  p[2] * MyBase::mm[ 2] + MyBase::mm[ 3]) / w),
1077  static_cast<T0>((p[0] * MyBase::mm[ 4] + p[1] * MyBase::mm[ 5] +
1078  p[2] * MyBase::mm[ 6] + MyBase::mm[ 7]) / w),
1079  static_cast<T0>((p[0] * MyBase::mm[ 8] + p[1] * MyBase::mm[ 9] +
1080  p[2] * MyBase::mm[10] + MyBase::mm[11]) / w));
1081  }
1082 
1083  return Vec3<T0>(0, 0, 0);
1084  }
1085 
1087  template<typename T0>
1089  {
1090  return Vec3<T0>(
1091  static_cast<T0>(v[0] * MyBase::mm[ 0] + v[1] * MyBase::mm[ 4] + v[2] * MyBase::mm[ 8]),
1092  static_cast<T0>(v[0] * MyBase::mm[ 1] + v[1] * MyBase::mm[ 5] + v[2] * MyBase::mm[ 9]),
1093  static_cast<T0>(v[0] * MyBase::mm[ 2] + v[1] * MyBase::mm[ 6] + v[2] * MyBase::mm[10]));
1094  }
1095 
1096 
1097 private:
1098  bool invert(Mat4<T> &inverse, T tolerance) const;
1099 
1100  T det2(const Mat4<T> &a, int i0, int i1, int j0, int j1) const {
1101  int i0row = i0 * 4;
1102  int i1row = i1 * 4;
1103  return a.mm[i0row+j0]*a.mm[i1row+j1] - a.mm[i0row+j1]*a.mm[i1row+j0];
1104  }
1105 
1106  T det3(const Mat4<T> &a, int i0, int i1, int i2,
1107  int j0, int j1, int j2) const {
1108  int i0row = i0 * 4;
1109  return a.mm[i0row+j0]*det2(a, i1,i2, j1,j2) +
1110  a.mm[i0row+j1]*det2(a, i1,i2, j2,j0) +
1111  a.mm[i0row+j2]*det2(a, i1,i2, j0,j1);
1112  }
1113 }; // class Mat4
1114 
1115 
1118 template <typename T0, typename T1>
1119 bool operator==(const Mat4<T0> &m0, const Mat4<T1> &m1)
1120 {
1121  const T0 *t0 = m0.asPointer();
1122  const T1 *t1 = m1.asPointer();
1123 
1124  for (int i=0; i<16; ++i) if (!isExactlyEqual(t0[i], t1[i])) return false;
1125  return true;
1126 }
1127 
1130 template <typename T0, typename T1>
1131 bool operator!=(const Mat4<T0> &m0, const Mat4<T1> &m1) { return !(m0 == m1); }
1132 
1135 template <typename S, typename T>
1137 {
1138  return m*scalar;
1139 }
1140 
1143 template <typename S, typename T>
1145 {
1147  result *= scalar;
1148  return result;
1149 }
1150 
1153 template<typename T, typename MT>
1156  const Vec4<T> &_v)
1157 {
1158  MT const *m = _m.asPointer();
1160  _v[0]*m[0] + _v[1]*m[1] + _v[2]*m[2] + _v[3]*m[3],
1161  _v[0]*m[4] + _v[1]*m[5] + _v[2]*m[6] + _v[3]*m[7],
1162  _v[0]*m[8] + _v[1]*m[9] + _v[2]*m[10] + _v[3]*m[11],
1163  _v[0]*m[12] + _v[1]*m[13] + _v[2]*m[14] + _v[3]*m[15]);
1164 }
1165 
1168 template<typename T, typename MT>
1170 operator*(const Vec4<T> &_v,
1171  const Mat4<MT> &_m)
1172 {
1173  MT const *m = _m.asPointer();
1175  _v[0]*m[0] + _v[1]*m[4] + _v[2]*m[8] + _v[3]*m[12],
1176  _v[0]*m[1] + _v[1]*m[5] + _v[2]*m[9] + _v[3]*m[13],
1177  _v[0]*m[2] + _v[1]*m[6] + _v[2]*m[10] + _v[3]*m[14],
1178  _v[0]*m[3] + _v[1]*m[7] + _v[2]*m[11] + _v[3]*m[15]);
1179 }
1180 
1183 template<typename T, typename MT>
1185 operator*(const Mat4<MT> &_m, const Vec3<T> &_v)
1186 {
1187  MT const *m = _m.asPointer();
1189  _v[0]*m[0] + _v[1]*m[1] + _v[2]*m[2] + m[3],
1190  _v[0]*m[4] + _v[1]*m[5] + _v[2]*m[6] + m[7],
1191  _v[0]*m[8] + _v[1]*m[9] + _v[2]*m[10] + m[11]);
1192 }
1193 
1196 template<typename T, typename MT>
1198 operator*(const Vec3<T> &_v, const Mat4<MT> &_m)
1199 {
1200  MT const *m = _m.asPointer();
1202  _v[0]*m[0] + _v[1]*m[4] + _v[2]*m[8] + m[12],
1203  _v[0]*m[1] + _v[1]*m[5] + _v[2]*m[9] + m[13],
1204  _v[0]*m[2] + _v[1]*m[6] + _v[2]*m[10] + m[14]);
1205 }
1206 
1209 template <typename T0, typename T1>
1211 operator+(const Mat4<T0> &m0, const Mat4<T1> &m1)
1212 {
1214  result += m1;
1215  return result;
1216 }
1217 
1220 template <typename T0, typename T1>
1222 operator-(const Mat4<T0> &m0, const Mat4<T1> &m1)
1223 {
1225  result -= m1;
1226  return result;
1227 }
1228 
1231 template <typename T0, typename T1>
1233 operator*(const Mat4<T0> &m0, const Mat4<T1> &m1)
1234 {
1236  result *= m1;
1237  return result;
1238 }
1239 
1240 
1244 template<typename T0, typename T1>
1246 {
1247  return Vec3<T1>(
1248  static_cast<T1>(m[0][0]*n[0] + m[0][1]*n[1] + m[0][2]*n[2]),
1249  static_cast<T1>(m[1][0]*n[0] + m[1][1]*n[1] + m[1][2]*n[2]),
1250  static_cast<T1>(m[2][0]*n[0] + m[2][1]*n[1] + m[2][2]*n[2]));
1251 }
1252 
1253 
1255 template<typename T>
1256 bool Mat4<T>::invert(Mat4<T> &inverse, T tolerance) const
1257 {
1258  Mat4<T> temp(*this);
1259  inverse.setIdentity();
1260 
1261  // Forward elimination step
1262  double det = 1.0;
1263  for (int i = 0; i < 4; ++i) {
1264  int row = i;
1265  double max = fabs(temp[i][i]);
1266 
1267  for (int k = i+1; k < 4; ++k) {
1268  if (fabs(temp[k][i]) > max) {
1269  row = k;
1270  max = fabs(temp[k][i]);
1271  }
1272  }
1273 
1274  if (isExactlyEqual(max, 0.0)) return false;
1275 
1276  // must move pivot to row i
1277  if (row != i) {
1278  det = -det;
1279  for (int k = 0; k < 4; ++k) {
1280  std::swap(temp[row][k], temp[i][k]);
1281  std::swap(inverse[row][k], inverse[i][k]);
1282  }
1283  }
1284 
1285  double pivot = temp[i][i];
1286  det *= pivot;
1287 
1288  // scale row i
1289  for (int k = 0; k < 4; ++k) {
1290  temp[i][k] /= pivot;
1291  inverse[i][k] /= pivot;
1292  }
1293 
1294  // eliminate in rows below i
1295  for (int j = i+1; j < 4; ++j) {
1296  double t = temp[j][i];
1297  if (!isExactlyEqual(t, 0.0)) {
1298  // subtract scaled row i from row j
1299  for (int k = 0; k < 4; ++k) {
1300  temp[j][k] -= temp[i][k] * t;
1301  inverse[j][k] -= inverse[i][k] * t;
1302  }
1303  }
1304  }
1305  }
1306 
1307  // Backward elimination step
1308  for (int i = 3; i > 0; --i) {
1309  for (int j = 0; j < i; ++j) {
1310  double t = temp[j][i];
1311 
1312  if (!isExactlyEqual(t, 0.0)) {
1313  for (int k = 0; k < 4; ++k) {
1314  inverse[j][k] -= inverse[i][k]*t;
1315  }
1316  }
1317  }
1318  }
1319  return det*det >= tolerance*tolerance;
1320 }
1321 
1322 template <typename T>
1323 inline bool isAffine(const Mat4<T>& m) {
1324  return (m.col(3) == Vec4<T>(0, 0, 0, 1));
1325 }
1326 
1327 template <typename T>
1328 inline bool hasTranslation(const Mat4<T>& m) {
1329  return (m.row(3) != Vec4<T>(0, 0, 0, 1));
1330 }
1331 
1332 
1335 using Mat4f = Mat4d;
1336 
1337 } // namespace math
1338 
1339 
1340 template<> inline math::Mat4s zeroVal<math::Mat4s>() { return math::Mat4s::zero(); }
1341 template<> inline math::Mat4d zeroVal<math::Mat4d>() { return math::Mat4d::zero(); }
1342 
1343 } // namespace OPENVDB_VERSION_NAME
1344 } // namespace openvdb
1345 
1346 #endif // OPENVDB_UTIL_MAT4_H_HAS_BEEN_INCLUDED
bool hasTranslation(const Mat4< T > &m)
Definition: Mat4.h:1328
T det() const
Determinant of matrix.
Definition: Mat4.h:656
bool operator!=(const Mat4< T0 > &m0, const Mat4< T1 > &m1)
Inequality operator, does exact floating point comparisons.
Definition: Mat4.h:1131
void preScale(const Vec3< T0 > &v)
Definition: Mat4.h:755
Vec3< T0 > pretransformH(const Vec3< T0 > &p) const
Transform a Vec3 by pre-multiplication, doing homogenous division.
Definition: Mat4.h:1067
void setRows(const Vec4< T > &v1, const Vec4< T > &v2, const Vec4< T > &v3, const Vec4< T > &v4)
Set the rows of this matrix to the vectors v1, v2, v3, v4.
Definition: Mat4.h:213
Mat4(const Mat4< Source > &m)
Conversion constructor.
Definition: Mat4.h:115
void pivot(int i, int j, Mat3< T > &S, Vec3< T > &D, Mat3< T > &Q)
Definition: Mat3.h:689
T & z()
Definition: Vec3.h:85
static Mat4 translation(const Vec3d &v)
Sets the matrix to a matrix that translates by v.
Definition: Mat4.h:686
bool isApproxEqual(const Type &a, const Type &b)
Return true if a is equal to b to within the default floating-point comparison tolerance.
Definition: Math.h:351
Mat4< typename promote< S, T >::type > operator*(S scalar, const Mat4< T > &m)
Multiply each element of the given matrix by scalar and return the result.
Definition: Mat4.h:1136
Mat4< double > Mat4d
Definition: Mat4.h:1334
General-purpose arithmetic and comparison routines, most of which accept arbitrary value types (or at...
static const Mat4< T > & zero()
Predefined constant for zero matrix.
Definition: Mat4.h:136
Vec3< typename promote< T, MT >::type > operator*(const Vec3< T > &_v, const Mat4< MT > &_m)
Multiply _v by _m and return the resulting vector.
Definition: Mat4.h:1198
Real ValueType
Definition: Mat.h:30
const T * operator[](int i) const
Definition: Mat4.h:185
void preShear(Axis axis0, Axis axis1, T shear)
Left multiplies a shearing transformation into the matrix.
Definition: Mat4.h:987
#define OPENVDB_THROW(exception, message)
Definition: Exceptions.h:82
static const Mat4< T > & identity()
Predefined constant for identity matrix.
Definition: Mat4.h:125
Definition: Mat.h:26
Vec4< T > row(int i) const
Get ith row, e.g. Vec4f v = m.row(1);.
Definition: Mat4.h:158
void setMat3(const Mat3< T > &m)
Set upper left to a Mat3.
Definition: Mat4.h:309
Vec4< T0 > pretransform(const Vec4< T0 > &v) const
Transform a Vec4 by pre-multiplication.
Definition: Mat4.h:1031
void postScale(const Vec3< T0 > &v)
Definition: Mat4.h:777
bool operator==(const Mat4< T0 > &m0, const Mat4< T1 > &m1)
Equality operator, does exact floating point comparisons.
Definition: Mat4.h:1119
void setToShear(Axis axis0, Axis axis1, T shearby)
Sets the matrix to a shear along axis0 by a fraction of axis1.
Definition: Mat4.h:979
void preRotate(Axis axis, T angle)
Left multiplies by a rotation clock-wiseabout the given axis into this matrix.
Definition: Mat4.h:817
Definition: Math.h:851
T mm[SIZE *SIZE]
Definition: Mat.h:165
void setCol(int j, const Vec4< T > &v)
Set jth column to vector v.
Definition: Mat4.h:165
void setToRotation(const Vec3< T > &v1, const Vec3< T > &v2)
Sets the matrix to a rotation that maps v1 onto v2 about the cross product of v1 and v2...
Definition: Mat4.h:811
T & y()
Definition: Vec3.h:84
T & operator()(int i, int j)
Definition: Mat4.h:195
Vec3< T1 > transformNormal(const Mat4< T0 > &m, const Vec3< T1 > &n)
Definition: Mat4.h:1245
Mat4(const Mat< 4, T > &m)
Copy constructor.
Definition: Mat4.h:104
Axis
Definition: Math.h:849
Vec4< typename promote< T, MT >::type > operator*(const Vec4< T > &_v, const Mat4< MT > &_m)
Multiply _v by _m and return the resulting vector.
Definition: Mat4.h:1170
Mat3< T > getMat3() const
Definition: Mat4.h:316
void setZero()
Definition: Mat4.h:263
Vec4< typename promote< T, MT >::type > operator*(const Mat4< MT > &_m, const Vec4< T > &_v)
Multiply _m by _v and return the resulting vector.
Definition: Mat4.h:1155
Mat4< typename promote< T0, T1 >::type > operator-(const Mat4< T0 > &m0, const Mat4< T1 > &m1)
Subtract corresponding elements of m0 and m1 and return the result.
Definition: Mat4.h:1222
Mat4(Source a, Source b, Source c, Source d, Source e, Source f, Source g, Source h, Source i, Source j, Source k, Source l, Source m, Source n, Source o, Source p)
Constructor given array of elements, the ordering is in row major form:
Definition: Mat4.h:64
Definition: Math.h:850
void setRow(int i, const Vec4< T > &v)
Set ith row to vector v.
Definition: Mat4.h:147
void setColumns(const Vec4< T > &v1, const Vec4< T > &v2, const Vec4< T > &v3, const Vec4< T > &v4)
Set the columns of this matrix to the vectors v1, v2, v3, v4.
Definition: Mat4.h:238
Vec4< T0 > transform(const Vec4< T0 > &v) const
Transform a Vec4 by post-multiplication.
Definition: Mat4.h:1017
T dot(const Vec3< T > &v) const
Dot product.
Definition: Vec3.h:189
#define OPENVDB_VERSION_NAME
The version namespace name for this library version.
Definition: version.h:102
Vec3< typename promote< T, MT >::type > operator*(const Mat4< MT > &_m, const Vec3< T > &_v)
Multiply _m by _v and return the resulting vector.
Definition: Mat4.h:1185
Mat4(Source *a)
Constructor given array of elements, the ordering is in row major form:
Definition: Mat4.h:49
Mat4< T > operator-() const
Negation operator, for e.g. m1 = -m2;.
Definition: Mat4.h:362
const Mat4< T > & operator*=(const Mat4< S > &m1)
Multiply this matrix by the given matrix.
Definition: Mat4.h:458
T * asPointer()
Definition: Mat3.h:192
void postTranslate(const Vec3< T0 > &tr)
Right multiplies by the specified translation matrix, i.e. (*this) * Trans.
Definition: Mat4.h:733
4x4 -matrix class.
Definition: Mat3.h:22
const Mat4< T > & operator+=(const Mat4< S > &m1)
Add each element of the given matrix to the corresponding element of this matrix. ...
Definition: Mat4.h:400
Mat4 transpose() const
Definition: Mat4.h:491
Definition: Mat4.h:24
Definition: Exceptions.h:13
Definition: Exceptions.h:56
const T * asPointer() const
Definition: Mat4.h:190
Vec4< T > col(int j) const
Get jth column, e.g. Vec4f v = m.col(0);.
Definition: Mat4.h:175
Vec3< T0 > transformH(const Vec3< T0 > &p) const
Transform a Vec3 by post-multiplication, doing homogenous divison.
Definition: Mat4.h:1045
void setToRotation(Axis axis, T angle)
Sets the matrix to a rotation about the given axis.
Definition: Mat4.h:802
void setToRotation(const Vec3< T > &axis, T angle)
Sets the matrix to a rotation about an arbitrary axis.
Definition: Mat4.h:807
Vec3< T0 > transform3x3(const Vec3< T0 > &v) const
Transform a Vec3 by post-multiplication, without translation.
Definition: Mat4.h:1088
T det() const
Determinant of matrix.
Definition: Mat3.h:500
bool eq(const Mat4 &m, T eps=1.0e-8) const
Return true if this matrix is equivalent to m within a tolerance of eps.
Definition: Mat4.h:352
T & x()
Reference to the component, e.g. v.x() = 4.5f;.
Definition: Vec3.h:83
Mat4< typename promote< T0, T1 >::type > operator+(const Mat4< T0 > &m0, const Mat4< T1 > &m1)
Add corresponding elements of m0 and m1 and return the result.
Definition: Mat4.h:1211
Mat4< typename promote< T0, T1 >::type > operator*(const Mat4< T0 > &m0, const Mat4< T1 > &m1)
Multiply m0 by m1 and return the resulting matrix.
Definition: Mat4.h:1233
Vec3< T0 > transform(const Vec3< T0 > &v) const
Transform a Vec3 by post-multiplication, without homogenous division.
Definition: Mat4.h:1024
T * asPointer()
Direct access to the internal data.
Definition: Mat4.h:189
void setToTranslation(const Vec3< T0 > &v)
Sets the matrix to a matrix that translates by v.
Definition: Mat4.h:697
3x3 matrix class.
Definition: Mat3.h:28
void preTranslate(const Vec3< T0 > &tr)
Left multiples by the specified translation, i.e. Trans * (*this)
Definition: Mat4.h:722
Definition: Mat.h:170
void postRotate(Axis axis, T angle)
Right multiplies by a rotation clock-wiseabout the given axis into this matrix.
Definition: Mat4.h:897
Vec3< T0 > pretransform(const Vec3< T0 > &v) const
Transform a Vec3 by pre-multiplication, without homogenous division.
Definition: Mat4.h:1038
Real value_type
Definition: Mat.h:29
MatType shear(Axis axis0, Axis axis1, typename MatType::value_type shear)
Set the matrix to a shear along axis0 by a fraction of axis1.
Definition: Mat.h:693
T angle(const Vec2< T > &v1, const Vec2< T > &v2)
Definition: Vec2.h:445
Vec3< T > getTranslation() const
Return the translation component.
Definition: Mat4.h:328
bool isExactlyEqual(const T0 &a, const T1 &b)
Return true if a is exactly equal to b.
Definition: Math.h:388
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h:154
const Mat4< T > & operator*=(S scalar)
Multiply each element of this matrix by scalar.
Definition: Mat4.h:374
const Mat4< T > & operator-=(const Mat4< S > &m1)
Subtract each element of the given matrix from the corresponding element of this matrix.
Definition: Mat4.h:429
Mat4 inverse(T tolerance=0) const
Definition: Mat4.h:504
void postShear(Axis axis0, Axis axis1, T shear)
Right multiplies a shearing transformation into the matrix.
Definition: Mat4.h:1002
void setIdentity()
Set this matrix to identity.
Definition: Mat4.h:284
T operator()(int i, int j) const
Definition: Mat4.h:205
Mat4< typename promote< S, T >::type > operator*(const Mat4< T > &m, S scalar)
Multiply each element of the given matrix by scalar and return the result.
Definition: Mat4.h:1144
Mat4(const Vec4< Source > &v1, const Vec4< Source > &v2, const Vec4< Source > &v3, const Vec4< Source > &v4, bool rows=true)
Definition: Mat4.h:93
Definition: Math.h:852
void setToScale(const Vec3< T0 > &v)
Sets the matrix to a matrix that scales by v.
Definition: Mat4.h:745
const Mat4 & operator=(const Mat4< Source > &m)
Assignment operator.
Definition: Mat4.h:342
Mat4()
Trivial constructor, the matrix is NOT initialized.
Definition: Mat4.h:39
bool isAffine(const Mat4< T > &m)
Definition: Mat4.h:1323
const std::enable_if<!VecTraits< T >::IsVec, T >::type & max(const T &a, const T &b)
Definition: Composite.h:106
void setTranslation(const Vec3< T > &t)
Definition: Mat4.h:333