31 #ifndef OPENVDB_MATH_MAT4_H_HAS_BEEN_INCLUDED
32 #define OPENVDB_MATH_MAT4_H_HAS_BEEN_INCLUDED
34 #include <openvdb/Exceptions.h>
35 #include <openvdb/Platform.h>
51 template<
typename T>
class Vec4;
75 template<
typename Source>
80 for (i = 0; i < 16; i++) {
92 template<
typename Source>
93 Mat4(Source a, Source b, Source c, Source d,
94 Source e, Source f, Source g, Source h,
95 Source i, Source j, Source k, Source l,
96 Source m, Source n, Source o, Source p)
120 template<
typename Source>
124 setBasis(v1, v2, v3, v4);
130 for (
int i = 0; i < 4; ++i) {
131 for (
int j = 0; j < 4; ++j) {
132 MyBase::mm[i*4 + j] = m[i][j];
138 template<
typename Source>
143 for (
int i=0; i<16; ++i) {
144 MyBase::mm[i] =
static_cast<T
>(src[i]);
163 MyBase::mm[i4+0] = v[0];
164 MyBase::mm[i4+1] = v[1];
165 MyBase::mm[i4+2] = v[2];
166 MyBase::mm[i4+3] = v[3];
173 return Vec4<T>((*this)(i,0), (*
this)(i,1), (*
this)(i,2), (*
this)(i,3));
180 MyBase::mm[ 0+j] = v[0];
181 MyBase::mm[ 4+j] = v[1];
182 MyBase::mm[ 8+j] = v[2];
183 MyBase::mm[12+j] = v[3];
190 return Vec4<T>((*this)(0,j), (*
this)(1,j), (*
this)(2,j), (*
this)(3,j));
197 const T*
operator[](
int i)
const {
return &(MyBase::mm[i<<2]); }
207 T& operator()(
int i,
int j)
211 return MyBase::mm[4*i+j];
217 T operator()(
int i,
int j)
const
221 return MyBase::mm[4*i+j];
228 MyBase::mm[ 0] = v1[0];
229 MyBase::mm[ 1] = v1[1];
230 MyBase::mm[ 2] = v1[2];
231 MyBase::mm[ 3] = v1[3];
233 MyBase::mm[ 4] = v2[0];
234 MyBase::mm[ 5] = v2[1];
235 MyBase::mm[ 6] = v2[2];
236 MyBase::mm[ 7] = v2[3];
238 MyBase::mm[ 8] = v3[0];
239 MyBase::mm[ 9] = v3[1];
240 MyBase::mm[10] = v3[2];
241 MyBase::mm[11] = v3[3];
243 MyBase::mm[12] = v4[0];
244 MyBase::mm[13] = v4[1];
245 MyBase::mm[14] = v4[2];
246 MyBase::mm[15] = v4[3];
299 for (
int i = 0; i < 3; i++)
300 for (
int j=0; j < 3; j++)
301 MyBase::mm[i*4+j] = m[i][j];
308 for (
int i = 0; i < 3; i++)
309 for (
int j = 0; j < 3; j++)
310 m[i][j] = MyBase::mm[i*4+j];
318 return Vec3<T>(MyBase::mm[12], MyBase::mm[13], MyBase::mm[14]);
323 MyBase::mm[12] = t[0];
324 MyBase::mm[13] = t[1];
325 MyBase::mm[14] = t[2];
329 template<
typename Source>
335 std::copy(src, (src + this->numElements()), MyBase::mm);
340 bool eq(
const Mat4 &m, T eps=1.0e-8)
const
342 for (
int i = 0; i < 16; i++) {
353 -MyBase::mm[ 0], -MyBase::mm[ 1], -MyBase::mm[ 2], -MyBase::mm[ 3],
354 -MyBase::mm[ 4], -MyBase::mm[ 5], -MyBase::mm[ 6], -MyBase::mm[ 7],
355 -MyBase::mm[ 8], -MyBase::mm[ 9], -MyBase::mm[10], -MyBase::mm[11],
356 -MyBase::mm[12], -MyBase::mm[13], -MyBase::mm[14], -MyBase::mm[15]
361 template <
typename S>
364 MyBase::mm[ 0] *= scalar;
365 MyBase::mm[ 1] *= scalar;
366 MyBase::mm[ 2] *= scalar;
367 MyBase::mm[ 3] *= scalar;
369 MyBase::mm[ 4] *= scalar;
370 MyBase::mm[ 5] *= scalar;
371 MyBase::mm[ 6] *= scalar;
372 MyBase::mm[ 7] *= scalar;
374 MyBase::mm[ 8] *= scalar;
375 MyBase::mm[ 9] *= scalar;
376 MyBase::mm[10] *= scalar;
377 MyBase::mm[11] *= scalar;
379 MyBase::mm[12] *= scalar;
380 MyBase::mm[13] *= scalar;
381 MyBase::mm[14] *= scalar;
382 MyBase::mm[15] *= scalar;
387 template <
typename S>
392 MyBase::mm[ 0] += s[ 0];
393 MyBase::mm[ 1] += s[ 1];
394 MyBase::mm[ 2] += s[ 2];
395 MyBase::mm[ 3] += s[ 3];
397 MyBase::mm[ 4] += s[ 4];
398 MyBase::mm[ 5] += s[ 5];
399 MyBase::mm[ 6] += s[ 6];
400 MyBase::mm[ 7] += s[ 7];
402 MyBase::mm[ 8] += s[ 8];
403 MyBase::mm[ 9] += s[ 9];
404 MyBase::mm[10] += s[10];
405 MyBase::mm[11] += s[11];
407 MyBase::mm[12] += s[12];
408 MyBase::mm[13] += s[13];
409 MyBase::mm[14] += s[14];
410 MyBase::mm[15] += s[15];
416 template <
typename S>
421 MyBase::mm[ 0] -= s[ 0];
422 MyBase::mm[ 1] -= s[ 1];
423 MyBase::mm[ 2] -= s[ 2];
424 MyBase::mm[ 3] -= s[ 3];
426 MyBase::mm[ 4] -= s[ 4];
427 MyBase::mm[ 5] -= s[ 5];
428 MyBase::mm[ 6] -= s[ 6];
429 MyBase::mm[ 7] -= s[ 7];
431 MyBase::mm[ 8] -= s[ 8];
432 MyBase::mm[ 9] -= s[ 9];
433 MyBase::mm[10] -= s[10];
434 MyBase::mm[11] -= s[11];
436 MyBase::mm[12] -= s[12];
437 MyBase::mm[13] -= s[13];
438 MyBase::mm[14] -= s[14];
439 MyBase::mm[15] -= s[15];
445 template <
typename S>
453 for (
int i = 0; i < 4; i++) {
455 MyBase::mm[i4+0] =
static_cast<T
>(s0[i4+0] * s1[ 0] +
460 MyBase::mm[i4+1] =
static_cast<T
>(s0[i4+0] * s1[ 1] +
465 MyBase::mm[i4+2] =
static_cast<T
>(s0[i4+0] * s1[ 2] +
470 MyBase::mm[i4+3] =
static_cast<T
>(s0[i4+0] * s1[ 3] +
482 MyBase::mm[ 0], MyBase::mm[ 4], MyBase::mm[ 8], MyBase::mm[12],
483 MyBase::mm[ 1], MyBase::mm[ 5], MyBase::mm[ 9], MyBase::mm[13],
484 MyBase::mm[ 2], MyBase::mm[ 6], MyBase::mm[10], MyBase::mm[14],
485 MyBase::mm[ 3], MyBase::mm[ 7], MyBase::mm[11], MyBase::mm[15]
515 T m0011 = m[0][0] * m[1][1];
516 T m0012 = m[0][0] * m[1][2];
517 T m0110 = m[0][1] * m[1][0];
518 T m0210 = m[0][2] * m[1][0];
519 T m0120 = m[0][1] * m[2][0];
520 T m0220 = m[0][2] * m[2][0];
522 T detA = m0011 * m[2][2] - m0012 * m[2][1] - m0110 * m[2][2]
523 + m0210 * m[2][1] + m0120 * m[1][2] - m0220 * m[1][1];
525 bool hasPerspective =
532 if (hasPerspective) {
533 det = m[0][3] * det3(m, 1,2,3, 0,2,1)
534 + m[1][3] * det3(m, 2,0,3, 0,2,1)
535 + m[2][3] * det3(m, 3,0,1, 0,2,1)
538 det = detA * m[3][3];
558 inv[0][0] = detA * ( m[1][1] * m[2][2] - m[1][2] * m[2][1]);
559 inv[0][1] = detA * (-m[0][1] * m[2][2] + m[0][2] * m[2][1]);
560 inv[0][2] = detA * ( m[0][1] * m[1][2] - m[0][2] * m[1][1]);
562 inv[1][0] = detA * (-m[1][0] * m[2][2] + m[1][2] * m[2][0]);
563 inv[1][1] = detA * ( m[0][0] * m[2][2] - m0220);
564 inv[1][2] = detA * ( m0210 - m0012);
566 inv[2][0] = detA * ( m[1][0] * m[2][1] - m[1][1] * m[2][0]);
567 inv[2][1] = detA * ( m0120 - m[0][0] * m[2][1]);
568 inv[2][2] = detA * ( m0011 - m0110);
570 if (hasPerspective) {
575 r[0] = m[3][0] * inv[0][0] + m[3][1] * inv[1][0]
576 + m[3][2] * inv[2][0];
577 r[1] = m[3][0] * inv[0][1] + m[3][1] * inv[1][1]
578 + m[3][2] * inv[2][1];
579 r[2] = m[3][0] * inv[0][2] + m[3][1] * inv[1][2]
580 + m[3][2] * inv[2][2];
583 p[0] = inv[0][0] * m[0][3] + inv[0][1] * m[1][3]
584 + inv[0][2] * m[2][3];
585 p[1] = inv[1][0] * m[0][3] + inv[1][1] * m[1][3]
586 + inv[1][2] * m[2][3];
587 p[2] = inv[2][0] * m[0][3] + inv[2][1] * m[1][3]
588 + inv[2][2] * m[2][3];
590 T h = m[3][3] - p.
dot(
Vec3<T>(m[3][0],m[3][1],m[3][2]));
601 inv[3][0] = -h * r[0];
602 inv[3][1] = -h * r[1];
603 inv[3][2] = -h * r[2];
605 inv[0][3] = -h * p[0];
606 inv[1][3] = -h * p[1];
607 inv[2][3] = -h * p[2];
613 inv[0][0] += p[0] * r[0];
614 inv[0][1] += p[0] * r[1];
615 inv[0][2] += p[0] * r[2];
616 inv[1][0] += p[1] * r[0];
617 inv[1][1] += p[1] * r[1];
618 inv[1][2] += p[1] * r[2];
619 inv[2][0] += p[2] * r[0];
620 inv[2][1] += p[2] * r[1];
621 inv[2][2] += p[2] * r[2];
625 inv[3][0] = - (m[3][0] * inv[0][0] + m[3][1] * inv[1][0]
626 + m[3][2] * inv[2][0]);
627 inv[3][1] = - (m[3][0] * inv[0][1] + m[3][1] * inv[1][1]
628 + m[3][2] * inv[2][1]);
629 inv[3][2] = - (m[3][0] * inv[0][2] + m[3][1] * inv[1][2]
630 + m[3][2] * inv[2][2]);
654 for (i = 0; i < 4; i++) {
655 ap = &MyBase::mm[ 0];
657 for (j = 0; j < 4; j++) {
658 for (k = 0; k < 4; k++) {
659 if ((k != i) && (j != 0)) {
666 det += sign * MyBase::mm[i] * submat.
det();
678 {
return snapBasis(*
this, axis, direction);}
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.
x();
777 MyBase::mm[ 2] *= v.
x();
779 MyBase::mm[ 4] *= v.
y();
780 MyBase::mm[ 5] *= v.
y();
781 MyBase::mm[ 6] *= v.
y();
783 MyBase::mm[ 8] *= v.
z();
784 MyBase::mm[ 9] *= v.
z();
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 = p[0] * MyBase::mm[ 3] + p[1] * MyBase::mm[ 7] + p[2] * MyBase::mm[11] + MyBase::mm[15];
1048 return Vec3<T0>(
static_cast<T0
>((p[0] * MyBase::mm[ 0] + p[1] * MyBase::mm[ 4] +
1049 p[2] * MyBase::mm[ 8] + MyBase::mm[12]) / w),
1050 static_cast<T0
>((p[0] * MyBase::mm[ 1] + p[1] * MyBase::mm[ 5] +
1051 p[2] * MyBase::mm[ 9] + MyBase::mm[13]) / w),
1052 static_cast<T0
>((p[0] * MyBase::mm[ 2] + p[1] * MyBase::mm[ 6] +
1053 p[2] * MyBase::mm[10] + MyBase::mm[14]) / w));
1060 template<
typename T0>
1066 w = p[0] * MyBase::mm[12] + p[1] * MyBase::mm[13] + p[2] * MyBase::mm[14] + MyBase::mm[15];
1069 return Vec3<T0>(
static_cast<T0
>((p[0] * MyBase::mm[ 0] + p[1] * MyBase::mm[ 1] +
1070 p[2] * MyBase::mm[ 2] + MyBase::mm[ 3]) / w),
1071 static_cast<T0
>((p[0] * MyBase::mm[ 4] + p[1] * MyBase::mm[ 5] +
1072 p[2] * MyBase::mm[ 6] + MyBase::mm[ 7]) / w),
1073 static_cast<T0
>((p[0] * MyBase::mm[ 8] + p[1] * MyBase::mm[ 9] +
1074 p[2] * MyBase::mm[10] + MyBase::mm[11]) / w));
1081 template<
typename T0>
1085 static_cast<T0
>(v[0] * MyBase::mm[ 0] + v[1] * MyBase::mm[ 4] + v[2] * MyBase::mm[ 8]),
1086 static_cast<T0>(v[0] * MyBase::mm[ 1] + v[1] * MyBase::mm[ 5] + v[2] * MyBase::mm[ 9]),
1087 static_cast<T0
>(v[0] * MyBase::mm[ 2] + v[1] * MyBase::mm[ 6] + v[2] * MyBase::mm[10]));
1096 #ifdef ALLOW_CAST_TO_POINTER
1099 #
if (__GNUC__ > 3 || (__GNUC__ == 3 && __GNUC_MINOR__ >= 4))
1100 __attribute__ ((deprecated))
1102 {
return MyBase::mm;}
1104 #
if (__GNUC__ > 3 || (__GNUC__ == 3 && __GNUC_MINOR__ >= 4))
1105 __attribute__ ((deprecated))
1107 {
return MyBase::mm;}
1116 int index0 =
static_cast<int>(axis0);
1117 int index1 =
static_cast<int>(axis1);
1119 MyBase::mm[index1 * 4 + 0] += shear * MyBase::mm[index0 * 4 ];
1120 MyBase::mm[index1 * 4 + 1] += shear * MyBase::mm[index0 * 4 + 1];
1121 MyBase::mm[index1 * 4 + 2] += shear * MyBase::mm[index0 * 4 + 2];
1122 MyBase::mm[index1 * 4 + 3] += shear * MyBase::mm[index0 * 4 + 3];
1127 template <
typename T0>
1136 template <
typename T0>
1139 this->mm[12] += v.
dot(
Vec3<T>(this->mm[ 0], this->mm[ 4], this->mm[ 8]));
1140 this->mm[13] += v.
dot(
Vec3<T>(this->mm[ 1], this->mm[ 5], this->mm[ 9]));
1141 this->mm[14] += v.
dot(
Vec3<T>(this->mm[ 2], this->mm[ 6], this->mm[10]));
1142 this->mm[15] += v.
dot(
Vec3<T>(this->mm[ 3], this->mm[ 7], this->mm[11]));
1151 T c =
static_cast<T
>(cos(angle));
1152 T s =
static_cast<T
>(sin(angle));
1159 a10 = c*this->mm[ 4] + s*this->mm[ 8];
1160 a11 = c*this->mm[ 5] + s*this->mm[ 9];
1161 a12 = c*this->mm[ 6] + s*this->mm[10];
1163 this->mm[ 8] = -s*this->mm[ 4] + c*this->mm[ 8];
1164 this->mm[ 9] = -s*this->mm[ 5] + c*this->mm[ 9];
1165 this->mm[10] = -s*this->mm[ 6] + c*this->mm[10];
1174 a00 = c*this->mm[ 0] - s*this->mm[ 8];
1175 a01 = c*this->mm[ 1] - s*this->mm[ 9];
1176 a02 = c*this->mm[ 2] - s*this->mm[10];
1178 this->mm[ 8] = s*this->mm[ 0] + c*this->mm[ 8];
1179 this->mm[ 9] = s*this->mm[ 1] + c*this->mm[ 9];
1180 this->mm[10] = s*this->mm[ 2] + c*this->mm[10];
1189 a00 = c*this->mm[ 0] + s*this->mm[ 4];
1190 a01 = c*this->mm[ 1] + s*this->mm[ 5];
1191 a02 = c*this->mm[ 2] + s*this->mm[ 6];
1193 this->mm[ 4] = -s*this->mm[ 0] + c*this->mm[ 4];
1194 this->mm[ 5] = -s*this->mm[ 1] + c*this->mm[ 5];
1195 this->mm[ 6] = -s*this->mm[ 2] + c*this->mm[ 6];
1211 template <
typename T0,
typename T1>
1216 this->mm[ 0] = src1[0] + src2[0];
1217 this->mm[ 1] = src1[1] + src2[1];
1218 this->mm[ 2] = src1[2] + src2[2];
1219 this->mm[ 3] = src1[3] + src2[3];
1221 this->mm[ 4] = src1[4] + src2[4];
1222 this->mm[ 5] = src1[5] + src2[5];
1223 this->mm[ 6] = src1[6] + src2[6];
1224 this->mm[ 7] = src1[7] + src2[7];
1226 this->mm[ 8] = src1[8] + src2[8];
1227 this->mm[ 9] = src1[9] + src2[9];
1228 this->mm[10] = src1[10] + src2[10];
1229 this->mm[11] = src1[11] + src2[11];
1231 this->mm[12] = src1[12] + src2[12];
1232 this->mm[13] = src1[13] + src2[13];
1233 this->mm[14] = src1[14] + src2[14];
1234 this->mm[15] = src1[15] + src2[15];
1240 template <
typename T0,
typename T1>
1245 this->mm[ 0] = src1[0] - src2[0];
1246 this->mm[ 1] = src1[1] - src2[1];
1247 this->mm[ 2] = src1[2] - src2[2];
1248 this->mm[ 3] = src1[3] - src2[3];
1250 this->mm[ 4] = src1[4] - src2[4];
1251 this->mm[ 5] = src1[5] - src2[5];
1252 this->mm[ 6] = src1[6] - src2[6];
1253 this->mm[ 7] = src1[7] - src2[7];
1255 this->mm[ 8] = src1[8] - src2[8];
1256 this->mm[ 9] = src1[9] - src2[9];
1257 this->mm[10] = src1[10] - src2[10];
1258 this->mm[11] = src1[11] - src2[11];
1260 this->mm[12] = src1[12] - src2[12];
1261 this->mm[13] = src1[13] - src2[13];
1262 this->mm[14] = src1[14] - src2[14];
1263 this->mm[15] = src1[15] - src2[15];
1268 template <
typename T0,
typename T1>
1272 this->mm[ 0] = scalar * src[0];
1273 this->mm[ 1] = scalar * src[1];
1274 this->mm[ 2] = scalar * src[2];
1275 this->mm[ 3] = scalar * src[3];
1277 this->mm[ 4] = scalar * src[4];
1278 this->mm[ 5] = scalar * src[5];
1279 this->mm[ 6] = scalar * src[6];
1280 this->mm[ 7] = scalar * src[7];
1282 this->mm[ 8] = scalar * src[8];
1283 this->mm[ 9] = scalar * src[9];
1284 this->mm[10] = scalar * src[10];
1285 this->mm[11] = scalar * src[11];
1287 this->mm[12] = scalar * src[12];
1288 this->mm[13] = scalar * src[13];
1289 this->mm[14] = scalar * src[14];
1290 this->mm[15] = scalar * src[15];
1296 template <
typename T0,
typename T1>
1303 for (i = 0; i < 4; i++) {
1305 this->mm[i4+0] =
static_cast<T
>(src1[i4+0]*src2[4*0+0] +
1306 src1[i4+1]*src2[4*1+0] +
1307 src1[i4+2]*src2[4*2+0] +
1308 src1[i4+3]*src2[4*3+0]);
1310 this->mm[i4+1] =
static_cast<T
>(src1[i4+0]*src2[4*0+1] +
1311 src1[i4+1]*src2[4*1+1] +
1312 src1[i4+2]*src2[4*2+1] +
1313 src1[i4+3]*src2[4*3+1]);
1315 this->mm[i4+2] =
static_cast<T
>(src1[i4+0]*src2[4*0+2] +
1316 src1[i4+1]*src2[4*1+2] +
1317 src1[i4+2]*src2[4*2+2] +
1318 src1[i4+3]*src2[4*3+2]);
1320 this->mm[i4+3] =
static_cast<T
>(src1[i4+0]*src2[4*0+3] +
1321 src1[i4+1]*src2[4*1+3] +
1322 src1[i4+2]*src2[4*2+3] +
1323 src1[i4+3]*src2[4*3+3]);
1334 T det2(
const Mat4<T> &a,
int i0,
int i1,
int j0,
int j1)
const {
1337 return a.
mm[i0row+j0]*a.
mm[i1row+j1] - a.
mm[i0row+j1]*a.
mm[i1row+j0];
1340 T det3(
const Mat4<T> &a,
int i0,
int i1,
int i2,
1341 int j0,
int j1,
int j2)
const {
1343 return a.mm[i0row+j0]*det2(a, i1,i2, j1,j2) +
1344 a.mm[i0row+j1]*det2(a, i1,i2, j2,j0) +
1345 a.mm[i0row+j2]*det2(a, i1,i2, j0,j1);
1348 static const Mat4<T> sIdentity;
1349 static const Mat4<T> sZero;
1353 template <
typename T>
1354 const Mat4<T> Mat4<T>::sIdentity = Mat4<T>(1, 0, 0, 0,
1359 template <
typename T>
1360 const Mat4<T> Mat4<T>::sZero = Mat4<T>(0, 0, 0, 0,
1367 template <
typename T0,
typename T1>
1373 for (
int i=0; i<16; ++i)
if (!
isExactlyEqual(t0[i], t1[i]))
return false;
1379 template <
typename T0,
typename T1>
1384 template <
typename S,
typename T>
1392 template <
typename S,
typename T>
1402 template<
typename T,
typename MT>
1409 _v[0]*m[0] + _v[1]*m[1] + _v[2]*m[2] + _v[3]*m[3],
1410 _v[0]*m[4] + _v[1]*m[5] + _v[2]*m[6] + _v[3]*m[7],
1411 _v[0]*m[8] + _v[1]*m[9] + _v[2]*m[10] + _v[3]*m[11],
1412 _v[0]*m[12] + _v[1]*m[13] + _v[2]*m[14] + _v[3]*m[15]);
1417 template<
typename T,
typename MT>
1424 _v[0]*m[0] + _v[1]*m[4] + _v[2]*m[8] + _v[3]*m[12],
1425 _v[0]*m[1] + _v[1]*m[5] + _v[2]*m[9] + _v[3]*m[13],
1426 _v[0]*m[2] + _v[1]*m[6] + _v[2]*m[10] + _v[3]*m[14],
1427 _v[0]*m[3] + _v[1]*m[7] + _v[2]*m[11] + _v[3]*m[15]);
1433 template<
typename T,
typename MT>
1440 _v[0]*m[0] + _v[1]*m[1] + _v[2]*m[2] + m[3],
1441 _v[0]*m[4] + _v[1]*m[5] + _v[2]*m[6] + m[7],
1442 _v[0]*m[8] + _v[1]*m[9] + _v[2]*m[10] + m[11]);
1448 template<
typename T,
typename MT>
1455 _v[0]*m[0] + _v[1]*m[4] + _v[2]*m[8] + m[12],
1456 _v[0]*m[1] + _v[1]*m[5] + _v[2]*m[9] + m[13],
1457 _v[0]*m[2] + _v[1]*m[6] + _v[2]*m[10] + m[14]);
1462 template <
typename T0,
typename T1>
1473 template <
typename T0,
typename T1>
1485 template <
typename T0,
typename T1>
1498 template<
typename T0,
typename T1>
1502 static_cast<T1
>(m[0][0]*n[0] + m[0][1]*n[1] + m[0][2]*n[2]),
1503 static_cast<T1>(m[1][0]*n[0] + m[1][1]*n[1] + m[1][2]*n[2]),
1504 static_cast<T1
>(m[2][0]*n[0] + m[2][1]*n[1] + m[2][2]*n[2]));
1509 template<
typename T>
1510 bool Mat4<T>::invert(Mat4<T> &inverse, T tolerance)
const
1512 Mat4<T> temp(*
this);
1513 inverse.setIdentity();
1517 for (
int i = 0; i < 4; ++i) {
1519 double max = fabs(temp[i][i]);
1521 for (
int k = i+1; k < 4; ++k) {
1522 if (fabs(temp[k][i]) > max) {
1524 max = fabs(temp[k][i]);
1533 for (
int k = 0; k < 4; ++k) {
1534 std::swap(temp[row][k], temp[i][k]);
1535 std::swap(inverse[row][k], inverse[i][k]);
1539 double pivot = temp[i][i];
1543 for (
int k = 0; k < 4; ++k) {
1544 temp[i][k] /= pivot;
1545 inverse[i][k] /= pivot;
1549 for (
int j = i+1; j < 4; ++j) {
1550 double t = temp[j][i];
1553 for (
int k = 0; k < 4; ++k) {
1554 temp[j][k] -= temp[i][k] * t;
1555 inverse[j][k] -= inverse[i][k] * t;
1562 for (
int i = 3; i > 0; --i) {
1563 for (
int j = 0; j < i; ++j) {
1564 double t = temp[j][i];
1567 for (
int k = 0; k < 4; ++k) {
1568 inverse[j][k] -= inverse[i][k]*t;
1573 return det*det >= tolerance*tolerance;
1576 template <
typename T>
1581 template <
typename T>
1590 #if DWREAL_IS_DOUBLE == 1
1594 #endif // DWREAL_IS_DOUBLE
1601 return math::Mat4s::identity();
1605 return math::Mat4d::identity();
1611 #endif // OPENVDB_UTIL_MAT4_H_HAS_BEEN_INCLUDED