1 /*
2  * Copyright 2013 The Android Open Source Project
3  *
4  * Licensed under the Apache License, Version 2.0 (the "License");
5  * you may not use this file except in compliance with the License.
6  * You may obtain a copy of the License at
7  *
8  *      http://www.apache.org/licenses/LICENSE-2.0
9  *
10  * Unless required by applicable law or agreed to in writing, software
11  * distributed under the License is distributed on an "AS IS" BASIS,
12  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13  * See the License for the specific language governing permissions and
14  * limitations under the License.
15  */
16 
17 #pragma once
18 
19 #include <math/mat3.h>
20 #include <math/quat.h>
21 #include <math/TMatHelpers.h>
22 #include <math/vec3.h>
23 #include <math/vec4.h>
24 
25 #include <stdint.h>
26 #include <sys/types.h>
27 #include <limits>
28 
29 #define PURE __attribute__((pure))
30 
31 #if __cplusplus >= 201402L
32 #define CONSTEXPR constexpr
33 #else
34 #define CONSTEXPR
35 #endif
36 
37 namespace android {
38 // -------------------------------------------------------------------------------------
39 namespace details {
40 
41 template<typename T>
42 class TQuaternion;
43 
44 /**
45  * A 4x4 column-major matrix class.
46  *
47  * Conceptually a 4x4 matrix is a an array of 4 column double4:
48  *
49  * mat4 m =
50  *      \f$
51  *      \left(
52  *      \begin{array}{cccc}
53  *      m[0] & m[1] & m[2] & m[3] \\
54  *      \end{array}
55  *      \right)
56  *      \f$
57  *      =
58  *      \f$
59  *      \left(
60  *      \begin{array}{cccc}
61  *      m[0][0] & m[1][0] & m[2][0] & m[3][0] \\
62  *      m[0][1] & m[1][1] & m[2][1] & m[3][1] \\
63  *      m[0][2] & m[1][2] & m[2][2] & m[3][2] \\
64  *      m[0][3] & m[1][3] & m[2][3] & m[3][3] \\
65  *      \end{array}
66  *      \right)
67  *      \f$
68  *      =
69  *      \f$
70  *      \left(
71  *      \begin{array}{cccc}
72  *      m(0,0) & m(0,1) & m(0,2) & m(0,3) \\
73  *      m(1,0) & m(1,1) & m(1,2) & m(1,3) \\
74  *      m(2,0) & m(2,1) & m(2,2) & m(2,3) \\
75  *      m(3,0) & m(3,1) & m(3,2) & m(3,3) \\
76  *      \end{array}
77  *      \right)
78  *      \f$
79  *
80  * m[n] is the \f$ n^{th} \f$ column of the matrix and is a double4.
81  *
82  */
83 template <typename T>
84 class TMat44 :  public TVecUnaryOperators<TMat44, T>,
85                 public TVecComparisonOperators<TMat44, T>,
86                 public TVecAddOperators<TMat44, T>,
87                 public TMatProductOperators<TMat44, T>,
88                 public TMatSquareFunctions<TMat44, T>,
89                 public TMatTransform<TMat44, T>,
90                 public TMatHelpers<TMat44, T>,
91                 public TMatDebug<TMat44, T> {
92 public:
93     enum no_init { NO_INIT };
94     typedef T value_type;
95     typedef T& reference;
96     typedef T const& const_reference;
97     typedef size_t size_type;
98     typedef TVec4<T> col_type;
99     typedef TVec4<T> row_type;
100 
101     static constexpr size_t COL_SIZE = col_type::SIZE;  // size of a column (i.e.: number of rows)
102     static constexpr size_t ROW_SIZE = row_type::SIZE;  // size of a row (i.e.: number of columns)
103     static constexpr size_t NUM_ROWS = COL_SIZE;
104     static constexpr size_t NUM_COLS = ROW_SIZE;
105 
106 private:
107     /*
108      *  <--  N columns  -->
109      *
110      *  a[0][0] a[1][0] a[2][0] ... a[N][0]    ^
111      *  a[0][1] a[1][1] a[2][1] ... a[N][1]    |
112      *  a[0][2] a[1][2] a[2][2] ... a[N][2]  M rows
113      *  ...                                    |
114      *  a[0][M] a[1][M] a[2][M] ... a[N][M]    v
115      *
116      *  COL_SIZE = M
117      *  ROW_SIZE = N
118      *  m[0] = [ a[0][0] a[0][1] a[0][2] ... a[0][M] ]
119      */
120 
121     col_type m_value[NUM_COLS];
122 
123 public:
124     // array access
125     inline constexpr col_type const& operator[](size_t column) const {
126 #if __cplusplus >= 201402L
127         // only possible in C++0x14 with constexpr
128         assert(column < NUM_COLS);
129 #endif
130         return m_value[column];
131     }
132 
133     inline col_type& operator[](size_t column) {
134         assert(column < NUM_COLS);
135         return m_value[column];
136     }
137 
138     // -----------------------------------------------------------------------
139     // we want the compiler generated versions for these...
140     TMat44(const TMat44&) = default;
141     ~TMat44() = default;
142     TMat44& operator = (const TMat44&) = default;
143 
144     /*
145      *  constructors
146      */
147 
148     // leaves object uninitialized. use with caution.
TMat44(no_init)149     explicit constexpr TMat44(no_init)
150             : m_value{ col_type(col_type::NO_INIT),
151                        col_type(col_type::NO_INIT),
152                        col_type(col_type::NO_INIT),
153                        col_type(col_type::NO_INIT) } {}
154 
155     /** initialize to identity.
156      *
157      *      \f$
158      *      \left(
159      *      \begin{array}{cccc}
160      *      1 & 0 & 0 & 0 \\
161      *      0 & 1 & 0 & 0 \\
162      *      0 & 0 & 1 & 0 \\
163      *      0 & 0 & 0 & 1 \\
164      *      \end{array}
165      *      \right)
166      *      \f$
167      */
168     CONSTEXPR TMat44();
169 
170     /** initialize to Identity*scalar.
171      *
172      *      \f$
173      *      \left(
174      *      \begin{array}{cccc}
175      *      v & 0 & 0 & 0 \\
176      *      0 & v & 0 & 0 \\
177      *      0 & 0 & v & 0 \\
178      *      0 & 0 & 0 & v \\
179      *      \end{array}
180      *      \right)
181      *      \f$
182      */
183     template<typename U>
184     explicit CONSTEXPR TMat44(U v);
185 
186     /** sets the diagonal to a vector.
187      *
188      *      \f$
189      *      \left(
190      *      \begin{array}{cccc}
191      *      v[0] & 0 & 0 & 0 \\
192      *      0 & v[1] & 0 & 0 \\
193      *      0 & 0 & v[2] & 0 \\
194      *      0 & 0 & 0 & v[3] \\
195      *      \end{array}
196      *      \right)
197      *      \f$
198      */
199     template <typename U>
200     explicit CONSTEXPR TMat44(const TVec4<U>& v);
201 
202     // construct from another matrix of the same size
203     template <typename U>
204     explicit CONSTEXPR TMat44(const TMat44<U>& rhs);
205 
206     /** construct from 4 column vectors.
207      *
208      *      \f$
209      *      \left(
210      *      \begin{array}{cccc}
211      *      v0 & v1 & v2 & v3 \\
212      *      \end{array}
213      *      \right)
214      *      \f$
215      */
216     template <typename A, typename B, typename C, typename D>
217     CONSTEXPR TMat44(const TVec4<A>& v0, const TVec4<B>& v1, const TVec4<C>& v2, const TVec4<D>& v3);
218 
219     /** construct from 16 elements in column-major form.
220      *
221      *      \f$
222      *      \left(
223      *      \begin{array}{cccc}
224      *      m[0][0] & m[1][0] & m[2][0] & m[3][0] \\
225      *      m[0][1] & m[1][1] & m[2][1] & m[3][1] \\
226      *      m[0][2] & m[1][2] & m[2][2] & m[3][2] \\
227      *      m[0][3] & m[1][3] & m[2][3] & m[3][3] \\
228      *      \end{array}
229      *      \right)
230      *      \f$
231      */
232     template <
233         typename A, typename B, typename C, typename D,
234         typename E, typename F, typename G, typename H,
235         typename I, typename J, typename K, typename L,
236         typename M, typename N, typename O, typename P>
237     CONSTEXPR TMat44(
238             A m00, B m01, C m02, D m03,
239             E m10, F m11, G m12, H m13,
240             I m20, J m21, K m22, L m23,
241             M m30, N m31, O m32, P m33);
242 
243     /**
244      * construct from a quaternion
245      */
246     template <typename U>
247     explicit CONSTEXPR TMat44(const TQuaternion<U>& q);
248 
249     /**
250      * construct from a C array in column major form.
251      */
252     template <typename U>
253     explicit CONSTEXPR TMat44(U const* rawArray);
254 
255     /**
256      * construct from a 3x3 matrix
257      */
258     template <typename U>
259     explicit CONSTEXPR TMat44(const TMat33<U>& matrix);
260 
261     /**
262      * construct from a 3x3 matrix and 3d translation
263      */
264     template <typename U, typename V>
265     CONSTEXPR TMat44(const TMat33<U>& matrix, const TVec3<V>& translation);
266 
267     /**
268      * construct from a 3x3 matrix and 4d last column.
269      */
270     template <typename U, typename V>
271     CONSTEXPR TMat44(const TMat33<U>& matrix, const TVec4<V>& column3);
272 
273     /*
274      *  helpers
275      */
276 
277     static CONSTEXPR TMat44 ortho(T left, T right, T bottom, T top, T near, T far);
278 
279     static CONSTEXPR TMat44 frustum(T left, T right, T bottom, T top, T near, T far);
280 
281     enum class Fov {
282         HORIZONTAL,
283         VERTICAL
284     };
285     static CONSTEXPR TMat44 perspective(T fov, T aspect, T near, T far, Fov direction = Fov::VERTICAL);
286 
287     template <typename A, typename B, typename C>
288     static CONSTEXPR TMat44 lookAt(const TVec3<A>& eye, const TVec3<B>& center, const TVec3<C>& up);
289 
290     template <typename A>
project(const TMat44 & projectionMatrix,TVec3<A> vertice)291     static CONSTEXPR TVec3<A> project(const TMat44& projectionMatrix, TVec3<A> vertice) {
292         TVec4<A> r = projectionMatrix * TVec4<A>{ vertice, 1 };
293         return r.xyz / r.w;
294     }
295 
296     template <typename A>
project(const TMat44 & projectionMatrix,TVec4<A> vertice)297     static CONSTEXPR TVec4<A> project(const TMat44& projectionMatrix, TVec4<A> vertice) {
298         vertice = projectionMatrix * vertice;
299         return { vertice.xyz / vertice.w, 1 };
300     }
301 
302     /**
303      * Constructs a 3x3 matrix from the upper-left corner of this 4x4 matrix
304      */
upperLeft()305     inline constexpr TMat33<T> upperLeft() const {
306         return TMat33<T>(m_value[0].xyz, m_value[1].xyz, m_value[2].xyz);
307     }
308 };
309 
310 // ----------------------------------------------------------------------------------------
311 // Constructors
312 // ----------------------------------------------------------------------------------------
313 
314 // Since the matrix code could become pretty big quickly, we don't inline most
315 // operations.
316 
317 template <typename T>
TMat44()318 CONSTEXPR TMat44<T>::TMat44() {
319     m_value[0] = col_type(1, 0, 0, 0);
320     m_value[1] = col_type(0, 1, 0, 0);
321     m_value[2] = col_type(0, 0, 1, 0);
322     m_value[3] = col_type(0, 0, 0, 1);
323 }
324 
325 template <typename T>
326 template <typename U>
TMat44(U v)327 CONSTEXPR TMat44<T>::TMat44(U v) {
328     m_value[0] = col_type(v, 0, 0, 0);
329     m_value[1] = col_type(0, v, 0, 0);
330     m_value[2] = col_type(0, 0, v, 0);
331     m_value[3] = col_type(0, 0, 0, v);
332 }
333 
334 template<typename T>
335 template<typename U>
TMat44(const TVec4<U> & v)336 CONSTEXPR TMat44<T>::TMat44(const TVec4<U>& v) {
337     m_value[0] = col_type(v.x, 0, 0, 0);
338     m_value[1] = col_type(0, v.y, 0, 0);
339     m_value[2] = col_type(0, 0, v.z, 0);
340     m_value[3] = col_type(0, 0, 0, v.w);
341 }
342 
343 // construct from 16 scalars
344 template<typename T>
345 template <
346     typename A, typename B, typename C, typename D,
347     typename E, typename F, typename G, typename H,
348     typename I, typename J, typename K, typename L,
349     typename M, typename N, typename O, typename P>
TMat44(A m00,B m01,C m02,D m03,E m10,F m11,G m12,H m13,I m20,J m21,K m22,L m23,M m30,N m31,O m32,P m33)350 CONSTEXPR TMat44<T>::TMat44(
351         A m00, B m01, C m02, D m03,
352         E m10, F m11, G m12, H m13,
353         I m20, J m21, K m22, L m23,
354         M m30, N m31, O m32, P m33) {
355     m_value[0] = col_type(m00, m01, m02, m03);
356     m_value[1] = col_type(m10, m11, m12, m13);
357     m_value[2] = col_type(m20, m21, m22, m23);
358     m_value[3] = col_type(m30, m31, m32, m33);
359 }
360 
361 template <typename T>
362 template <typename U>
TMat44(const TMat44<U> & rhs)363 CONSTEXPR TMat44<T>::TMat44(const TMat44<U>& rhs) {
364     for (size_t col = 0; col < NUM_COLS; ++col) {
365         m_value[col] = col_type(rhs[col]);
366     }
367 }
368 
369 // Construct from 4 column vectors.
370 template <typename T>
371 template <typename A, typename B, typename C, typename D>
TMat44(const TVec4<A> & v0,const TVec4<B> & v1,const TVec4<C> & v2,const TVec4<D> & v3)372 CONSTEXPR TMat44<T>::TMat44(
373         const TVec4<A>& v0, const TVec4<B>& v1,
374         const TVec4<C>& v2, const TVec4<D>& v3) {
375     m_value[0] = col_type(v0);
376     m_value[1] = col_type(v1);
377     m_value[2] = col_type(v2);
378     m_value[3] = col_type(v3);
379 }
380 
381 // Construct from raw array, in column-major form.
382 template <typename T>
383 template <typename U>
TMat44(U const * rawArray)384 CONSTEXPR TMat44<T>::TMat44(U const* rawArray) {
385     for (size_t col = 0; col < NUM_COLS; ++col) {
386         for (size_t row = 0; row < NUM_ROWS; ++row) {
387             m_value[col][row] = *rawArray++;
388         }
389     }
390 }
391 
392 template <typename T>
393 template <typename U>
TMat44(const TQuaternion<U> & q)394 CONSTEXPR TMat44<T>::TMat44(const TQuaternion<U>& q) {
395     const U n = q.x*q.x + q.y*q.y + q.z*q.z + q.w*q.w;
396     const U s = n > 0 ? 2/n : 0;
397     const U x = s*q.x;
398     const U y = s*q.y;
399     const U z = s*q.z;
400     const U xx = x*q.x;
401     const U xy = x*q.y;
402     const U xz = x*q.z;
403     const U xw = x*q.w;
404     const U yy = y*q.y;
405     const U yz = y*q.z;
406     const U yw = y*q.w;
407     const U zz = z*q.z;
408     const U zw = z*q.w;
409     m_value[0] = col_type(1-yy-zz,    xy+zw,    xz-yw,   0);
410     m_value[1] = col_type(  xy-zw,  1-xx-zz,    yz+xw,   0);  // NOLINT
411     m_value[2] = col_type(  xz+yw,    yz-xw,  1-xx-yy,   0);  // NOLINT
412     m_value[3] = col_type(      0,        0,        0,   1);  // NOLINT
413 }
414 
415 template <typename T>
416 template <typename U>
TMat44(const TMat33<U> & m)417 CONSTEXPR TMat44<T>::TMat44(const TMat33<U>& m) {
418     m_value[0] = col_type(m[0][0], m[0][1], m[0][2], 0);
419     m_value[1] = col_type(m[1][0], m[1][1], m[1][2], 0);
420     m_value[2] = col_type(m[2][0], m[2][1], m[2][2], 0);
421     m_value[3] = col_type(      0,       0,       0, 1);  // NOLINT
422 }
423 
424 template <typename T>
425 template <typename U, typename V>
TMat44(const TMat33<U> & m,const TVec3<V> & v)426 CONSTEXPR TMat44<T>::TMat44(const TMat33<U>& m, const TVec3<V>& v) {
427     m_value[0] = col_type(m[0][0], m[0][1], m[0][2], 0);
428     m_value[1] = col_type(m[1][0], m[1][1], m[1][2], 0);
429     m_value[2] = col_type(m[2][0], m[2][1], m[2][2], 0);
430     m_value[3] = col_type(   v[0],    v[1],    v[2], 1);  // NOLINT
431 }
432 
433 template <typename T>
434 template <typename U, typename V>
TMat44(const TMat33<U> & m,const TVec4<V> & v)435 CONSTEXPR TMat44<T>::TMat44(const TMat33<U>& m, const TVec4<V>& v) {
436     m_value[0] = col_type(m[0][0], m[0][1], m[0][2],    0);  // NOLINT
437     m_value[1] = col_type(m[1][0], m[1][1], m[1][2],    0);  // NOLINT
438     m_value[2] = col_type(m[2][0], m[2][1], m[2][2],    0);  // NOLINT
439     m_value[3] = col_type(   v[0],    v[1],    v[2], v[3]);  // NOLINT
440 }
441 
442 // ----------------------------------------------------------------------------------------
443 // Helpers
444 // ----------------------------------------------------------------------------------------
445 
446 template <typename T>
ortho(T left,T right,T bottom,T top,T near,T far)447 CONSTEXPR TMat44<T> TMat44<T>::ortho(T left, T right, T bottom, T top, T near, T far) {
448     TMat44<T> m;
449     m[0][0] =  2 / (right - left);
450     m[1][1] =  2 / (top   - bottom);
451     m[2][2] = -2 / (far   - near);
452     m[3][0] = -(right + left)   / (right - left);
453     m[3][1] = -(top   + bottom) / (top   - bottom);
454     m[3][2] = -(far   + near)   / (far   - near);
455     return m;
456 }
457 
458 template <typename T>
frustum(T left,T right,T bottom,T top,T near,T far)459 CONSTEXPR TMat44<T> TMat44<T>::frustum(T left, T right, T bottom, T top, T near, T far) {
460     TMat44<T> m;
461     m[0][0] =  (2 * near) / (right - left);
462     m[1][1] =  (2 * near) / (top   - bottom);
463     m[2][0] =  (right + left)   / (right - left);
464     m[2][1] =  (top   + bottom) / (top   - bottom);
465     m[2][2] = -(far   + near)   / (far   - near);
466     m[2][3] = -1;
467     m[3][2] = -(2 * far * near) / (far   - near);
468     m[3][3] =  0;
469     return m;
470 }
471 
472 template <typename T>
perspective(T fov,T aspect,T near,T far,TMat44::Fov direction)473 CONSTEXPR TMat44<T> TMat44<T>::perspective(T fov, T aspect, T near, T far, TMat44::Fov direction) {
474     T h;
475     T w;
476 
477     if (direction == TMat44::Fov::VERTICAL) {
478         h = std::tan(fov * M_PI / 360.0f) * near;
479         w = h * aspect;
480     } else {
481         w = std::tan(fov * M_PI / 360.0f) * near;
482         h = w / aspect;
483     }
484     return frustum(-w, w, -h, h, near, far);
485 }
486 
487 /*
488  * Returns a matrix representing the pose of a virtual camera looking towards -Z in its
489  * local Y-up coordinate system. "eye" is where the camera is located, "center" is the points its
490  * looking at and "up" defines where the Y axis of the camera's local coordinate system is.
491  */
492 template <typename T>
493 template <typename A, typename B, typename C>
lookAt(const TVec3<A> & eye,const TVec3<B> & center,const TVec3<C> & up)494 CONSTEXPR TMat44<T> TMat44<T>::lookAt(const TVec3<A>& eye, const TVec3<B>& center, const TVec3<C>& up) {
495     TVec3<T> z_axis(normalize(center - eye));
496     TVec3<T> norm_up(normalize(up));
497     if (std::abs(dot(z_axis, norm_up)) > 0.999) {
498         // Fix up vector if we're degenerate (looking straight up, basically)
499         norm_up = { norm_up.z, norm_up.x, norm_up.y };
500     }
501     TVec3<T> x_axis(normalize(cross(z_axis, norm_up)));
502     TVec3<T> y_axis(cross(x_axis, z_axis));
503     return TMat44<T>(
504             TVec4<T>(x_axis, 0),
505             TVec4<T>(y_axis, 0),
506             TVec4<T>(-z_axis, 0),
507             TVec4<T>(eye, 1));
508 }
509 
510 // ----------------------------------------------------------------------------------------
511 // Arithmetic operators outside of class
512 // ----------------------------------------------------------------------------------------
513 
514 /* We use non-friend functions here to prevent the compiler from using
515  * implicit conversions, for instance of a scalar to a vector. The result would
516  * not be what the caller expects.
517  *
518  * Also note that the order of the arguments in the inner loop is important since
519  * it determines the output type (only relevant when T != U).
520  */
521 
522 // matrix * column-vector, result is a vector of the same type than the input vector
523 template <typename T, typename U>
524 CONSTEXPR typename TMat44<T>::col_type PURE operator *(const TMat44<T>& lhs, const TVec4<U>& rhs) {
525     // Result is initialized to zero.
526     typename TMat44<T>::col_type result;
527     for (size_t col = 0; col < TMat44<T>::NUM_COLS; ++col) {
528         result += lhs[col] * rhs[col];
529     }
530     return result;
531 }
532 
533 // mat44 * vec3, result is vec3( mat44 * {vec3, 1} )
534 template <typename T, typename U>
535 CONSTEXPR typename TMat44<T>::col_type PURE operator *(const TMat44<T>& lhs, const TVec3<U>& rhs) {
536     return lhs * TVec4<U>{ rhs, 1 };
537 }
538 
539 
540 // row-vector * matrix, result is a vector of the same type than the input vector
541 template <typename T, typename U>
542 CONSTEXPR typename TMat44<U>::row_type PURE operator *(const TVec4<U>& lhs, const TMat44<T>& rhs) {
543     typename TMat44<U>::row_type result(TMat44<U>::row_type::NO_INIT);
544     for (size_t col = 0; col < TMat44<T>::NUM_COLS; ++col) {
545         result[col] = dot(lhs, rhs[col]);
546     }
547     return result;
548 }
549 
550 // matrix * scalar, result is a matrix of the same type than the input matrix
551 template <typename T, typename U>
552 constexpr typename std::enable_if<std::is_arithmetic<U>::value, TMat44<T>>::type PURE
553 operator *(TMat44<T> lhs, U rhs) {
554     return lhs *= rhs;
555 }
556 
557 // scalar * matrix, result is a matrix of the same type than the input matrix
558 template <typename T, typename U>
559 constexpr typename std::enable_if<std::is_arithmetic<U>::value, TMat44<T>>::type PURE
560 operator *(U lhs, const TMat44<T>& rhs) {
561     return rhs * lhs;
562 }
563 
564 // ----------------------------------------------------------------------------------------
565 
566 /* FIXME: this should go into TMatSquareFunctions<> but for some reason
567  * BASE<T>::col_type is not accessible from there (???)
568  */
569 template<typename T>
diag(const TMat44<T> & m)570 typename TMat44<T>::col_type PURE diag(const TMat44<T>& m) {
571     return matrix::diag(m);
572 }
573 
574 } // namespace details
575 
576 // ----------------------------------------------------------------------------------------
577 
578 typedef details::TMat44<double> mat4d;
579 typedef details::TMat44<float> mat4;
580 typedef details::TMat44<float> mat4f;
581 
582 // ----------------------------------------------------------------------------------------
583 }  // namespace android
584 
585 #undef PURE
586 #undef CONSTEXPR
587