4 #ifndef OPENVDB_MATH_MAT4_H_HAS_BEEN_INCLUDED
5 #define OPENVDB_MATH_MAT4_H_HAS_BEEN_INCLUDED
24 template<
typename T>
class Vec4;
30 class Mat4:
public Mat<4, T>
39 #if OPENVDB_ABI_VERSION_NUMBER >= 8
49 for (
int i = 0; i < 4; ++i) {
50 for (
int j = 0; j < 4; ++j) {
51 MyBase::mm[i*4 + j] = m[i][j];
64 template<
typename Source>
67 for (
int i = 0; i < 16; i++) {
68 MyBase::mm[i] =
static_cast<T
>(a[i]);
79 template<
typename Source>
80 Mat4(Source a, Source b, Source c, Source d,
81 Source e, Source f, Source g, Source h,
82 Source i, Source j, Source k, Source l,
83 Source m, Source n, Source o, Source p)
85 MyBase::mm[ 0] =
static_cast<T
>(a);
86 MyBase::mm[ 1] =
static_cast<T
>(b);
87 MyBase::mm[ 2] =
static_cast<T
>(c);
88 MyBase::mm[ 3] =
static_cast<T
>(d);
90 MyBase::mm[ 4] =
static_cast<T
>(e);
91 MyBase::mm[ 5] =
static_cast<T
>(f);
92 MyBase::mm[ 6] =
static_cast<T
>(g);
93 MyBase::mm[ 7] =
static_cast<T
>(h);
95 MyBase::mm[ 8] =
static_cast<T
>(i);
96 MyBase::mm[ 9] =
static_cast<T
>(j);
97 MyBase::mm[10] =
static_cast<T
>(k);
98 MyBase::mm[11] =
static_cast<T
>(l);
100 MyBase::mm[12] =
static_cast<T
>(m);
101 MyBase::mm[13] =
static_cast<T
>(n);
102 MyBase::mm[14] =
static_cast<T
>(o);
103 MyBase::mm[15] =
static_cast<T
>(p);
108 template<
typename Source>
113 this->setRows(v1, v2, v3, v4);
115 this->setColumns(v1, v2, v3, v4);
120 template<
typename Source>
123 const Source *src = m.asPointer();
125 for (
int i=0; i<16; ++i) {
126 MyBase::mm[i] =
static_cast<T
>(src[i]);
157 MyBase::mm[i4+0] = v[0];
158 MyBase::mm[i4+1] = v[1];
159 MyBase::mm[i4+2] = v[2];
160 MyBase::mm[i4+3] = v[3];
167 return Vec4<T>((*
this)(i,0), (*
this)(i,1), (*
this)(i,2), (*
this)(i,3));
174 MyBase::mm[ 0+j] = v[0];
175 MyBase::mm[ 4+j] = v[1];
176 MyBase::mm[ 8+j] = v[2];
177 MyBase::mm[12+j] = v[3];
184 return Vec4<T>((*
this)(0,j), (*
this)(1,j), (*
this)(2,j), (*
this)(3,j));
194 return MyBase::mm[4*i+j];
204 return MyBase::mm[4*i+j];
211 MyBase::mm[ 0] = v1[0];
212 MyBase::mm[ 1] = v1[1];
213 MyBase::mm[ 2] = v1[2];
214 MyBase::mm[ 3] = v1[3];
216 MyBase::mm[ 4] = v2[0];
217 MyBase::mm[ 5] = v2[1];
218 MyBase::mm[ 6] = v2[2];
219 MyBase::mm[ 7] = v2[3];
221 MyBase::mm[ 8] = v3[0];
222 MyBase::mm[ 9] = v3[1];
223 MyBase::mm[10] = v3[2];
224 MyBase::mm[11] = v3[3];
226 MyBase::mm[12] = v4[0];
227 MyBase::mm[13] = v4[1];
228 MyBase::mm[14] = v4[2];
229 MyBase::mm[15] = v4[3];
236 MyBase::mm[ 0] = v1[0];
237 MyBase::mm[ 1] = v2[0];
238 MyBase::mm[ 2] = v3[0];
239 MyBase::mm[ 3] = v4[0];
241 MyBase::mm[ 4] = v1[1];
242 MyBase::mm[ 5] = v2[1];
243 MyBase::mm[ 6] = v3[1];
244 MyBase::mm[ 7] = v4[1];
246 MyBase::mm[ 8] = v1[2];
247 MyBase::mm[ 9] = v2[2];
248 MyBase::mm[10] = v3[2];
249 MyBase::mm[11] = v4[2];
251 MyBase::mm[12] = v1[3];
252 MyBase::mm[13] = v2[3];
253 MyBase::mm[14] = v3[3];
254 MyBase::mm[15] = v4[3];
306 for (
int i = 0; i < 3; i++)
307 for (
int j=0; j < 3; j++)
308 MyBase::mm[i*4+j] = m[i][j];
315 for (
int i = 0; i < 3; i++)
316 for (
int j = 0; j < 3; j++)
317 m[i][j] = MyBase::mm[i*4+j];
325 return Vec3<T>(MyBase::mm[12], MyBase::mm[13], MyBase::mm[14]);
330 MyBase::mm[12] = t[0];
331 MyBase::mm[13] = t[1];
332 MyBase::mm[14] = t[2];
336 template<
typename Source>
339 const Source *src = m.asPointer();
342 std::copy(src, (src + this->numElements()), MyBase::mm);
347 bool eq(
const Mat4 &m, T eps=1.0e-8)
const
349 for (
int i = 0; i < 16; i++) {
360 -MyBase::mm[ 0], -MyBase::mm[ 1], -MyBase::mm[ 2], -MyBase::mm[ 3],
361 -MyBase::mm[ 4], -MyBase::mm[ 5], -MyBase::mm[ 6], -MyBase::mm[ 7],
362 -MyBase::mm[ 8], -MyBase::mm[ 9], -MyBase::mm[10], -MyBase::mm[11],
363 -MyBase::mm[12], -MyBase::mm[13], -MyBase::mm[14], -MyBase::mm[15]
368 template <
typename S>
371 MyBase::mm[ 0] *= scalar;
372 MyBase::mm[ 1] *= scalar;
373 MyBase::mm[ 2] *= scalar;
374 MyBase::mm[ 3] *= scalar;
376 MyBase::mm[ 4] *= scalar;
377 MyBase::mm[ 5] *= scalar;
378 MyBase::mm[ 6] *= scalar;
379 MyBase::mm[ 7] *= scalar;
381 MyBase::mm[ 8] *= scalar;
382 MyBase::mm[ 9] *= scalar;
383 MyBase::mm[10] *= scalar;
384 MyBase::mm[11] *= scalar;
386 MyBase::mm[12] *= scalar;
387 MyBase::mm[13] *= scalar;
388 MyBase::mm[14] *= scalar;
389 MyBase::mm[15] *= scalar;
394 template <
typename S>
397 const S* s = m1.asPointer();
399 MyBase::mm[ 0] += s[ 0];
400 MyBase::mm[ 1] += s[ 1];
401 MyBase::mm[ 2] += s[ 2];
402 MyBase::mm[ 3] += s[ 3];
404 MyBase::mm[ 4] += s[ 4];
405 MyBase::mm[ 5] += s[ 5];
406 MyBase::mm[ 6] += s[ 6];
407 MyBase::mm[ 7] += s[ 7];
409 MyBase::mm[ 8] += s[ 8];
410 MyBase::mm[ 9] += s[ 9];
411 MyBase::mm[10] += s[10];
412 MyBase::mm[11] += s[11];
414 MyBase::mm[12] += s[12];
415 MyBase::mm[13] += s[13];
416 MyBase::mm[14] += s[14];
417 MyBase::mm[15] += s[15];
423 template <
typename S>
426 const S* s = m1.asPointer();
428 MyBase::mm[ 0] -= s[ 0];
429 MyBase::mm[ 1] -= s[ 1];
430 MyBase::mm[ 2] -= s[ 2];
431 MyBase::mm[ 3] -= s[ 3];
433 MyBase::mm[ 4] -= s[ 4];
434 MyBase::mm[ 5] -= s[ 5];
435 MyBase::mm[ 6] -= s[ 6];
436 MyBase::mm[ 7] -= s[ 7];
438 MyBase::mm[ 8] -= s[ 8];
439 MyBase::mm[ 9] -= s[ 9];
440 MyBase::mm[10] -= s[10];
441 MyBase::mm[11] -= s[11];
443 MyBase::mm[12] -= s[12];
444 MyBase::mm[13] -= s[13];
445 MyBase::mm[14] -= s[14];
446 MyBase::mm[15] -= s[15];
452 template <
typename S>
457 const T* s0 = m0.asPointer();
458 const S* s1 = m1.asPointer();
460 for (
int i = 0; i < 4; i++) {
462 MyBase::mm[i4+0] =
static_cast<T
>(s0[i4+0] * s1[ 0] +
467 MyBase::mm[i4+1] =
static_cast<T
>(s0[i4+0] * s1[ 1] +
472 MyBase::mm[i4+2] =
static_cast<T
>(s0[i4+0] * s1[ 2] +
477 MyBase::mm[i4+3] =
static_cast<T
>(s0[i4+0] * s1[ 3] +
489 MyBase::mm[ 0], MyBase::mm[ 4], MyBase::mm[ 8], MyBase::mm[12],
490 MyBase::mm[ 1], MyBase::mm[ 5], MyBase::mm[ 9], MyBase::mm[13],
491 MyBase::mm[ 2], MyBase::mm[ 6], MyBase::mm[10], MyBase::mm[14],
492 MyBase::mm[ 3], MyBase::mm[ 7], MyBase::mm[11], MyBase::mm[15]
522 T m0011 = m[0][0] * m[1][1];
523 T m0012 = m[0][0] * m[1][2];
524 T m0110 = m[0][1] * m[1][0];
525 T m0210 = m[0][2] * m[1][0];
526 T m0120 = m[0][1] * m[2][0];
527 T m0220 = m[0][2] * m[2][0];
529 T detA = m0011 * m[2][2] - m0012 * m[2][1] - m0110 * m[2][2]
530 + m0210 * m[2][1] + m0120 * m[1][2] - m0220 * m[1][1];
532 bool hasPerspective =
539 if (hasPerspective) {
540 det = m[0][3] * det3(m, 1,2,3, 0,2,1)
541 + m[1][3] * det3(m, 2,0,3, 0,2,1)
542 + m[2][3] * det3(m, 3,0,1, 0,2,1)
545 det = detA * m[3][3];
556 invertible = m.invert(inv, tolerance);
565 inv[0][0] = detA * ( m[1][1] * m[2][2] - m[1][2] * m[2][1]);
566 inv[0][1] = detA * (-m[0][1] * m[2][2] + m[0][2] * m[2][1]);
567 inv[0][2] = detA * ( m[0][1] * m[1][2] - m[0][2] * m[1][1]);
569 inv[1][0] = detA * (-m[1][0] * m[2][2] + m[1][2] * m[2][0]);
570 inv[1][1] = detA * ( m[0][0] * m[2][2] - m0220);
571 inv[1][2] = detA * ( m0210 - m0012);
573 inv[2][0] = detA * ( m[1][0] * m[2][1] - m[1][1] * m[2][0]);
574 inv[2][1] = detA * ( m0120 - m[0][0] * m[2][1]);
575 inv[2][2] = detA * ( m0011 - m0110);
577 if (hasPerspective) {
582 r[0] = m[3][0] * inv[0][0] + m[3][1] * inv[1][0]
583 + m[3][2] * inv[2][0];
584 r[1] = m[3][0] * inv[0][1] + m[3][1] * inv[1][1]
585 + m[3][2] * inv[2][1];
586 r[2] = m[3][0] * inv[0][2] + m[3][1] * inv[1][2]
587 + m[3][2] * inv[2][2];
590 p[0] = inv[0][0] * m[0][3] + inv[0][1] * m[1][3]
591 + inv[0][2] * m[2][3];
592 p[1] = inv[1][0] * m[0][3] + inv[1][1] * m[1][3]
593 + inv[1][2] * m[2][3];
594 p[2] = inv[2][0] * m[0][3] + inv[2][1] * m[1][3]
595 + inv[2][2] * m[2][3];
597 T h = m[3][3] - p.
dot(
Vec3<T>(m[3][0],m[3][1],m[3][2]));
608 inv[3][0] = -h * r[0];
609 inv[3][1] = -h * r[1];
610 inv[3][2] = -h * r[2];
612 inv[0][3] = -h * p[0];
613 inv[1][3] = -h * p[1];
614 inv[2][3] = -h * p[2];
620 inv[0][0] += p[0] * r[0];
621 inv[0][1] += p[0] * r[1];
622 inv[0][2] += p[0] * r[2];
623 inv[1][0] += p[1] * r[0];
624 inv[1][1] += p[1] * r[1];
625 inv[1][2] += p[1] * r[2];
626 inv[2][0] += p[2] * r[0];
627 inv[2][1] += p[2] * r[1];
628 inv[2][2] += p[2] * r[2];
632 inv[3][0] = - (m[3][0] * inv[0][0] + m[3][1] * inv[1][0]
633 + m[3][2] * inv[2][0]);
634 inv[3][1] = - (m[3][0] * inv[0][1] + m[3][1] * inv[1][1]
635 + m[3][2] * inv[2][1]);
636 inv[3][2] = - (m[3][0] * inv[0][2] + m[3][1] * inv[1][2]
637 + m[3][2] * inv[2][2]);
661 for (i = 0; i < 4; i++) {
662 ap = &MyBase::mm[ 0];
664 for (j = 0; j < 4; j++) {
665 for (k = 0; k < 4; k++) {
666 if ((k != i) && (j != 0)) {
673 det += T(sign) * MyBase::mm[i] * submat.
det();
684 T(1), T(0), T(0), T(0),
685 T(0), T(1), T(0), T(0),
686 T(0), T(0), T(1), T(0),
687 T(v.
x()), T(v.
y()),T(v.
z()), T(1));
691 template <
typename T0>
709 MyBase::mm[12] = v.
x();
710 MyBase::mm[13] = v.
y();
711 MyBase::mm[14] = v.
z();
716 template <
typename T0>
722 *
this = Tr * (*this);
727 template <
typename T0>
733 *
this = (*this) * Tr;
739 template <
typename T0>
743 MyBase::mm[ 0] = v.
x();
744 MyBase::mm[ 5] = v.
y();
745 MyBase::mm[10] = v.
z();
749 template <
typename T0>
752 MyBase::mm[ 0] *= v.
x();
753 MyBase::mm[ 1] *= v.
x();
754 MyBase::mm[ 2] *= v.
x();
755 MyBase::mm[ 3] *= v.
x();
757 MyBase::mm[ 4] *= v.
y();
758 MyBase::mm[ 5] *= v.
y();
759 MyBase::mm[ 6] *= v.
y();
760 MyBase::mm[ 7] *= v.
y();
762 MyBase::mm[ 8] *= v.
z();
763 MyBase::mm[ 9] *= v.
z();
764 MyBase::mm[10] *= v.
z();
765 MyBase::mm[11] *= v.
z();
771 template <
typename T0>
775 MyBase::mm[ 0] *= v.
x();
776 MyBase::mm[ 1] *= v.
y();
777 MyBase::mm[ 2] *= v.
z();
779 MyBase::mm[ 4] *= v.
x();
780 MyBase::mm[ 5] *= v.
y();
781 MyBase::mm[ 6] *= v.
z();
783 MyBase::mm[ 8] *= v.
x();
784 MyBase::mm[ 9] *= v.
y();
785 MyBase::mm[10] *= v.
z();
787 MyBase::mm[12] *= v.
x();
788 MyBase::mm[13] *= v.
y();
789 MyBase::mm[14] *= v.
z();
814 T c =
static_cast<T
>(cos(
angle));
815 T s = -
static_cast<T
>(sin(
angle));
822 a4 = c * MyBase::mm[ 4] - s * MyBase::mm[ 8];
823 a5 = c * MyBase::mm[ 5] - s * MyBase::mm[ 9];
824 a6 = c * MyBase::mm[ 6] - s * MyBase::mm[10];
825 a7 = c * MyBase::mm[ 7] - s * MyBase::mm[11];
828 MyBase::mm[ 8] = s * MyBase::mm[ 4] + c * MyBase::mm[ 8];
829 MyBase::mm[ 9] = s * MyBase::mm[ 5] + c * MyBase::mm[ 9];
830 MyBase::mm[10] = s * MyBase::mm[ 6] + c * MyBase::mm[10];
831 MyBase::mm[11] = s * MyBase::mm[ 7] + c * MyBase::mm[11];
844 a0 = c * MyBase::mm[ 0] + s * MyBase::mm[ 8];
845 a1 = c * MyBase::mm[ 1] + s * MyBase::mm[ 9];
846 a2 = c * MyBase::mm[ 2] + s * MyBase::mm[10];
847 a3 = c * MyBase::mm[ 3] + s * MyBase::mm[11];
849 MyBase::mm[ 8] = -s * MyBase::mm[ 0] + c * MyBase::mm[ 8];
850 MyBase::mm[ 9] = -s * MyBase::mm[ 1] + c * MyBase::mm[ 9];
851 MyBase::mm[10] = -s * MyBase::mm[ 2] + c * MyBase::mm[10];
852 MyBase::mm[11] = -s * MyBase::mm[ 3] + c * MyBase::mm[11];
866 a0 = c * MyBase::mm[ 0] - s * MyBase::mm[ 4];
867 a1 = c * MyBase::mm[ 1] - s * MyBase::mm[ 5];
868 a2 = c * MyBase::mm[ 2] - s * MyBase::mm[ 6];
869 a3 = c * MyBase::mm[ 3] - s * MyBase::mm[ 7];
871 MyBase::mm[ 4] = s * MyBase::mm[ 0] + c * MyBase::mm[ 4];
872 MyBase::mm[ 5] = s * MyBase::mm[ 1] + c * MyBase::mm[ 5];
873 MyBase::mm[ 6] = s * MyBase::mm[ 2] + c * MyBase::mm[ 6];
874 MyBase::mm[ 7] = s * MyBase::mm[ 3] + c * MyBase::mm[ 7];
894 T c =
static_cast<T
>(cos(
angle));
895 T s = -
static_cast<T
>(sin(
angle));
904 a2 = c * MyBase::mm[ 2] - s * MyBase::mm[ 1];
905 a6 = c * MyBase::mm[ 6] - s * MyBase::mm[ 5];
906 a10 = c * MyBase::mm[10] - s * MyBase::mm[ 9];
907 a14 = c * MyBase::mm[14] - s * MyBase::mm[13];
910 MyBase::mm[ 1] = c * MyBase::mm[ 1] + s * MyBase::mm[ 2];
911 MyBase::mm[ 5] = c * MyBase::mm[ 5] + s * MyBase::mm[ 6];
912 MyBase::mm[ 9] = c * MyBase::mm[ 9] + s * MyBase::mm[10];
913 MyBase::mm[13] = c * MyBase::mm[13] + s * MyBase::mm[14];
917 MyBase::mm[10] = a10;
918 MyBase::mm[14] = a14;
926 a2 = c * MyBase::mm[ 2] + s * MyBase::mm[ 0];
927 a6 = c * MyBase::mm[ 6] + s * MyBase::mm[ 4];
928 a10 = c * MyBase::mm[10] + s * MyBase::mm[ 8];
929 a14 = c * MyBase::mm[14] + s * MyBase::mm[12];
931 MyBase::mm[ 0] = c * MyBase::mm[ 0] - s * MyBase::mm[ 2];
932 MyBase::mm[ 4] = c * MyBase::mm[ 4] - s * MyBase::mm[ 6];
933 MyBase::mm[ 8] = c * MyBase::mm[ 8] - s * MyBase::mm[10];
934 MyBase::mm[12] = c * MyBase::mm[12] - s * MyBase::mm[14];
938 MyBase::mm[10] = a10;
939 MyBase::mm[14] = a14;
947 a1 = c * MyBase::mm[ 1] - s * MyBase::mm[ 0];
948 a5 = c * MyBase::mm[ 5] - s * MyBase::mm[ 4];
949 a9 = c * MyBase::mm[ 9] - s * MyBase::mm[ 8];
950 a13 = c * MyBase::mm[13] - s * MyBase::mm[12];
952 MyBase::mm[ 0] = c * MyBase::mm[ 0] + s * MyBase::mm[ 1];
953 MyBase::mm[ 4] = c * MyBase::mm[ 4] + s * MyBase::mm[ 5];
954 MyBase::mm[ 8] = c * MyBase::mm[ 8] + s * MyBase::mm[ 9];
955 MyBase::mm[12] = c * MyBase::mm[12] + s * MyBase::mm[13];
960 MyBase::mm[13] = a13;
976 *
this = shear<Mat4<T> >(axis0, axis1, shearby);
984 int index0 =
static_cast<int>(axis0);
985 int index1 =
static_cast<int>(axis1);
988 MyBase::mm[index1 * 4 + 0] +=
shear * MyBase::mm[index0 * 4 + 0];
989 MyBase::mm[index1 * 4 + 1] +=
shear * MyBase::mm[index0 * 4 + 1];
990 MyBase::mm[index1 * 4 + 2] +=
shear * MyBase::mm[index0 * 4 + 2];
991 MyBase::mm[index1 * 4 + 3] +=
shear * MyBase::mm[index0 * 4 + 3];
999 int index0 =
static_cast<int>(axis0);
1000 int index1 =
static_cast<int>(axis1);
1003 MyBase::mm[index0 + 0] +=
shear * MyBase::mm[index1 + 0];
1004 MyBase::mm[index0 + 4] +=
shear * MyBase::mm[index1 + 4];
1005 MyBase::mm[index0 + 8] +=
shear * MyBase::mm[index1 + 8];
1006 MyBase::mm[index0 + 12] +=
shear * MyBase::mm[index1 + 12];
1011 template<
typename T0>
1014 return static_cast< Vec4<T0> >(v * *
this);
1018 template<
typename T0>
1021 return static_cast< Vec3<T0> >(v * *
this);
1025 template<
typename T0>
1028 return static_cast< Vec4<T0> >(*
this * v);
1032 template<
typename T0>
1035 return static_cast< Vec3<T0> >(*
this * v);
1039 template<
typename T0>
1045 w =
static_cast<T0
>(p[0] * MyBase::mm[ 3] + p[1] * MyBase::mm[ 7]
1046 + p[2] * MyBase::mm[11] + MyBase::mm[15]);
1049 return Vec3<T0>(
static_cast<T0
>((p[0] * MyBase::mm[ 0] + p[1] * MyBase::mm[ 4] +
1050 p[2] * MyBase::mm[ 8] + MyBase::mm[12]) / w),
1051 static_cast<T0
>((p[0] * MyBase::mm[ 1] + p[1] * MyBase::mm[ 5] +
1052 p[2] * MyBase::mm[ 9] + MyBase::mm[13]) / w),
1053 static_cast<T0
>((p[0] * MyBase::mm[ 2] + p[1] * MyBase::mm[ 6] +
1054 p[2] * MyBase::mm[10] + MyBase::mm[14]) / w));
1061 template<
typename T0>
1067 w = p[0] * MyBase::mm[12] + p[1] * MyBase::mm[13] + p[2] * MyBase::mm[14] + MyBase::mm[15];
1070 return Vec3<T0>(
static_cast<T0
>((p[0] * MyBase::mm[ 0] + p[1] * MyBase::mm[ 1] +
1071 p[2] * MyBase::mm[ 2] + MyBase::mm[ 3]) / w),
1072 static_cast<T0
>((p[0] * MyBase::mm[ 4] + p[1] * MyBase::mm[ 5] +
1073 p[2] * MyBase::mm[ 6] + MyBase::mm[ 7]) / w),
1074 static_cast<T0
>((p[0] * MyBase::mm[ 8] + p[1] * MyBase::mm[ 9] +
1075 p[2] * MyBase::mm[10] + MyBase::mm[11]) / w));
1082 template<
typename T0>
1086 static_cast<T0
>(v[0] * MyBase::mm[ 0] + v[1] * MyBase::mm[ 4] + v[2] * MyBase::mm[ 8]),
1087 static_cast<T0
>(v[0] * MyBase::mm[ 1] + v[1] * MyBase::mm[ 5] + v[2] * MyBase::mm[ 9]),
1088 static_cast<T0
>(v[0] * MyBase::mm[ 2] + v[1] * MyBase::mm[ 6] + v[2] * MyBase::mm[10]));
1093 bool invert(
Mat4<T> &inverse, T tolerance)
const;
1095 T det2(
const Mat4<T> &a,
int i0,
int i1,
int j0,
int j1)
const {
1098 return a.mm[i0row+j0]*a.mm[i1row+j1] - a.mm[i0row+j1]*a.mm[i1row+j0];
1101 T det3(
const Mat4<T> &a,
int i0,
int i1,
int i2,
1102 int j0,
int j1,
int j2)
const {
1104 return a.mm[i0row+j0]*det2(a, i1,i2, j1,j2) +
1105 a.mm[i0row+j1]*det2(a, i1,i2, j2,j0) +
1106 a.mm[i0row+j2]*det2(a, i1,i2, j0,j1);
1113 template <
typename T0,
typename T1>
1116 const T0 *t0 = m0.asPointer();
1117 const T1 *t1 = m1.asPointer();
1119 for (
int i=0; i<16; ++i)
if (!
isExactlyEqual(t0[i], t1[i]))
return false;
1125 template <
typename T0,
typename T1>
1130 template <
typename S,
typename T>
1138 template <
typename S,
typename T>
1148 template<
typename T,
typename MT>
1153 MT
const *m = _m.asPointer();
1155 _v[0]*m[0] + _v[1]*m[1] + _v[2]*m[2] + _v[3]*m[3],
1156 _v[0]*m[4] + _v[1]*m[5] + _v[2]*m[6] + _v[3]*m[7],
1157 _v[0]*m[8] + _v[1]*m[9] + _v[2]*m[10] + _v[3]*m[11],
1158 _v[0]*m[12] + _v[1]*m[13] + _v[2]*m[14] + _v[3]*m[15]);
1163 template<
typename T,
typename MT>
1168 MT
const *m = _m.asPointer();
1170 _v[0]*m[0] + _v[1]*m[4] + _v[2]*m[8] + _v[3]*m[12],
1171 _v[0]*m[1] + _v[1]*m[5] + _v[2]*m[9] + _v[3]*m[13],
1172 _v[0]*m[2] + _v[1]*m[6] + _v[2]*m[10] + _v[3]*m[14],
1173 _v[0]*m[3] + _v[1]*m[7] + _v[2]*m[11] + _v[3]*m[15]);
1178 template<
typename T,
typename MT>
1182 MT
const *m = _m.asPointer();
1184 _v[0]*m[0] + _v[1]*m[1] + _v[2]*m[2] + m[3],
1185 _v[0]*m[4] + _v[1]*m[5] + _v[2]*m[6] + m[7],
1186 _v[0]*m[8] + _v[1]*m[9] + _v[2]*m[10] + m[11]);
1191 template<
typename T,
typename MT>
1195 MT
const *m = _m.asPointer();
1197 _v[0]*m[0] + _v[1]*m[4] + _v[2]*m[8] + m[12],
1198 _v[0]*m[1] + _v[1]*m[5] + _v[2]*m[9] + m[13],
1199 _v[0]*m[2] + _v[1]*m[6] + _v[2]*m[10] + m[14]);
1204 template <
typename T0,
typename T1>
1215 template <
typename T0,
typename T1>
1226 template <
typename T0,
typename T1>
1239 template<
typename T0,
typename T1>
1243 static_cast<T1
>(m[0][0]*n[0] + m[0][1]*n[1] + m[0][2]*n[2]),
1244 static_cast<T1
>(m[1][0]*n[0] + m[1][1]*n[1] + m[1][2]*n[2]),
1245 static_cast<T1
>(m[2][0]*n[0] + m[2][1]*n[1] + m[2][2]*n[2]));
1250 template<
typename T>
1254 inverse.setIdentity();
1258 for (
int i = 0; i < 4; ++i) {
1260 double max = fabs(temp[i][i]);
1262 for (
int k = i+1; k < 4; ++k) {
1263 if (fabs(temp[k][i]) >
max) {
1265 max = fabs(temp[k][i]);
1274 for (
int k = 0; k < 4; ++k) {
1275 std::swap(temp[row][k], temp[i][k]);
1276 std::swap(inverse[row][k], inverse[i][k]);
1280 double pivot = temp[i][i];
1284 for (
int k = 0; k < 4; ++k) {
1285 temp[i][k] /=
pivot;
1286 inverse[i][k] /=
pivot;
1290 for (
int j = i+1; j < 4; ++j) {
1291 double t = temp[j][i];
1294 for (
int k = 0; k < 4; ++k) {
1295 temp[j][k] -= temp[i][k] * t;
1296 inverse[j][k] -= inverse[i][k] * t;
1303 for (
int i = 3; i > 0; --i) {
1304 for (
int j = 0; j < i; ++j) {
1305 double t = temp[j][i];
1308 for (
int k = 0; k < 4; ++k) {
1309 inverse[j][k] -= inverse[i][k]*t;
1314 return det*det >= tolerance*tolerance;
1317 template <
typename T>
1319 return (m.col(3) ==
Vec4<T>(0, 0, 0, 1));
1322 template <
typename T>
1324 return (m.row(3) !=
Vec4<T>(0, 0, 0, 1));
1327 template<
typename T>
1332 const T* ip = m.asPointer();
1333 T* op = out.asPointer();
1334 for (
unsigned i = 0; i < 16; ++i, ++op, ++ip) *op =
math::Abs(*ip);
1338 template<
typename Type1,
typename Type2>
1343 const Type1* ip = m.asPointer();
1344 Type1* op = out.asPointer();
1345 for (
unsigned i = 0; i < 16; ++i, ++op, ++ip) {
1353 template<
typename T>
1357 return cwiseLessThan<4, T>(m0, m1);
1360 template<
typename T>
1364 return cwiseGreaterThan<4, T>(m0, m1);
1371 #if OPENVDB_ABI_VERSION_NUMBER >= 8
1379 template<>
inline math::Mat4s zeroVal<math::Mat4s>() {
return math::Mat4s::zero(); }
1380 template<>
inline math::Mat4d zeroVal<math::Mat4d>() {
return math::Mat4d::zero(); }
1385 #endif // OPENVDB_UTIL_MAT4_H_HAS_BEEN_INCLUDED