Direct-BT v3.3.0-1-gc2d430c
Direct-BT - Direct Bluetooth Programming.
mat4f.hpp
Go to the documentation of this file.
1/*
2 * Author: Sven Gothel <sgothel@jausoft.com>
3 * Copyright (c) 2014-2024 Gothel Software e.K.
4 *
5 * Permission is hereby granted, free of charge, to any person obtaining
6 * a copy of this software and associated documentation files (the
7 * "Software"), to deal in the Software without restriction, including
8 * without limitation the rights to use, copy, modify, merge, publish,
9 * distribute, sublicense, and/or sell copies of the Software, and to
10 * permit persons to whom the Software is furnished to do so, subject to
11 * the following conditions:
12 *
13 * The above copyright notice and this permission notice shall be
14 * included in all copies or substantial portions of the Software.
15 *
16 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
17 * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
18 * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
19 * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE
20 * LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION
21 * OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION
22 * WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
23 */
24#ifndef JAU_MAT4f_HPP_
25#define JAU_MAT4f_HPP_
26
27#include <cmath>
28#include <cstdarg>
29#include <cstdint>
30#include <cassert>
31#include <limits>
32#include <string>
33#include <vector>
34#include <initializer_list>
35#include <iostream>
36
37#include <jau/float_math.hpp>
39#include <jau/math/vec3f.hpp>
40#include <jau/math/vec4f.hpp>
41#include <jau/math/recti.hpp>
43
44namespace jau::math::geom {
45 class Frustum; // forward
46}
47
48namespace jau::math {
49
50 /** \addtogroup Math
51 *
52 * @{
53 */
54
55 template<typename Value_type,
56 std::enable_if_t<std::is_floating_point_v<Value_type>, bool>>
57 class Quaternion; // forward
58
59/**
60 * Basic 4x4 value_type matrix implementation using fields for intensive use-cases (host operations).
61 * <p>
62 * Implementation covers {@link FloatUtil} matrix functionality, exposed in an object oriented manner.
63 * </p>
64 * <p>
65 * Unlike {@link com.jogamp.math.util.PMVmat4f PMVmat4f}, this class only represents one single matrix.
66 * </p>
67 * <p>
68 * For array operations the layout is expected in column-major order
69 * matching OpenGL's implementation, illustration:
70 * <pre>
71 Row-Major Column-Major (OpenGL):
72
73 | 0 1 2 tx |
74 | |
75 | 4 5 6 ty |
76 M = | |
77 | 8 9 10 tz |
78 | |
79 | 12 13 14 15 |
80
81 R C R C
82 m[0*4+3] = tx; m[0+4*3] = tx;
83 m[1*4+3] = ty; m[1+4*3] = ty;
84 m[2*4+3] = tz; m[2+4*3] = tz;
85
86 RC (std subscript order) RC (std subscript order)
87 m03 = tx; m03 = tx;
88 m13 = ty; m13 = ty;
89 m23 = tz; m23 = tz;
90
91 * </pre>
92 * </p>
93 * <p>
94 * <ul>
95 * <li><a href="http://web.archive.org/web/20041029003853/http://www.j3d.org/matrix_faq/matrfaq_latest.html">Matrix-FAQ</a></li>
96 * <li><a href="https://en.wikipedia.org/wiki/Matrix_%28mathematics%29">Wikipedia-Matrix</a></li>
97 * <li><a href="http://www.euclideanspace.com/maths/algebra/matrix/index.htm">euclideanspace.com-Matrix</a></li>
98 * </ul>
99 * </p>
100 * <p>
101 * Implementation utilizes unrolling of small vertices and matrices wherever possible
102 * while trying to access memory in a linear fashion for performance reasons, see:
103 * <ul>
104 * <li><a href="https://lessthanoptimal.github.io/Java-Matrix-Benchmark/">java-matrix-benchmark</a></li>
105 * <li><a href="https://github.com/lessthanoptimal/ejml">EJML Efficient Java Matrix Library</a></li>
106 * </ul>
107 * </p>
108 */
109
110template<typename Value_type,
111 std::enable_if_t<std::is_floating_point_v<Value_type>, bool> = true>
112class alignas(Value_type) Matrix4 {
113 public:
116 typedef const value_type* const_pointer;
121
125
126 constexpr static const value_type zero = value_type(0);
127 constexpr static const value_type one = value_type(1);
128 constexpr static const value_type two = value_type(2);
129 constexpr static const value_type half = one/two;
130
131 /**
132 * Inversion Epsilon, used with equals method to determine if two inverted matrices are close enough to be considered equal.
133 * <p>
134 * Using {@value}, which is ~84 times `std::numeric_limits<value_type>::epsilon()`.
135 * </p>
136 */
137 constexpr static const value_type inv_deviation = value_type(84) * std::numeric_limits<value_type>::epsilon(); // 84 * EPSILON(1.1920929E-7f) = 1.0E-5f
138
139 private:
140 // RC
141 value_type m00, m10, m20, m30; // column 0
142 value_type m01, m11, m21, m31; // column 1
143 value_type m02, m12, m22, m32; // column 2
144 value_type m03, m13, m23, m33; // column 3
145
146 friend geom::Frustum;
148
149 public:
150
151 /**
152 * Creates a new identity matrix.
153 */
154 constexpr Matrix4() noexcept
155 : m00(one), m10(zero), m20(zero), m30(zero),
156 m01(zero), m11(one), m21(zero), m31(zero),
157 m02(zero), m12(zero), m22(one), m32(zero),
158 m03(zero), m13(zero), m23(zero), m33(one)
159 { }
160
161 /**
162 * Creates a new matrix based on given value_type[4*4] column major order.
163 * @param m 4x4 matrix in column-major order
164 */
165 constexpr Matrix4(const_iterator m) noexcept
166 : m00(*m), m10(*(++m)), m20(*(++m)), m30(*(++m)), // column 0
167 m01(*(++m)), m11(*(++m)), m21(*(++m)), m31(*(++m)), // column 1
168 m02(*(++m)), m12(*(++m)), m22(*(++m)), m32(*(++m)), // column 2
169 m03(*(++m)), m13(*(++m)), m23(*(++m)), m33(*(++m)) // column 3
170 {}
171
172 /**
173 * Creates a new matrix based on given value_type initializer list in column major order.
174 * @param m source initializer list value_type data to be copied into this new instance, implied size must be >= 16
175 */
176 constexpr Matrix4(std::initializer_list<value_type> m) noexcept
177 : Matrix4( m.begin() )
178 {
179 assert(m.size() >= 16 );
180 }
181
182 /**
183 * Creates a new matrix copying the values of the given {@code src} matrix.
184 */
185 constexpr Matrix4(const Matrix4& o) noexcept
186 : Matrix4( o.cbegin() )
187 { }
188
189 /**
190 * Copy assignment using the the values of the given {@code src} matrix.
191 */
192 constexpr Matrix4& operator=(const Matrix4& o) noexcept { return load(o); }
193
194 constexpr bool equals(const Matrix4& o, const value_type epsilon=std::numeric_limits<value_type>::epsilon()) const noexcept {
195 if( this == &o ) {
196 return true;
197 } else {
198 return jau::equals(m00, o.m00, epsilon) &&
199 jau::equals(m01, o.m01, epsilon) &&
200 jau::equals(m02, o.m02, epsilon) &&
201 jau::equals(m03, o.m03, epsilon) &&
202 jau::equals(m10, o.m10, epsilon) &&
203 jau::equals(m11, o.m11, epsilon) &&
204 jau::equals(m12, o.m12, epsilon) &&
205 jau::equals(m13, o.m13, epsilon) &&
206 jau::equals(m20, o.m20, epsilon) &&
207 jau::equals(m21, o.m21, epsilon) &&
208 jau::equals(m22, o.m22, epsilon) &&
209 jau::equals(m23, o.m23, epsilon) &&
210 jau::equals(m30, o.m30, epsilon) &&
211 jau::equals(m31, o.m31, epsilon) &&
212 jau::equals(m32, o.m32, epsilon) &&
213 jau::equals(m33, o.m33, epsilon);
214 }
215 }
216 constexpr bool operator==(const Matrix4& rhs) const noexcept { return equals(rhs); }
217
218 //
219 // Write to Matrix via set(..) or load(..)
220 //
221
222 /**
223 * Returns writable reference to the {@code i}th component of this column-major order matrix, 0 <= i < 16 w/o boundary check
224 */
225 constexpr reference operator[](size_t i) noexcept {
226 assert( i < 16 );
227 return (&m00)[i];
228 }
229
230 /** Sets the {@code i}th component of this column-major order matrix with value_type {@code v}, 0 <= i < 16 w/o boundary check*/
231 constexpr void set(const jau::nsize_t i, const value_type v) noexcept {
232 assert( i < 16 );
233 (&m00)[i] = v;
234 }
235
236 explicit operator pointer() noexcept { return &m00; }
237 constexpr iterator begin() noexcept { return &m00; }
238
239 /**
240 * Set this matrix to identity.
241 * <pre>
242 Translation matrix (Column Order):
243 1 0 0 0
244 0 1 0 0
245 0 0 1 0
246 0 0 0 1
247 * </pre>
248 * @return this matrix for chaining
249 */
250 constexpr Matrix4& loadIdentity() noexcept {
251 m00 = m11 = m22 = m33 = one;
252 m01 = m02 = m03 =
253 m10 = m12 = m13 =
254 m20 = m21 = m23 =
255 m30 = m31 = m32 = zero;
256 return *this;
257 }
258
259 /**
260 * Load the values of the given matrix {@code src} to this matrix w/o boundary check.
261 * @param src 4x4 matrix value_type[16] in column-major order
262 * @return this matrix for chaining
263 */
264 constexpr Matrix4& load(const_iterator src) noexcept {
265 // RC
266 m00 = *src; // column 0
267 m10 = *(++src);
268 m20 = *(++src);
269 m30 = *(++src);
270 m01 = *(++src); // column 1
271 m11 = *(++src);
272 m21 = *(++src);
273 m31 = *(++src);
274 m02 = *(++src); // column 2
275 m12 = *(++src);
276 m22 = *(++src);
277 m32 = *(++src);
278 m03 = *(++src); // column 3
279 m13 = *(++src);
280 m23 = *(++src);
281 m33 = *(++src);
282 return *this;
283 }
284 /**
285 * Load the values of the given matrix {@code src} to this matrix w/o boundary check
286 * @param src the source values
287 * @return this matrix for chaining
288 */
289 constexpr Matrix4& load(const Matrix4& src) noexcept {
290 return load( src.cbegin() );
291 }
292
293 //
294 // Read out Matrix via get(..)
295 //
296
297 /**
298 * Returns read-only {@code i}th component of the given column-major order matrix, 0 <= i < 16 w/o boundary check
299 */
300 constexpr value_type operator[](size_t i) const noexcept {
301 assert( i < 16 );
302 return (&m00)[i];
303 }
304
305 /** Returns the {@code i}th component of the given column-major order matrix, 0 <= i < 16, w/o boundary check */
306 constexpr value_type get(const jau::nsize_t i) const noexcept {
307 assert( i < 16 );
308 return (&m00)[i];
309 }
310
311 explicit operator const_pointer() const noexcept { return &m00; }
312 constexpr const_iterator cbegin() const noexcept { return &m00; }
313
314 /**
315 * Get the named column of the given column-major matrix to v_out w/o boundary check.
316 * @param column named column to copy
317 * @param v_out the column-vector storage
318 * @return given result vector <i>v_out</i> for chaining
319 */
320 constexpr Vec4& getColumn(const jau::nsize_t column, Vec4& v_out) const noexcept {
321 assert( column < 4 );
322 return v_out.set( get(0+column*4),
323 get(1+column*4),
324 get(2+column*4),
325 get(3+column*4) );
326 }
327
328 /**
329 * Get the named column of the given column-major matrix to v_out w/o boundary check.
330 * @param column named column to copy
331 * @return result vector holding the requested column
332 */
333 constexpr Vec4 getColumn(const jau::nsize_t column) const noexcept {
334 assert( column < 4 );
335 return Vec4( get(0+column*4),
336 get(1+column*4),
337 get(2+column*4),
338 get(3+column*4) );
339 }
340
341 /**
342 * Get the named column of the given column-major matrix to v_out w/o boundary check.
343 * @param column named column to copy
344 * @param v_out the column-vector storage
345 * @return given result vector <i>v_out</i> for chaining
346 */
347 constexpr Vec3& getColumn(const jau::nsize_t column, Vec3& v_out) const noexcept {
348 return v_out.set( get(0+column*4),
349 get(1+column*4),
350 get(2+column*4) );
351 }
352
353 /**
354 * Get the named row of the given column-major matrix to v_out w/ boundary check.
355 * @param row named row to copy
356 * @param v_out the row-vector storage
357 * @return given result vector <i>v_out</i> for chaining
358 */
359 constexpr Vec4& getRow(const jau::nsize_t row, Vec4& v_out) const noexcept {
360 return v_out.set( get(row+0*4),
361 get(row+1*4),
362 get(row+2*4),
363 get(row+3*4) );
364 }
365 /**
366 * Get the named column of the given column-major matrix to v_out w/o boundary check.
367 * @param row named row to copy
368 * @return result vector holding the requested row
369 */
370 constexpr Vec4 getRow(const jau::nsize_t row) const noexcept {
371 return Vec4( get(row+0*4),
372 get(row+1*4),
373 get(row+2*4),
374 get(row+3*4) );
375 }
376
377 /**
378 * Get the named row of the given column-major matrix to v_out w/o boundary check.
379 * @param row named row to copy
380 * @param v_out the row-vector assert( i < 16 )e
381 * @return given result vector <i>v_out</i> for chaining
382 */
383 constexpr Vec3& getRow(const jau::nsize_t row, Vec3& v_out) const noexcept {
384 assert( row <= 2 );
385 return v_out.set( get(row+0*4),
386 get(row+1*4),
387 get(row+2*4) );
388 }
389
390 /**
391 * Get this matrix into the given value_type[16] array in column major order w/o boundary check.
392 *
393 * @param dst value_type[16] array storage in column major order
394 * @return {@code dst} for chaining
395 */
396 constexpr iterator get(iterator dst) const noexcept {
397 iterator dst_i = dst;
398 *dst_i = m00; // column 0
399 *(++dst_i) = m10;
400 *(++dst_i) = m20;
401 *(++dst_i) = m30;
402 *(++dst_i) = m01; // column 1
403 *(++dst_i) = m11;
404 *(++dst_i) = m21;
405 *(++dst_i) = m31;
406 *(++dst_i) = m02; // column 2
407 *(++dst_i) = m12;
408 *(++dst_i) = m22;
409 *(++dst_i) = m32;
410 *(++dst_i) = m03; // column 3
411 *(++dst_i) = m13;
412 *(++dst_i) = m23;
413 *(++dst_i) = m33;
414 return dst;
415 }
416
417 /**
418 * Get this matrix into the given {@link FloatBuffer} in column major order.
419 *
420 * @param dst 4x4 matrix std::vector in column-major order starting at {@code dst_off}
421 * @param dst_off offset for matrix {@code dst}
422 * @return {@code dst} for chaining
423 */
424 constexpr std::vector<value_type>& get(std::vector<value_type>& dst, size_t dst_off) const noexcept {
425 assert( dst.size() >= dst_off+16 && dst_off <= std::numeric_limits<size_t>::max() - 15 );
426 get( &dst[dst_off++] );
427 return dst;
428 }
429
430 //
431 // Basic matrix operations
432 //
433
434 /**
435 * Returns the determinant of this matrix
436 * @return the matrix determinant
437 */
438 value_type determinant() const noexcept {
439 value_type ret = 0;
440 ret += m00 * ( + m11*(m22*m33 - m23*m32) - m12*(m21*m33 - m23*m31) + m13*(m21*m32 - m22*m31));
441 ret -= m01 * ( + m10*(m22*m33 - m23*m32) - m12*(m20*m33 - m23*m30) + m13*(m20*m32 - m22*m30));
442 ret += m02 * ( + m10*(m21*m33 - m23*m31) - m11*(m20*m33 - m23*m30) + m13*(m20*m31 - m21*m30));
443 ret -= m03 * ( + m10*(m21*m32 - m22*m31) - m11*(m20*m32 - m22*m30) + m12*(m20*m31 - m21*m30));
444 return ret;
445 }
446
447 /**
448 * Transpose this matrix.
449 *
450 * @return this matrix for chaining
451 */
452 Matrix4& transpose() noexcept {
453 value_type tmp;
454
455 tmp = m10;
456 m10 = m01;
457 m01 = tmp;
458
459 tmp = m20;
460 m20 = m02;
461 m02 = tmp;
462
463 tmp = m30;
464 m30 = m03;
465 m03 = tmp;
466
467 tmp = m21;
468 m21 = m12;
469 m12 = tmp;
470
471 tmp = m31;
472 m31 = m13;
473 m13 = tmp;
474
475 tmp = m32;
476 m32 = m23;
477 m23 = tmp;
478
479 return *this;
480 }
481
482 /**
483 * Transpose the given {@code src} matrix into this matrix.
484 *
485 * @param src source 4x4 matrix
486 * @return this matrix (result) for chaining
487 */
488 Matrix4& transpose(const Matrix4& src) noexcept {
489 if( &src == this ) {
490 return transpose();
491 }
492 m00 = src.m00;
493 m10 = src.m01;
494 m20 = src.m02;
495 m30 = src.m03;
496
497 m01 = src.m10;
498 m11 = src.m11;
499 m21 = src.m12;
500 m31 = src.m13;
501
502 m02 = src.m20;
503 m12 = src.m21;
504 m22 = src.m22;
505 m32 = src.m23;
506
507 m03 = src.m30;
508 m13 = src.m31;
509 m23 = src.m32;
510 m33 = src.m33;
511 return *this;
512 }
513
514 /**
515 * Invert this matrix.
516 * @return false if this matrix is singular and inversion not possible, otherwise true
517 */
518 bool invert() noexcept {
519 const value_type amax = absMax();
520 if( zero == amax ) {
521 return false;
522 }
523 const value_type scale = one/amax;
524 const value_type a00 = m00*scale;
525 const value_type a10 = m10*scale;
526 const value_type a20 = m20*scale;
527 const value_type a30 = m30*scale;
528
529 const value_type a01 = m01*scale;
530 const value_type a11 = m11*scale;
531 const value_type a21 = m21*scale;
532 const value_type a31 = m31*scale;
533
534 const value_type a02 = m02*scale;
535 const value_type a12 = m12*scale;
536 const value_type a22 = m22*scale;
537 const value_type a32 = m32*scale;
538
539 const value_type a03 = m03*scale;
540 const value_type a13 = m13*scale;
541 const value_type a23 = m23*scale;
542 const value_type a33 = m33*scale;
543
544 const value_type b00 = + a11*(a22*a33 - a23*a32) - a12*(a21*a33 - a23*a31) + a13*(a21*a32 - a22*a31);
545 const value_type b01 = -( + a10*(a22*a33 - a23*a32) - a12*(a20*a33 - a23*a30) + a13*(a20*a32 - a22*a30));
546 const value_type b02 = + a10*(a21*a33 - a23*a31) - a11*(a20*a33 - a23*a30) + a13*(a20*a31 - a21*a30);
547 const value_type b03 = -( + a10*(a21*a32 - a22*a31) - a11*(a20*a32 - a22*a30) + a12*(a20*a31 - a21*a30));
548
549 const value_type b10 = -( + a01*(a22*a33 - a23*a32) - a02*(a21*a33 - a23*a31) + a03*(a21*a32 - a22*a31));
550 const value_type b11 = + a00*(a22*a33 - a23*a32) - a02*(a20*a33 - a23*a30) + a03*(a20*a32 - a22*a30);
551 const value_type b12 = -( + a00*(a21*a33 - a23*a31) - a01*(a20*a33 - a23*a30) + a03*(a20*a31 - a21*a30));
552 const value_type b13 = + a00*(a21*a32 - a22*a31) - a01*(a20*a32 - a22*a30) + a02*(a20*a31 - a21*a30);
553
554 const value_type b20 = + a01*(a12*a33 - a13*a32) - a02*(a11*a33 - a13*a31) + a03*(a11*a32 - a12*a31);
555 const value_type b21 = -( + a00*(a12*a33 - a13*a32) - a02*(a10*a33 - a13*a30) + a03*(a10*a32 - a12*a30));
556 const value_type b22 = + a00*(a11*a33 - a13*a31) - a01*(a10*a33 - a13*a30) + a03*(a10*a31 - a11*a30);
557 const value_type b23 = -( + a00*(a11*a32 - a12*a31) - a01*(a10*a32 - a12*a30) + a02*(a10*a31 - a11*a30));
558
559 const value_type b30 = -( + a01*(a12*a23 - a13*a22) - a02*(a11*a23 - a13*a21) + a03*(a11*a22 - a12*a21));
560 const value_type b31 = + a00*(a12*a23 - a13*a22) - a02*(a10*a23 - a13*a20) + a03*(a10*a22 - a12*a20);
561 const value_type b32 = -( + a00*(a11*a23 - a13*a21) - a01*(a10*a23 - a13*a20) + a03*(a10*a21 - a11*a20));
562 const value_type b33 = + a00*(a11*a22 - a12*a21) - a01*(a10*a22 - a12*a20) + a02*(a10*a21 - a11*a20);
563
564 const value_type det = (a00*b00 + a01*b01 + a02*b02 + a03*b03) / scale;
565 if( 0 == det ) {
566 return false;
567 }
568 const value_type invdet = one / det;
569
570 m00 = b00 * invdet;
571 m10 = b01 * invdet;
572 m20 = b02 * invdet;
573 m30 = b03 * invdet;
574
575 m01 = b10 * invdet;
576 m11 = b11 * invdet;
577 m21 = b12 * invdet;
578 m31 = b13 * invdet;
579
580 m02 = b20 * invdet;
581 m12 = b21 * invdet;
582 m22 = b22 * invdet;
583 m32 = b23 * invdet;
584
585 m03 = b30 * invdet;
586 m13 = b31 * invdet;
587 m23 = b32 * invdet;
588 m33 = b33 * invdet;
589 return true;
590 }
591
592 /**
593 * Invert the {@code src} matrix values into this matrix
594 * @param src the source matrix, which values are to be inverted
595 * @return false if {@code src} matrix is singular and inversion not possible, otherwise true
596 */
597 bool invert(const Matrix4& src) noexcept {
598 const value_type amax = src.absMax();
599 if( zero == amax ) {
600 return false;
601 }
602 const value_type scale = one/amax;
603 const value_type a00 = src.m00*scale;
604 const value_type a10 = src.m10*scale;
605 const value_type a20 = src.m20*scale;
606 const value_type a30 = src.m30*scale;
607
608 const value_type a01 = src.m01*scale;
609 const value_type a11 = src.m11*scale;
610 const value_type a21 = src.m21*scale;
611 const value_type a31 = src.m31*scale;
612
613 const value_type a02 = src.m02*scale;
614 const value_type a12 = src.m12*scale;
615 const value_type a22 = src.m22*scale;
616 const value_type a32 = src.m32*scale;
617
618 const value_type a03 = src.m03*scale;
619 const value_type a13 = src.m13*scale;
620 const value_type a23 = src.m23*scale;
621 const value_type a33 = src.m33*scale;
622
623 const value_type b00 = + a11*(a22*a33 - a23*a32) - a12*(a21*a33 - a23*a31) + a13*(a21*a32 - a22*a31);
624 const value_type b01 = -( + a10*(a22*a33 - a23*a32) - a12*(a20*a33 - a23*a30) + a13*(a20*a32 - a22*a30));
625 const value_type b02 = + a10*(a21*a33 - a23*a31) - a11*(a20*a33 - a23*a30) + a13*(a20*a31 - a21*a30);
626 const value_type b03 = -( + a10*(a21*a32 - a22*a31) - a11*(a20*a32 - a22*a30) + a12*(a20*a31 - a21*a30));
627
628 const value_type b10 = -( + a01*(a22*a33 - a23*a32) - a02*(a21*a33 - a23*a31) + a03*(a21*a32 - a22*a31));
629 const value_type b11 = + a00*(a22*a33 - a23*a32) - a02*(a20*a33 - a23*a30) + a03*(a20*a32 - a22*a30);
630 const value_type b12 = -( + a00*(a21*a33 - a23*a31) - a01*(a20*a33 - a23*a30) + a03*(a20*a31 - a21*a30));
631 const value_type b13 = + a00*(a21*a32 - a22*a31) - a01*(a20*a32 - a22*a30) + a02*(a20*a31 - a21*a30);
632
633 const value_type b20 = + a01*(a12*a33 - a13*a32) - a02*(a11*a33 - a13*a31) + a03*(a11*a32 - a12*a31);
634 const value_type b21 = -( + a00*(a12*a33 - a13*a32) - a02*(a10*a33 - a13*a30) + a03*(a10*a32 - a12*a30));
635 const value_type b22 = + a00*(a11*a33 - a13*a31) - a01*(a10*a33 - a13*a30) + a03*(a10*a31 - a11*a30);
636 const value_type b23 = -( + a00*(a11*a32 - a12*a31) - a01*(a10*a32 - a12*a30) + a02*(a10*a31 - a11*a30));
637
638 const value_type b30 = -( + a01*(a12*a23 - a13*a22) - a02*(a11*a23 - a13*a21) + a03*(a11*a22 - a12*a21));
639 const value_type b31 = + a00*(a12*a23 - a13*a22) - a02*(a10*a23 - a13*a20) + a03*(a10*a22 - a12*a20);
640 const value_type b32 = -( + a00*(a11*a23 - a13*a21) - a01*(a10*a23 - a13*a20) + a03*(a10*a21 - a11*a20));
641 const value_type b33 = + a00*(a11*a22 - a12*a21) - a01*(a10*a22 - a12*a20) + a02*(a10*a21 - a11*a20);
642
643 const value_type det = (a00*b00 + a01*b01 + a02*b02 + a03*b03) / scale;
644
645 if( 0 == det ) {
646 return false;
647 }
648 const value_type invdet = one / det;
649
650 m00 = b00 * invdet;
651 m10 = b01 * invdet;
652 m20 = b02 * invdet;
653 m30 = b03 * invdet;
654
655 m01 = b10 * invdet;
656 m11 = b11 * invdet;
657 m21 = b12 * invdet;
658 m31 = b13 * invdet;
659
660 m02 = b20 * invdet;
661 m12 = b21 * invdet;
662 m22 = b22 * invdet;
663 m32 = b23 * invdet;
664
665 m03 = b30 * invdet;
666 m13 = b31 * invdet;
667 m23 = b32 * invdet;
668 m33 = b33 * invdet;
669 return true;
670 }
671
672 private:
673 /** Returns the maximum abs(mxy) field */
674 value_type absMax() const noexcept {
675 value_type max = std::abs(m00);
676 max = std::max(max, std::abs(m01));
677 max = std::max(max, std::abs(m02));
678 max = std::max(max, std::abs(m03));
679
680 max = std::max(max, std::abs(m10));
681 max = std::max(max, std::abs(m11));
682 max = std::max(max, std::abs(m12));
683 max = std::max(max, std::abs(m13));
684
685 max = std::max(max, std::abs(m20));
686 max = std::max(max, std::abs(m21));
687 max = std::max(max, std::abs(m22));
688 max = std::max(max, std::abs(m23));
689
690 max = std::max(max, std::abs(m30));
691 max = std::max(max, std::abs(m31));
692 max = std::max(max, std::abs(m32));
693 max = std::max(max, std::abs(m33));
694 return max;
695 }
696
697 public:
698 /**
699 * Multiply matrix with scalar: [this] = [this] x [s]
700 * @param s a scalar
701 * @return this matrix for chaining
702 */
703 constexpr Matrix4& operator*=( const value_type s ) noexcept {
704 m00 *= s; m10 *= s; m20 *= s; m30 *= s;
705 m01 *= s; m11 *= s; m21 *= s; m31 *= s;
706 m02 *= s; m12 *= s; m22 *= s; m32 *= s;
707 m03 *= s; m13 *= s; m23 *= s; m33 *= s;
708 return *this;
709 }
710
711 /**
712 * Multiply matrix: [this] = [this] x [b]
713 * @param b 4x4 matrix
714 * @return this matrix for chaining
715 * @see #mul(mat4f, mat4f)
716 */
717 constexpr Matrix4& mul(const Matrix4& b) noexcept {
718 // return mul(new mat4f(this), b); // <- roughly half speed
719 value_type ai0=m00; // row-0, m[0+0*4]
720 value_type ai1=m01;
721 value_type ai2=m02;
722 value_type ai3=m03;
723 m00 = ai0 * b.m00 + ai1 * b.m10 + ai2 * b.m20 + ai3 * b.m30 ;
724 m01 = ai0 * b.m01 + ai1 * b.m11 + ai2 * b.m21 + ai3 * b.m31 ;
725 m02 = ai0 * b.m02 + ai1 * b.m12 + ai2 * b.m22 + ai3 * b.m32 ;
726 m03 = ai0 * b.m03 + ai1 * b.m13 + ai2 * b.m23 + ai3 * b.m33 ;
727
728 ai0=m10; //row-1, m[1+0*4]
729 ai1=m11;
730 ai2=m12;
731 ai3=m13;
732 m10 = ai0 * b.m00 + ai1 * b.m10 + ai2 * b.m20 + ai3 * b.m30 ;
733 m11 = ai0 * b.m01 + ai1 * b.m11 + ai2 * b.m21 + ai3 * b.m31 ;
734 m12 = ai0 * b.m02 + ai1 * b.m12 + ai2 * b.m22 + ai3 * b.m32 ;
735 m13 = ai0 * b.m03 + ai1 * b.m13 + ai2 * b.m23 + ai3 * b.m33 ;
736
737 ai0=m20; // row-2, m[2+0*4]
738 ai1=m21;
739 ai2=m22;
740 ai3=m23;
741 m20 = ai0 * b.m00 + ai1 * b.m10 + ai2 * b.m20 + ai3 * b.m30 ;
742 m21 = ai0 * b.m01 + ai1 * b.m11 + ai2 * b.m21 + ai3 * b.m31 ;
743 m22 = ai0 * b.m02 + ai1 * b.m12 + ai2 * b.m22 + ai3 * b.m32 ;
744 m23 = ai0 * b.m03 + ai1 * b.m13 + ai2 * b.m23 + ai3 * b.m33 ;
745
746 ai0=m30; // row-3, m[3+0*4]
747 ai1=m31;
748 ai2=m32;
749 ai3=m33;
750 m30 = ai0 * b.m00 + ai1 * b.m10 + ai2 * b.m20 + ai3 * b.m30 ;
751 m31 = ai0 * b.m01 + ai1 * b.m11 + ai2 * b.m21 + ai3 * b.m31 ;
752 m32 = ai0 * b.m02 + ai1 * b.m12 + ai2 * b.m22 + ai3 * b.m32 ;
753 m33 = ai0 * b.m03 + ai1 * b.m13 + ai2 * b.m23 + ai3 * b.m33 ;
754 return *this;
755 }
756 /**
757 * Multiply matrix: [this] = [this] x [b]
758 * @param b 4x4 matrix
759 * @return this matrix for chaining
760 * @see #mul(mat4f, mat4f)
761 */
762 constexpr Matrix4& operator*=( const Matrix4& rhs ) noexcept {
763 return mul( rhs );
764 }
765
766 /**
767 * Multiply matrix: [this] = [a] x [b]
768 * @param a 4x4 matrix, can't be this matrix
769 * @param b 4x4 matrix, can't be this matrix
770 * @return this matrix for chaining
771 * @see #mul(mat4f)
772 */
773 constexpr Matrix4& mul(const Matrix4& a, const Matrix4& b) noexcept {
774 // row-0, m[0+0*4]
775 m00 = a.m00 * b.m00 + a.m01 * b.m10 + a.m02 * b.m20 + a.m03 * b.m30 ;
776 m01 = a.m00 * b.m01 + a.m01 * b.m11 + a.m02 * b.m21 + a.m03 * b.m31 ;
777 m02 = a.m00 * b.m02 + a.m01 * b.m12 + a.m02 * b.m22 + a.m03 * b.m32 ;
778 m03 = a.m00 * b.m03 + a.m01 * b.m13 + a.m02 * b.m23 + a.m03 * b.m33 ;
779
780 //row-1, m[1+0*4]
781 m10 = a.m10 * b.m00 + a.m11 * b.m10 + a.m12 * b.m20 + a.m13 * b.m30 ;
782 m11 = a.m10 * b.m01 + a.m11 * b.m11 + a.m12 * b.m21 + a.m13 * b.m31 ;
783 m12 = a.m10 * b.m02 + a.m11 * b.m12 + a.m12 * b.m22 + a.m13 * b.m32 ;
784 m13 = a.m10 * b.m03 + a.m11 * b.m13 + a.m12 * b.m23 + a.m13 * b.m33 ;
785
786 // row-2, m[2+0*4]
787 m20 = a.m20 * b.m00 + a.m21 * b.m10 + a.m22 * b.m20 + a.m23 * b.m30 ;
788 m21 = a.m20 * b.m01 + a.m21 * b.m11 + a.m22 * b.m21 + a.m23 * b.m31 ;
789 m22 = a.m20 * b.m02 + a.m21 * b.m12 + a.m22 * b.m22 + a.m23 * b.m32 ;
790 m23 = a.m20 * b.m03 + a.m21 * b.m13 + a.m22 * b.m23 + a.m23 * b.m33 ;
791
792 // row-3, m[3+0*4]
793 m30 = a.m30 * b.m00 + a.m31 * b.m10 + a.m32 * b.m20 + a.m33 * b.m30 ;
794 m31 = a.m30 * b.m01 + a.m31 * b.m11 + a.m32 * b.m21 + a.m33 * b.m31 ;
795 m32 = a.m30 * b.m02 + a.m31 * b.m12 + a.m32 * b.m22 + a.m33 * b.m32 ;
796 m33 = a.m30 * b.m03 + a.m31 * b.m13 + a.m32 * b.m23 + a.m33 * b.m33 ;
797
798 return *this;
799 }
800
801 /**
802 * @param v_in 4-component column-vector, can be v_out for in-place transformation
803 * @param v_out this * v_in
804 * @returns v_out for chaining
805 */
806 constexpr Vec4& mulVec4(const Vec4& v_in, Vec4& v_out) const noexcept {
807 // (one matrix row in column-major order) X (column vector)
808 const value_type x = v_in.x, y = v_in.y, z = v_in.z, w = v_in.w;
809 v_out.set( x * m00 + y * m01 + z * m02 + w * m03,
810 x * m10 + y * m11 + z * m12 + w * m13,
811 x * m20 + y * m21 + z * m22 + w * m23,
812 x * m30 + y * m31 + z * m32 + w * m33 );
813 return v_out;
814 }
815
816 /**
817 * Returns new Vec4, with this * v_in
818 * @param v_in 4-component column-vector
819 */
820 constexpr Vec4 operator*(const Vec4& rhs) const noexcept {
821 // (one matrix row in column-major order) X (column vector)
822 const value_type x = rhs.x, y = rhs.y, z = rhs.z, w = rhs.w;
823 return Vec4( x * m00 + y * m01 + z * m02 + w * m03,
824 x * m10 + y * m11 + z * m12 + w * m13,
825 x * m20 + y * m21 + z * m22 + w * m23,
826 x * m30 + y * m31 + z * m32 + w * m33 );
827 }
828
829 /**
830 * @param v_inout 4-component column-vector input and output, i.e. in-place transformation
831 * @returns v_inout for chaining
832 */
833 constexpr Vec4& mulVec4(Vec4& v_inout) const noexcept {
834 // (one matrix row in column-major order) X (column vector)
835 const value_type x = v_inout.x, y = v_inout.y, z = v_inout.z, w = v_inout.w;
836 v_inout.set( x * m00 + y * m01 + z * m02 + w * m03,
837 x * m10 + y * m11 + z * m12 + w * m13,
838 x * m20 + y * m21 + z * m22 + w * m23,
839 x * m30 + y * m31 + z * m32 + w * m33 );
840 return v_inout;
841 }
842
843 /**
844 * Affine 3f-vector transformation by 4x4 matrix
845 *
846 * 4x4 matrix multiplication with 3-component vector,
847 * using {@code 1} for for {@code v_in.w} and dropping {@code v_out.w},
848 * which shall be {@code 1}.
849 *
850 * @param v_in 3-component column-vector {@link vec3f}, can be v_out for in-place transformation
851 * @param v_out m_in * v_in, 3-component column-vector {@link vec3f}
852 * @returns v_out for chaining
853 */
854 constexpr Vec3& mulVec3(const Vec3& v_in, Vec3& v_out) const noexcept {
855 // (one matrix row in column-major order) X (column vector)
856 const value_type x = v_in.x, y = v_in.y, z = v_in.z;
857 v_out.set( x * m00 + y * m01 + z * m02 + one * m03,
858 x * m10 + y * m11 + z * m12 + one * m13,
859 x * m20 + y * m21 + z * m22 + one * m23 );
860 return v_out;
861 }
862 /**
863 * Returns new Vec3, with affine 3f-vector transformation by this 4x4 matrix: this * v_in
864 *
865 * 4x4 matrix multiplication with 3-component vector,
866 * using {@code 1} for for {@code v_in.w} and dropping {@code v_out.w},
867 * which shall be {@code 1}.
868 *
869 * @param v_in 3-component column-vector {@link vec3f}
870 */
871 constexpr Vec3 operator*(const Vec3& rhs) const noexcept {
872 // (one matrix row in column-major order) X (column vector)
873 const value_type x = rhs.x, y = rhs.y, z = rhs.z;
874 return Vec3( x * m00 + y * m01 + z * m02 + one * m03,
875 x * m10 + y * m11 + z * m12 + one * m13,
876 x * m20 + y * m21 + z * m22 + one * m23 );
877 }
878
879 /**
880 * Affine 3f-vector transformation by 4x4 matrix
881 *
882 * 4x4 matrix multiplication with 3-component vector,
883 * using {@code 1} for for {@code v_inout.w} and dropping {@code v_inout.w},
884 * which shall be {@code 1}.
885 *
886 * @param v_inout 3-component column-vector {@link vec3f} input and output, i.e. in-place transformation
887 * @returns v_inout for chaining
888 */
889 constexpr Vec3& mulVec3(Vec3& v_inout) const noexcept {
890 // (one matrix row in column-major order) X (column vector)
891 const value_type x = v_inout.x, y = v_inout.y, z = v_inout.z;
892 v_inout.set( x * m00 + y * m01 + z * m02 + one * m03,
893 x * m10 + y * m11 + z * m12 + one * m13,
894 x * m20 + y * m21 + z * m22 + one * m23 );
895 return v_inout;
896 }
897
898 //
899 // Matrix setTo...(), affine + basic
900 //
901
902 /**
903 * Set this matrix to translation.
904 * <pre>
905 Translation matrix (Column Order):
906 1 0 0 0
907 0 1 0 0
908 0 0 1 0
909 x y z 1
910 * </pre>
911 * @param x x-axis translate
912 * @param y y-axis translate
913 * @param z z-axis translate
914 * @return this matrix for chaining
915 */
916 constexpr Matrix4& setToTranslation(const value_type x, const value_type y, const value_type z) noexcept {
917 m00 = m11 = m22 = m33 = one;
918 m03 = x;
919 m13 = y;
920 m23 = z;
921 m01 = m02 =
922 m10 = m12 =
923 m20 = m21 =
924 m30 = m31 = m32 = zero;
925 return *this;
926 }
927
928 /**
929 * Set this matrix to translation.
930 * <pre>
931 Translation matrix (Column Order):
932 1 0 0 0
933 0 1 0 0
934 0 0 1 0
935 x y z 1
936 * </pre>
937 * @param t translate vec3f
938 * @return this matrix for chaining
939 */
940 constexpr Matrix4& setToTranslation(const Vec3& t) noexcept {
941 return setToTranslation(t.x, t.y, t.z);
942 }
943
944 /**
945 * Set this matrix to scale.
946 * <pre>
947 Scale matrix (Any Order):
948 x 0 0 0
949 0 y 0 0
950 0 0 z 0
951 0 0 0 1
952 * </pre>
953 * @param x x-axis scale
954 * @param y y-axis scale
955 * @param z z-axis scale
956 * @return this matrix for chaining
957 */
958 constexpr Matrix4& setToScale(const value_type x, const value_type y, const value_type z) noexcept {
959 m33 = one;
960 m00 = x;
961 m11 = y;
962 m22 = z;
963 m01 = m02 = m03 =
964 m10 = m12 = m13 =
965 m20 = m21 = m23 =
966 m30 = m31 = m32 = zero;
967 return *this;
968 }
969
970 /**
971 * Set this matrix to scale.
972 * <pre>
973 Scale matrix (Any Order):
974 x 0 0 0
975 0 y 0 0
976 0 0 z 0
977 0 0 0 1
978 * </pre>
979 * @param s scale vec3f
980 * @return this matrix for chaining
981 */
982 constexpr Matrix4& setToScale(const Vec3& s) noexcept {
983 return setToScale(s.x, s.y, s.z);
984 }
985
986 /**
987 * Set this matrix to rotation from the given axis and angle in radians.
988 * <pre>
989 Rotation matrix (Column Order):
990 xx(1-c)+c xy(1-c)+zs xz(1-c)-ys 0
991 xy(1-c)-zs yy(1-c)+c yz(1-c)+xs 0
992 xz(1-c)+ys yz(1-c)-xs zz(1-c)+c 0
993 0 0 0 1
994 * </pre>
995 * @see <a href="http://web.archive.org/web/20041029003853/http://www.j3d.org/matrix_faq/matrfaq_latest.html#Q38">Matrix-FAQ Q38</a>
996 * @param ang_rad angle in radians
997 * @param x x of rotation axis
998 * @param y y of rotation axis
999 * @param z z of rotation axis
1000 * @return this matrix for chaining
1001 */
1003 const value_type c = std::cos(ang_rad);
1004 const value_type ic= one - c;
1005 const value_type s = std::sin(ang_rad);
1006
1007 Vec3 tmp(x, y, z); tmp.normalize();
1008 x = tmp.x; y = tmp.y; z = tmp.z;
1009
1010 const value_type xy = x*y;
1011 const value_type xz = x*z;
1012 const value_type xs = x*s;
1013 const value_type ys = y*s;
1014 const value_type yz = y*z;
1015 const value_type zs = z*s;
1016 m00 = x*x*ic+c;
1017 m10 = xy*ic+zs;
1018 m20 = xz*ic-ys;
1019 m30 = zero;
1020
1021 m01 = xy*ic-zs;
1022 m11 = y*y*ic+c;
1023 m21 = yz*ic+xs;
1024 m31 = zero;
1025
1026 m02 = xz*ic+ys;
1027 m12 = yz*ic-xs;
1028 m22 = z*z*ic+c;
1029 m32 = zero;
1030
1031 m03 = zero;
1032 m13 = zero;
1033 m23 = zero;
1034 m33 = one;
1035
1036 return *this;
1037 }
1038
1039 /**
1040 * Set this matrix to rotation from the given axis and angle in radians.
1041 * <pre>
1042 Rotation matrix (Column Order):
1043 xx(1-c)+c xy(1-c)+zs xz(1-c)-ys 0
1044 xy(1-c)-zs yy(1-c)+c yz(1-c)+xs 0
1045 xz(1-c)+ys yz(1-c)-xs zz(1-c)+c 0
1046 0 0 0 1
1047 * </pre>
1048 * @see <a href="http://web.archive.org/web/20041029003853/http://www.j3d.org/matrix_faq/matrfaq_latest.html#Q38">Matrix-FAQ Q38</a>
1049 * @param ang_rad angle in radians
1050 * @param axis rotation axis
1051 * @return this matrix for chaining
1052 */
1053 constexpr_cxx26 Matrix4& setToRotationAxis(const value_type ang_rad, const Vec3& axis) noexcept {
1054 return setToRotationAxis(ang_rad, axis.x, axis.y, axis.z);
1055 }
1056
1057 /**
1058 * Set this matrix to rotation from the given Euler rotation angles in radians.
1059 * <p>
1060 * The rotations are applied in the given order:
1061 * <ul>
1062 * <li>y - heading</li>
1063 * <li>z - attitude</li>
1064 * <li>x - bank</li>
1065 * </ul>
1066 * </p>
1067 * @param bankX the Euler pitch angle in radians. (rotation about the X axis)
1068 * @param headingY the Euler yaw angle in radians. (rotation about the Y axis)
1069 * @param attitudeZ the Euler roll angle in radians. (rotation about the Z axis)
1070 * @return this matrix for chaining
1071 * <p>
1072 * Implementation does not use Quaternion and hence is exposed to
1073 * <a href="http://web.archive.org/web/20041029003853/http://www.j3d.org/matrix_faq/matrfaq_latest.html#Q34">Gimbal-Lock</a>,
1074 * consider using Quaternion::toMatrix().
1075 * </p>
1076 * @see <a href="http://web.archive.org/web/20041029003853/http://www.j3d.org/matrix_faq/matrfaq_latest.html#Q36">Matrix-FAQ Q36</a>
1077 * @see <a href="http://www.euclideanspace.com/maths/geometry/rotations/conversions/eulerToMatrix/index.htm">euclideanspace.com-eulerToMatrix</a>
1078 * @see Quaternion::toMatrix()
1079 */
1080 constexpr_cxx26 Matrix4& setToRotationEuler(const value_type bankX, const value_type headingY, const value_type attitudeZ) noexcept {
1081 // Assuming the angles are in radians.
1082 const value_type ch = std::cos(headingY);
1083 const value_type sh = std::sin(headingY);
1084 const value_type ca = std::cos(attitudeZ);
1085 const value_type sa = std::sin(attitudeZ);
1086 const value_type cb = std::cos(bankX);
1087 const value_type sb = std::sin(bankX);
1088
1089 m00 = ch*ca;
1090 m10 = sa;
1091 m20 = -sh*ca;
1092 m30 = zero;
1093
1094 m01 = sh*sb - ch*sa*cb;
1095 m11 = ca*cb;
1096 m21 = sh*sa*cb + ch*sb;
1097 m31 = zero;
1098
1099 m02 = ch*sa*sb + sh*cb;
1100 m12 = -ca*sb;
1101 m22 = -sh*sa*sb + ch*cb;
1102 m32 = zero;
1103
1104 m03 = zero;
1105 m13 = zero;
1106 m23 = zero;
1107 m33 = one;
1108
1109 return *this;
1110 }
1111
1112 /**
1113 * Set this matrix to rotation from the given Euler rotation angles in radians.
1114 * <p>
1115 * The rotations are applied in the given order:
1116 * <ul>
1117 * <li>y - heading</li>
1118 * <li>z - attitude</li>
1119 * <li>x - bank</li>
1120 * </ul>
1121 * </p>
1122 * @param angradXYZ euler angle vector in radians holding x-bank, y-heading and z-attitude
1123 * @return this quaternion for chaining.
1124 * <p>
1125 * Implementation does not use Quaternion and hence is exposed to
1126 * <a href="http://web.archive.org/web/20041029003853/http://www.j3d.org/matrix_faq/matrfaq_latest.html#Q34">Gimbal-Lock</a>,
1127 * consider using Quaternion::toMatrix().
1128 * </p>
1129 * @see <a href="http://web.archive.org/web/20041029003853/http://www.j3d.org/matrix_faq/matrfaq_latest.html#Q36">Matrix-FAQ Q36</a>
1130 * @see <a href="http://www.euclideanspace.com/maths/geometry/rotations/conversions/eulerToMatrix/index.htm">euclideanspace.com-eulerToMatrix</a>
1131 * @see Quaternion::toMatrix()
1132 */
1133 constexpr_cxx26 Matrix4& setToRotationEuler(const Vec3& angradXYZ) noexcept {
1134 return setToRotationEuler(angradXYZ.x, angradXYZ.y, angradXYZ.z);
1135 }
1136
1137 /**
1138 * Set this matrix to orthogonal projection.
1139 * <pre>
1140 Ortho matrix (Column Order):
1141 2/dx 0 0 0
1142 0 2/dy 0 0
1143 0 0 2/dz 0
1144 tx ty tz 1
1145 * </pre>
1146 * @param left
1147 * @param right
1148 * @param bottom
1149 * @param top
1150 * @param zNear
1151 * @param zFar
1152 * @return this matrix for chaining
1153 */
1154 constexpr Matrix4& setToOrtho(const value_type left, const value_type right,
1155 const value_type bottom, const value_type top,
1156 const value_type zNear, const value_type zFar) noexcept {
1157 {
1158 // m00 = m11 = m22 = m33 = one;
1159 m10 = m20 = m30 = zero;
1160 m01 = m21 = m31 = zero;
1161 m02 = m12 = m32 = zero;
1162 // m03 = m13 = m23 = zero;
1163 }
1164 const value_type dx=right-left;
1165 const value_type dy=top-bottom;
1166 const value_type dz=zFar-zNear;
1167 const value_type tx=-one*(right+left)/dx;
1168 const value_type ty=-one*(top+bottom)/dy;
1169 const value_type tz=-one*(zFar+zNear)/dz;
1170
1171 m00 = two/dx;
1172 m11 = two/dy;
1173 m22 = -two/dz;
1174
1175 m03 = tx;
1176 m13 = ty;
1177 m23 = tz;
1178 m33 = one;
1179
1180 return *this;
1181 }
1182
1183 /**
1184 * Set this matrix to frustum.
1185 * <pre>
1186 Frustum matrix (Column Order):
1187 2*zNear/dx 0 0 0
1188 0 2*zNear/dy 0 0
1189 A B C -1
1190 0 0 D 0
1191 * </pre>
1192 * @param left
1193 * @param right
1194 * @param bottom
1195 * @param top
1196 * @param zNear
1197 * @param zFar
1198 * @return this matrix for chaining
1199 * @throws IllegalArgumentException if {@code zNear <= 0} or {@code zFar <= zNear}
1200 * or {@code left == right}, or {@code bottom == top}.
1201 */
1202 Matrix4& setToFrustum(const value_type left, const value_type right,
1203 const value_type bottom, const value_type top,
1204 const value_type zNear, const value_type zFar) {
1205 if( zNear <= zero || zFar <= zNear ) {
1206 throw jau::IllegalArgumentException("Requirements zNear > 0 and zFar > zNear, but zNear "+std::to_string(zNear)+", zFar "+std::to_string(zFar), E_FILE_LINE);
1207 }
1208 if( left == right || top == bottom) {
1209 throw jau::IllegalArgumentException("GL_INVALID_VALUE: top,bottom and left,right must not be equal", E_FILE_LINE);
1210 }
1211 {
1212 // m00 = m11 = m22 = m33 = 1f;
1213 m10 = m20 = m30 = zero;
1214 m01 = m21 = m31 = zero;
1215 m03 = m13 = zero;
1216 }
1217 const value_type zNear2 = two*zNear;
1218 const value_type dx=right-left;
1219 const value_type dy=top-bottom;
1220 const value_type dz=zFar-zNear;
1221 const value_type A=(right+left)/dx;
1222 const value_type B=(top+bottom)/dy;
1223 const value_type C=-one*(zFar+zNear)/dz;
1224 const value_type D=-two*(zFar*zNear)/dz;
1225
1226 m00 = zNear2/dx;
1227 m11 = zNear2/dy;
1228
1229 m02 = A;
1230 m12 = B;
1231 m22 = C;
1232 m32 = -one;
1233
1234 m23 = D;
1235 m33 = zero;
1236
1237 return *this;
1238 }
1239
1240 /**
1241 * Set this matrix to perspective {@link #setToFrustum(value_type, value_type, value_type, value_type, value_type, value_type) frustum} projection.
1242 *
1243 * @param fovy_rad angle in radians
1244 * @param aspect aspect ratio width / height
1245 * @param zNear
1246 * @param zFar
1247 * @return this matrix for chaining
1248 * @throws IllegalArgumentException if {@code zNear <= 0} or {@code zFar <= zNear}
1249 * @see #setToFrustum(value_type, value_type, value_type, value_type, value_type, value_type)
1250 */
1251 Matrix4& setToPerspective(const value_type fovy_rad, const value_type aspect, const value_type zNear, const value_type zFar) {
1252 const value_type top = std::tan(fovy_rad/two) * zNear; // use tangent of half-fov !
1253 const value_type bottom = -one * top; // -1f * fovhvTan.top * zNear
1254 const value_type left = aspect * bottom; // aspect * -1f * fovhvTan.top * zNear
1255 const value_type right = aspect * top; // aspect * fovhvTan.top * zNear
1256 return setToFrustum(left, right, bottom, top, zNear, zFar);
1257 }
1258
1259 /**
1260 * Set this matrix to perspective {@link #setToFrustum(value_type, value_type, value_type, value_type, value_type, value_type) frustum} projection.
1261 *
1262 * @param fovhv {@link FovHVHalves} field of view in both directions, may not be centered, either in radians or tangent
1263 * @param zNear
1264 * @param zFar
1265 * @return this matrix for chaining
1266 * @throws IllegalArgumentException if {@code zNear <= 0} or {@code zFar <= zNear}
1267 * @see #setToFrustum(value_type, value_type, value_type, value_type, value_type, value_type)
1268 * @see Frustum#updateByFovDesc(mat4f, com.jogamp.math.geom.Frustum.FovDesc)
1269 */
1270 Matrix4& setToPerspective(const FovHVHalves& fovhv, const value_type zNear, const value_type zFar) {
1271 const FovHVHalves fovhvTan = fovhv.toTangents(); // use tangent of half-fov !
1272 const value_type top = fovhvTan.top * zNear;
1273 const value_type bottom = -one * fovhvTan.bottom * zNear;
1274 const value_type left = -one * fovhvTan.left * zNear;
1275 const value_type right = fovhvTan.right * zNear;
1276 return setToFrustum(left, right, bottom, top, zNear, zFar);
1277 }
1278
1279 /**
1280 * Set this matrix to the <i>look-at</i> matrix based on given parameters.
1281 * <p>
1282 * Consist out of two matrix multiplications:
1283 * <pre>
1284 * <b>R</b> = <b>L</b> x <b>T</b>,
1285 * with <b>L</b> for <i>look-at</i> matrix and
1286 * <b>T</b> for eye translation.
1287 *
1288 * Result <b>R</b> can be utilized for <i>projection or modelview</i> multiplication, i.e.
1289 * <b>M</b> = <b>M</b> x <b>R</b>,
1290 * with <b>M</b> being the <i>projection or modelview</i> matrix.
1291 * </pre>
1292 * </p>
1293 * @param eye 3 component eye vector
1294 * @param center 3 component center vector
1295 * @param up 3 component up vector
1296 * @param tmp temporary mat4f used for multiplication
1297 * @return this matrix for chaining
1298 */
1299 constexpr Matrix4& setToLookAt(const Vec3& eye, const Vec3& center, const Vec3& up, Matrix4& tmp) noexcept {
1300 // normalized forward!
1301 const Vec3 fwd = ( center - eye ).normalize();
1302
1303 /* Side = forward x up, normalized */
1304 const Vec3 side = fwd.cross(up).normalize();
1305
1306 /* Recompute up as: up = side x forward */
1307 const Vec3 up2 = side.cross(fwd);
1308
1309 m00 = side.x;
1310 m10 = up2.x;
1311 m20 = -fwd.x;
1312 m30 = 0;
1313
1314 m01 = side.y;
1315 m11 = up2.y;
1316 m21 = -fwd.y;
1317 m31 = 0;
1318
1319 m02 = side.z;
1320 m12 = up2.z;
1321 m22 = -fwd.z;
1322 m32 = 0;
1323
1324 m03 = 0;
1325 m13 = 0;
1326 m23 = 0;
1327 m33 = 1;
1328
1329 return mul( tmp.setToTranslation( -eye.x, -eye.y, -eye.z ) );
1330 }
1331
1332 /**
1333 * Set this matrix to the <i>pick</i> matrix based on given parameters.
1334 * <p>
1335 * Traditional <code>gluPickMatrix</code> implementation.
1336 * </p>
1337 * <p>
1338 * Consist out of two matrix multiplications:
1339 * <pre>
1340 * <b>R</b> = <b>T</b> x <b>S</b>,
1341 * with <b>T</b> for viewport translation matrix and
1342 * <b>S</b> for viewport scale matrix.
1343 *
1344 * Result <b>R</b> can be utilized for <i>projection</i> multiplication, i.e.
1345 * <b>P</b> = <b>P</b> x <b>R</b>,
1346 * with <b>P</b> being the <i>projection</i> matrix.
1347 * </pre>
1348 * </p>
1349 * <p>
1350 * To effectively use the generated pick matrix for picking,
1351 * call {@link #setToPick(value_type, value_type, value_type, value_type, Recti, mat4f) setToPick(..)}
1352 * and multiply a {@link #setToPerspective(value_type, value_type, value_type, value_type) custom perspective matrix}
1353 * by this pick matrix. Then you may load the result onto the perspective matrix stack.
1354 * </p>
1355 * @param x the center x-component of a picking region in window coordinates
1356 * @param y the center y-component of a picking region in window coordinates
1357 * @param deltaX the width of the picking region in window coordinates.
1358 * @param deltaY the height of the picking region in window coordinates.
1359 * @param viewport Rect4i viewport
1360 * @param mat4Tmp temp storage
1361 * @return true if successful or false if either delta value is <= zero.
1362 */
1363 constexpr bool setToPick(const value_type x, const value_type y, const value_type deltaX, const value_type deltaY,
1364 const Recti& viewport, Matrix4& mat4Tmp) noexcept {
1365 if (deltaX <= 0 || deltaY <= 0) {
1366 return false;
1367 }
1368 /* Translate and scale the picked region to the entire window */
1369 setToTranslation( ( viewport.width() - two * ( x - viewport.x() ) ) / deltaX,
1370 ( viewport.height() - two * ( y - viewport.y() ) ) / deltaY,
1371 0);
1372 mat4Tmp.setToScale( viewport.width() / deltaX, viewport.height() / deltaY, one );
1373 mul(mat4Tmp);
1374 return true;
1375 }
1376
1377 //
1378 // Matrix affine operations using setTo..()
1379 //
1380
1381 /**
1382 * Rotate this matrix about give axis and angle in radians, i.e. multiply by {@link #setToRotationAxis(value_type, value_type, value_type, value_type) axis-rotation matrix}.
1383 * @see <a href="http://web.archive.org/web/20041029003853/http://www.j3d.org/matrix_faq/matrfaq_latest.html#Q38">Matrix-FAQ Q38</a>
1384 * @param angrad angle in radians
1385 * @param x x of rotation axis
1386 * @param y y of rotation axis
1387 * @param z z of rotation axis
1388 * @param tmp temporary mat4f used for multiplication
1389 * @return this matrix for chaining
1390 */
1391 constexpr_cxx26 Matrix4& rotate(const value_type ang_rad, const value_type x, const value_type y, const value_type z, Matrix4& tmp) noexcept {
1392 return mul( tmp.setToRotationAxis(ang_rad, x, y, z) );
1393 }
1394
1395 /**
1396 * Rotate this matrix about give axis and angle in radians, i.e. multiply by {@link #setToRotationAxis(value_type, vec3f) axis-rotation matrix}.
1397 * @see <a href="http://web.archive.org/web/20041029003853/http://www.j3d.org/matrix_faq/matrfaq_latest.html#Q38">Matrix-FAQ Q38</a>
1398 * @param angrad angle in radians
1399 * @param axis rotation axis
1400 * @param tmp temporary mat4f used for multiplication
1401 * @return this matrix for chaining
1402 */
1403 constexpr_cxx26 Matrix4& rotate(const value_type ang_rad, const Vec3& axis, Matrix4& tmp) noexcept {
1404 return mul( tmp.setToRotationAxis(ang_rad, axis) );
1405 }
1406
1407 /**
1408 * Translate this matrix, i.e. multiply by {@link #setToTranslation(value_type, value_type, value_type) translation matrix}.
1409 * @param x x translation
1410 * @param y y translation
1411 * @param z z translation
1412 * @param tmp temporary mat4f used for multiplication
1413 * @return this matrix for chaining
1414 */
1415 constexpr Matrix4& translate(const value_type x, const value_type y, const value_type z, Matrix4& tmp) noexcept {
1416 return mul( tmp.setToTranslation(x, y, z) );
1417 }
1418
1419 /**
1420 * Translate this matrix, i.e. multiply by {@link #setToTranslation(vec3f) translation matrix}.
1421 * @param t translation vec3f
1422 * @param tmp temporary mat4f used for multiplication
1423 * @return this matrix for chaining
1424 */
1425 constexpr Matrix4& translate(const Vec3& t, Matrix4& tmp) noexcept {
1426 return mul( tmp.setToTranslation(t) );
1427 }
1428
1429 /**
1430 * Scale this matrix, i.e. multiply by {@link #setToScale(value_type, value_type, value_type) scale matrix}.
1431 * @param x x scale
1432 * @param y y scale
1433 * @param z z scale
1434 * @param tmp temporary mat4f used for multiplication
1435 * @return this matrix for chaining
1436 */
1437 constexpr Matrix4& scale(const value_type x, const value_type y, const value_type z, Matrix4& tmp) noexcept {
1438 return mul( tmp.setToScale(x, y, z) );
1439 }
1440
1441 /**
1442 * Scale this matrix, i.e. multiply by {@link #setToScale(value_type, value_type, value_type) scale matrix}.
1443 * @param s scale for x-, y- and z-axis
1444 * @param tmp temporary mat4f used for multiplication
1445 * @return this matrix for chaining
1446 */
1447 constexpr Matrix4& scale(const value_type s, Matrix4& tmp) noexcept {
1448 return mul( tmp.setToScale(s, s, s) );
1449 }
1450
1451 //
1452 // Static multi Matrix ops
1453 //
1454
1455 /**
1456 * Map object coordinates to window coordinates.
1457 * <p>
1458 * Traditional <code>gluProject</code> implementation.
1459 * </p>
1460 *
1461 * @param obj object position, 3 component vector
1462 * @param mMv modelview matrix
1463 * @param mP projection matrix
1464 * @param viewport Rect4i viewport
1465 * @param winPos 3 component window coordinate, the result
1466 * @return true if successful, otherwise false (z is 1)
1467 */
1468 static bool mapObjToWin(const Vec3& obj, const Matrix4& mMv, const Matrix4& mP,
1469 const Recti& viewport, Vec3& winPos) noexcept
1470 {
1471 // vec4Tmp2 = Mv * o
1472 // rawWinPos = P * vec4Tmp2
1473 // rawWinPos = P * ( Mv * o )
1474 // rawWinPos = P * Mv * o
1475 Vec4 vec4Tmp2 = mMv * Vec4(obj, 1.0f);
1476
1477 Vec4 rawWinPos = mP * vec4Tmp2;
1478
1479 if ( zero == rawWinPos.w ) {
1480 return false;
1481 }
1482
1483 const value_type s = ( one / rawWinPos.w ) * half;
1484
1485 // Map x, y and z to range 0-1 (w is ignored)
1486 rawWinPos.scale(s).add(half, half, half, 0.0f);
1487
1488 // Map x,y to viewport
1489 winPos.set( rawWinPos.x * viewport.width() + viewport.x(),
1490 rawWinPos.y * viewport.height() + viewport.y(),
1491 rawWinPos.z );
1492
1493 return true;
1494 }
1495
1496 /**
1497 * Map object coordinates to window coordinates.
1498 * <p>
1499 * Traditional <code>gluProject</code> implementation.
1500 * </p>
1501 *
1502 * @param obj object position, 3 component vector
1503 * @param mPMv [projection] x [modelview] matrix, i.e. P x Mv
1504 * @param viewport Rect4i viewport
1505 * @param winPos 3 component window coordinate, the result
1506 * @return true if successful, otherwise false (z is 1)
1507 */
1508 static bool mapObjToWin(const Vec3& obj, const Matrix4& mPMv,
1509 const Recti& viewport, Vec3& winPos) noexcept
1510 {
1511 // rawWinPos = P * Mv * o
1512 Vec4 rawWinPos = mPMv * Vec4(obj, 1);
1513
1514 if ( zero == rawWinPos.w ) {
1515 return false;
1516 }
1517
1518 const value_type s = ( one / rawWinPos.w ) * half;
1519
1520 // Map x, y and z to range 0-1 (w is ignored)
1521 rawWinPos.scale(s).add(half, half, half, 0.0f);
1522
1523 // Map x,y to viewport
1524 winPos.set( rawWinPos.x * viewport.width() + viewport.x(),
1525 rawWinPos.y * viewport.height() + viewport.y(),
1526 rawWinPos.z );
1527
1528 return true;
1529 }
1530
1531 /**
1532 * Map window coordinates to object coordinates.
1533 * <p>
1534 * Traditional <code>gluUnProject</code> implementation.
1535 * </p>
1536 *
1537 * @param winx
1538 * @param winy
1539 * @param winz
1540 * @param mMv 4x4 modelview matrix
1541 * @param mP 4x4 projection matrix
1542 * @param viewport Rect4i viewport
1543 * @param objPos 3 component object coordinate, the result
1544 * @param mat4Tmp 16 component matrix for temp storage
1545 * @return true if successful, otherwise false (failed to invert matrix, or becomes infinity due to zero z)
1546 */
1547 static bool mapWinToObj(const value_type winx, const value_type winy, const value_type winz,
1548 const Matrix4& mMv, const Matrix4& mP,
1549 const Recti& viewport,
1550 Vec3& objPos,
1551 Matrix4& mat4Tmp) noexcept
1552 {
1553 // invPMv = Inv(P x Mv)
1554 Matrix4& invPMv = mat4Tmp.mul(mP, mMv);
1555 if( !invPMv.invert() ) {
1556 return false;
1557 }
1558
1559 Vec4 winPos(winx, winy, winz, 1.0f);
1560
1561 // Map x and y from window coordinates
1562 winPos.add(-viewport.x(), -viewport.y(), 0.0f, 0.0f).mul(1.0f/viewport.width(), 1.0f/viewport.height(), 1.0f, 1.0f);
1563
1564 // Map to range -1 to 1
1565 winPos.mul(2.0f, 2.0f, 2.0f, 1.0f).add(-1.0f, -1.0f, -1.0f, 0.0f);
1566
1567 // rawObjPos = Inv(P x Mv) * winPos
1568 Vec4 rawObjPos = invPMv * winPos;
1569
1570 if ( zero == rawObjPos.w ) {
1571 return false;
1572 }
1573
1574 rawObjPos.scale(1.0f / rawObjPos.w).getVec3(objPos);
1575 return true;
1576 }
1577
1578 /**
1579 * Map window coordinates to object coordinates.
1580 * <p>
1581 * Traditional <code>gluUnProject</code> implementation.
1582 * </p>
1583 *
1584 * @param winx
1585 * @param winy
1586 * @param winz
1587 * @param invPMv inverse [projection] x [modelview] matrix, i.e. Inv(P x Mv), if null method returns false
1588 * @param viewport Rect4i viewport
1589 * @param objPos 3 component object coordinate, the result
1590 * @return true if successful, otherwise false (null invert matrix, or becomes infinity due to zero z)
1591 */
1592 static bool mapWinToObj(const value_type winx, const value_type winy, const value_type winz,
1593 const Matrix4& invPMv,
1594 const Recti& viewport,
1595 Vec3& objPos) noexcept
1596 {
1597 Vec4 winPos(winx, winy, winz, 1.0f);
1598
1599 // Map x and y from window coordinates
1600 winPos.add(-viewport.x(), -viewport.y(), 0.0f, 0.0f).mul(1.0f/viewport.width(), 1.0f/viewport.height(), 1.0f, 1.0f);
1601
1602 // Map to range -1 to 1
1603 winPos.mul(2.0f, 2.0f, 2.0f, 1.0f).add(-1.0f, -1.0f, -1.0f, 0.0f);
1604
1605 // rawObjPos = Inv(P x Mv) * winPos
1606 Vec4 rawObjPos = invPMv * winPos;
1607
1608 if ( zero == rawObjPos.w ) {
1609 return false;
1610 }
1611
1612 rawObjPos.scale(1.0f / rawObjPos.w).getVec3(objPos);
1613 return true;
1614 }
1615
1616 /**
1617 * Map two window coordinates to two object coordinates,
1618 * distinguished by their z component.
1619 * <p>
1620 * Traditional <code>gluUnProject</code> implementation.
1621 * </p>
1622 *
1623 * @param winx
1624 * @param winy
1625 * @param winz1
1626 * @param winz2
1627 * @param invPMv inverse [projection] x [modelview] matrix, i.e. Inv(P x Mv), if null method returns false
1628 * @param viewport Rect4i viewport vector
1629 * @param objPos1 3 component object coordinate, the result
1630 * @return true if successful, otherwise false (null invert matrix, or becomes infinity due to zero z)
1631 */
1632 static bool mapWinToObj(const value_type winx, const value_type winy, const value_type winz1, const value_type winz2,
1633 const Matrix4& invPMv,
1634 const Recti& viewport,
1635 Vec3& objPos1, Vec3& objPos2) noexcept
1636 {
1637 Vec4 winPos(winx, winy, winz1, 1.0f);
1638
1639 // Map x and y from window coordinates
1640 winPos.add(-viewport.x(), -viewport.y(), 0.0f, 0.0f).mul(1.0f/viewport.width(), 1.0f/viewport.height(), 1.0f, 1.0f);
1641
1642 // Map to range -1 to 1
1643 winPos.mul(2.0f, 2.0f, 2.0f, 1.0f).add(-1.0f, -1.0f, -1.0f, 0.0f);
1644
1645 // rawObjPos = Inv(P x Mv) * winPos1
1646 Vec4 rawObjPos = invPMv * winPos;
1647
1648 if ( zero == rawObjPos.w ) {
1649 return false;
1650 }
1651 rawObjPos.scale(1.0f / rawObjPos.w).getVec3(objPos1);
1652
1653 //
1654 // winz2
1655 //
1656 // Map Z to range -1 to 1
1657 winPos.z = winz2 * 2.0f - 1.0f;
1658
1659 // rawObjPos = Inv(P x Mv) * winPos2
1660 invPMv.mulVec4(winPos, rawObjPos);
1661
1662 if ( zero == rawObjPos.w ) {
1663 return false;
1664 }
1665 rawObjPos.scale(1.0f / rawObjPos.w).getVec3(objPos2);
1666
1667 return true;
1668 }
1669
1670 /**
1671 * Map window coordinates to object coordinates.
1672 * <p>
1673 * Traditional <code>gluUnProject4</code> implementation.
1674 * </p>
1675 *
1676 * @param winx
1677 * @param winy
1678 * @param winz
1679 * @param clipw
1680 * @param mMv 4x4 modelview matrix
1681 * @param mP 4x4 projection matrix
1682 * @param viewport Rect4i viewport vector
1683 * @param near
1684 * @param far
1685 * @param obj_pos 4 component object coordinate, the result
1686 * @param mat4Tmp 16 component matrix for temp storage
1687 * @return true if successful, otherwise false (failed to invert matrix, or becomes infinity due to zero z)
1688 */
1689 static bool mapWinToObj4(const value_type winx, const value_type winy, const value_type winz, const value_type clipw,
1690 const Matrix4& mMv, const Matrix4& mP,
1691 const Recti& viewport,
1692 const value_type near, const value_type far,
1693 Vec4& objPos,
1694 Matrix4& mat4Tmp) noexcept
1695 {
1696 // invPMv = Inv(P x Mv)
1697 Matrix4& invPMv = mat4Tmp.mul(mP, mMv);
1698 if( !invPMv.invert() ) {
1699 return false;
1700 }
1701 Vec4 winPos(winx, winy, winz, clipw);
1702
1703 // Map x and y from window coordinates
1704 winPos.add(-viewport.x(), -viewport.y(), -near, 0.0f).mul(1.0f/viewport.width(), 1.0f/viewport.height(), 1.0f/(far-near), 1.0f);
1705
1706 // Map to range -1 to 1
1707 winPos.mul(2.0f, 2.0f, 2.0f, 1.0f).add(-1.0f, -1.0f, -1.0f, 0.0f);
1708
1709 // objPos = Inv(P x Mv) * winPos
1710 invPMv.mulVec4(winPos, objPos);
1711
1712 if ( zero == objPos.w ) {
1713 return false;
1714 }
1715 return true;
1716 }
1717
1718 /**
1719 * Map window coordinates to object coordinates.
1720 * <p>
1721 * Traditional <code>gluUnProject4</code> implementation.
1722 * </p>
1723 *
1724 * @param winx
1725 * @param winy
1726 * @param winz
1727 * @param clipw
1728 * @param invPMv inverse [projection] x [modelview] matrix, i.e. Inv(P x Mv), if null method returns false
1729 * @param viewport Rect4i viewport vector
1730 * @param near
1731 * @param far
1732 * @param obj_pos 4 component object coordinate, the result
1733 * @return true if successful, otherwise false (null invert matrix, or becomes infinity due to zero z)
1734 */
1735 static bool mapWinToObj4(const value_type winx, const value_type winy, const value_type winz, const value_type clipw,
1736 const Matrix4& invPMv,
1737 const Recti& viewport,
1738 const value_type near, const value_type far,
1739 Vec4& objPos) noexcept
1740 {
1741 Vec4 winPos(winx, winy, winz, clipw);
1742
1743 // Map x and y from window coordinates
1744 winPos.add(-viewport.x(), -viewport.y(), -near, 0.0f).mul(1.0f/viewport.width(), 1.0f/viewport.height(), 1.0f/(far-near), 1.0f);
1745
1746 // Map to range -1 to 1
1747 winPos.mul(2.0f, 2.0f, 2.0f, 1.0f).add(-1.0f, -1.0f, -1.0f, 0.0f);
1748
1749 // objPos = Inv(P x Mv) * winPos
1750 invPMv.mulVec4(winPos, objPos);
1751
1752 if ( zero == objPos.w ) {
1753 return false;
1754 }
1755 return true;
1756 }
1757
1758 /**
1759 * Map two window coordinates w/ shared X/Y and distinctive Z
1760 * to a {@link Ray}. The resulting {@link Ray} maybe used for <i>picking</i>
1761 * using a {@link AABBox#getRayIntersection(vec3f, Ray, value_type, boolean)}.
1762 * <p>
1763 * Notes for picking <i>winz0</i> and <i>winz1</i>:
1764 * <ul>
1765 * <li>see {@link FloatUtil#getZBufferEpsilon(int, value_type, value_type)}</li>
1766 * <li>see {@link FloatUtil#getZBufferValue(int, value_type, value_type, value_type)}</li>
1767 * <li>see {@link FloatUtil#getOrthoWinZ(value_type, value_type, value_type)}</li>
1768 * </ul>
1769 * </p>
1770 * @param winx
1771 * @param winy
1772 * @param winz0
1773 * @param winz1
1774 * @param mMv 4x4 modelview matrix
1775 * @param mP 4x4 projection matrix
1776 * @param viewport Rect4i viewport
1777 * @param ray storage for the resulting {@link Ray}
1778 * @param mat4Tmp1 16 component matrix for temp storage
1779 * @param mat4Tmp2 16 component matrix for temp storage
1780 * @return true if successful, otherwise false (failed to invert matrix, or becomes z is infinity)
1781 */
1782 static bool mapWinToRay(const value_type winx, const value_type winy, const value_type winz0, const value_type winz1,
1783 const Matrix4& mMv, const Matrix4& mP,
1784 const Recti& viewport,
1785 Ray3& ray,
1786 Matrix4& mat4Tmp1) noexcept
1787 {
1788 // invPMv = Inv(P x Mv)
1789 const Matrix4 invPMv = mat4Tmp1.mul(mP, mMv);
1790 if( !invPMv.invert() ) {
1791 return false;
1792 }
1793
1794 if( mapWinToObj(winx, winy, winz0, winz1, invPMv, viewport, ray.orig, ray.dir) ) {
1795 ray.dir.sub(ray.orig).normalize();
1796 return true;
1797 } else {
1798 return false;
1799 }
1800 }
1801
1802 /**
1803 * Map two window coordinates w/ shared X/Y and distinctive Z
1804 * to a {@link Ray}. The resulting {@link Ray} maybe used for <i>picking</i>
1805 * using a {@link AABBox#getRayIntersection(vec3f, Ray, value_type, boolean)}.
1806 * <p>
1807 * Notes for picking <i>winz0</i> and <i>winz1</i>:
1808 * <ul>
1809 * <li>see {@link FloatUtil#getZBufferEpsilon(int, value_type, value_type)}</li>
1810 * <li>see {@link FloatUtil#getZBufferValue(int, value_type, value_type, value_type)}</li>
1811 * <li>see {@link FloatUtil#getOrthoWinZ(value_type, value_type, value_type)}</li>
1812 * </ul>
1813 * </p>
1814 * @param winx
1815 * @param winy
1816 * @param winz0
1817 * @param winz1
1818 * @param invPMv inverse [projection] x [modelview] matrix, i.e. Inv(P x Mv), if null method returns false
1819 * @param viewport Rect4i viewport
1820 * @param ray storage for the resulting {@link Ray}
1821 * @return true if successful, otherwise false (null invert matrix, or becomes z is infinity)
1822 */
1823 static bool mapWinToRay(const value_type winx, const value_type winy, const value_type winz0, const value_type winz1,
1824 const Matrix4& invPMv,
1825 const Recti& viewport,
1826 Ray3& ray) noexcept
1827 {
1828 if( mapWinToObj(winx, winy, winz0, winz1, invPMv, viewport, ray.orig, ray.dir) ) {
1829 (ray.dir -= ray.orig).normalize();
1830 return true;
1831 } else {
1832 return false;
1833 }
1834 }
1835
1836 /**
1837 * Returns a formatted string representation of this matrix
1838 * @param rowPrefix prefix for each row
1839 * @param f format string for each value_type element, e.g. "%10.5f"
1840 */
1841 std::string toString(const std::string& rowPrefix, const std::string& f) const noexcept {
1842 std::string sb;
1843 value_type tmp[16];
1844 get(tmp);
1845 return jau::mat_to_string(sb, rowPrefix, f, tmp, 4, 4, false /* rowMajorOrder */); // creates a copy-out!
1846 }
1847
1848 /**
1849 * Returns a formatted string representation of this matrix
1850 * @param rowPrefix prefix for each row
1851 */
1852 std::string toString(const std::string& rowPrefix) const noexcept { return toString(rowPrefix, "%13.9f"); }
1853
1854 std::string toString() const noexcept { return toString("", "%13.9f"); }
1855};
1856
1857template<typename T,
1858 std::enable_if_t<std::is_floating_point_v<T>, bool> = true>
1859constexpr Matrix4<T> operator*(const Matrix4<T>& lhs, const Matrix4<T>& rhs ) noexcept {
1860 Matrix4<T> r(lhs); r.mul(rhs); return r;
1861}
1862
1863template<typename T,
1864 std::enable_if_t<std::is_floating_point_v<T>, bool> = true>
1865constexpr Matrix4<T> operator*(const Matrix4<T>& lhs, const T s ) noexcept {
1866 Matrix4<T> r(lhs); r *= s; return r;
1867}
1868
1869template<typename T,
1870 std::enable_if_t<std::is_floating_point_v<T>, bool> = true>
1871constexpr Matrix4<T> operator*(const T s, const Matrix4<T>& rhs) noexcept {
1872 Matrix4<T> r(rhs); r *= s; return r;
1873}
1874
1875template<typename T,
1876 std::enable_if_t<std::is_floating_point_v<T>, bool> = true>
1877std::ostream& operator<<(std::ostream& out, const Matrix4<T>& v) noexcept {
1878 return out << v.toString();
1879}
1880
1882
1883static_assert(alignof(float) == alignof(Mat4f));
1884static_assert(sizeof(float)*16 == sizeof(Mat4f));
1885
1886/**@}*/
1887
1888} // namespace jau::math
1889
1890#endif // JAU_MAT4f_HPP_
#define E_FILE_LINE
Horizontal and vertical field of view (FOV) halves, allowing a non-centered projection.
float top
Half vertical FOV from center to top, either in inTangents or radians.
float right
Half horizontal FOV from center to right, either in inTangents or radians.
float left
Half horizontal FOV from center to left, either in inTangents or radians.
FovHVHalves toTangents() const noexcept
Returns this instance in tangent values.
float bottom
Half vertical FOV from center to bottom, either in inTangents or radians.
Basic 4x4 value_type matrix implementation using fields for intensive use-cases (host operations).
Definition: mat4f.hpp:112
constexpr Matrix4 & translate(const value_type x, const value_type y, const value_type z, Matrix4 &tmp) noexcept
Translate this matrix, i.e.
Definition: mat4f.hpp:1415
static bool mapWinToObj(const value_type winx, const value_type winy, const value_type winz1, const value_type winz2, const Matrix4 &invPMv, const Recti &viewport, Vec3 &objPos1, Vec3 &objPos2) noexcept
Map two window coordinates to two object coordinates, distinguished by their z component.
Definition: mat4f.hpp:1632
constexpr Matrix4 & setToScale(const value_type x, const value_type y, const value_type z) noexcept
Set this matrix to scale.
Definition: mat4f.hpp:958
constexpr bool operator==(const Matrix4 &rhs) const noexcept
Definition: mat4f.hpp:216
constexpr value_type get(const jau::nsize_t i) const noexcept
Returns the ith component of the given column-major order matrix, 0 <= i < 16, w/o boundary check.
Definition: mat4f.hpp:306
constexpr iterator get(iterator dst) const noexcept
Get this matrix into the given value_type[16] array in column major order w/o boundary check.
Definition: mat4f.hpp:396
Vector4F< value_type, std::is_floating_point_v< Value_type > > Vec4
Definition: mat4f.hpp:123
constexpr_cxx26 Matrix4 & setToRotationAxis(const value_type ang_rad, value_type x, value_type y, value_type z) noexcept
Set this matrix to rotation from the given axis and angle in radians.
Definition: mat4f.hpp:1002
value_type * pointer
Definition: mat4f.hpp:115
constexpr Vec4 getRow(const jau::nsize_t row) const noexcept
Get the named column of the given column-major matrix to v_out w/o boundary check.
Definition: mat4f.hpp:370
std::string toString() const noexcept
Definition: mat4f.hpp:1854
constexpr void set(const jau::nsize_t i, const value_type v) noexcept
Sets the ith component of this column-major order matrix with value_type v, 0 <= i < 16 w/o boundary ...
Definition: mat4f.hpp:231
constexpr value_type operator[](size_t i) const noexcept
Returns read-only ith component of the given column-major order matrix, 0 <= i < 16 w/o boundary chec...
Definition: mat4f.hpp:300
static bool mapWinToObj4(const value_type winx, const value_type winy, const value_type winz, const value_type clipw, const Matrix4 &mMv, const Matrix4 &mP, const Recti &viewport, const value_type near, const value_type far, Vec4 &objPos, Matrix4 &mat4Tmp) noexcept
Map window coordinates to object coordinates.
Definition: mat4f.hpp:1689
constexpr_cxx26 Matrix4 & setToRotationAxis(const value_type ang_rad, const Vec3 &axis) noexcept
Set this matrix to rotation from the given axis and angle in radians.
Definition: mat4f.hpp:1053
bool invert() noexcept
Invert this matrix.
Definition: mat4f.hpp:518
Matrix4 & setToPerspective(const FovHVHalves &fovhv, const value_type zNear, const value_type zFar)
Set this matrix to perspective frustum projection.
Definition: mat4f.hpp:1270
Value_type value_type
Definition: mat4f.hpp:114
static bool mapWinToRay(const value_type winx, const value_type winy, const value_type winz0, const value_type winz1, const Matrix4 &mMv, const Matrix4 &mP, const Recti &viewport, Ray3 &ray, Matrix4 &mat4Tmp1) noexcept
Map two window coordinates w/ shared X/Y and distinctive Z to a Ray.
Definition: mat4f.hpp:1782
constexpr Vec3 & getColumn(const jau::nsize_t column, Vec3 &v_out) const noexcept
Get the named column of the given column-major matrix to v_out w/o boundary check.
Definition: mat4f.hpp:347
constexpr Vec3 & mulVec3(const Vec3 &v_in, Vec3 &v_out) const noexcept
Affine 3f-vector transformation by 4x4 matrix.
Definition: mat4f.hpp:854
constexpr const_iterator cbegin() const noexcept
Definition: mat4f.hpp:312
constexpr reference operator[](size_t i) noexcept
Returns writable reference to the ith component of this column-major order matrix,...
Definition: mat4f.hpp:225
Ray3F< value_type, std::is_floating_point_v< Value_type > > Ray3
Definition: mat4f.hpp:124
Matrix4 & transpose() noexcept
Transpose this matrix.
Definition: mat4f.hpp:452
constexpr bool setToPick(const value_type x, const value_type y, const value_type deltaX, const value_type deltaY, const Recti &viewport, Matrix4 &mat4Tmp) noexcept
Set this matrix to the pick matrix based on given parameters.
Definition: mat4f.hpp:1363
constexpr Vec4 & mulVec4(const Vec4 &v_in, Vec4 &v_out) const noexcept
Definition: mat4f.hpp:806
constexpr Vec4 getColumn(const jau::nsize_t column) const noexcept
Get the named column of the given column-major matrix to v_out w/o boundary check.
Definition: mat4f.hpp:333
static bool mapWinToRay(const value_type winx, const value_type winy, const value_type winz0, const value_type winz1, const Matrix4 &invPMv, const Recti &viewport, Ray3 &ray) noexcept
Map two window coordinates w/ shared X/Y and distinctive Z to a Ray.
Definition: mat4f.hpp:1823
value_type determinant() const noexcept
Returns the determinant of this matrix.
Definition: mat4f.hpp:438
constexpr Matrix4 & scale(const value_type x, const value_type y, const value_type z, Matrix4 &tmp) noexcept
Scale this matrix, i.e.
Definition: mat4f.hpp:1437
Matrix4 & setToPerspective(const value_type fovy_rad, const value_type aspect, const value_type zNear, const value_type zFar)
Set this matrix to perspective frustum projection.
Definition: mat4f.hpp:1251
constexpr bool equals(const Matrix4 &o, const value_type epsilon=std::numeric_limits< value_type >::epsilon()) const noexcept
Definition: mat4f.hpp:194
const value_type & const_reference
Definition: mat4f.hpp:118
constexpr Vec3 & mulVec3(Vec3 &v_inout) const noexcept
Affine 3f-vector transformation by 4x4 matrix.
Definition: mat4f.hpp:889
constexpr iterator begin() noexcept
Definition: mat4f.hpp:237
static bool mapObjToWin(const Vec3 &obj, const Matrix4 &mPMv, const Recti &viewport, Vec3 &winPos) noexcept
Map object coordinates to window coordinates.
Definition: mat4f.hpp:1508
std::string toString(const std::string &rowPrefix, const std::string &f) const noexcept
Returns a formatted string representation of this matrix.
Definition: mat4f.hpp:1841
static bool mapWinToObj(const value_type winx, const value_type winy, const value_type winz, const Matrix4 &invPMv, const Recti &viewport, Vec3 &objPos) noexcept
Map window coordinates to object coordinates.
Definition: mat4f.hpp:1592
Matrix4 & setToFrustum(const value_type left, const value_type right, const value_type bottom, const value_type top, const value_type zNear, const value_type zFar)
Set this matrix to frustum.
Definition: mat4f.hpp:1202
constexpr Vec3 operator*(const Vec3 &rhs) const noexcept
Returns new Vec3, with affine 3f-vector transformation by this 4x4 matrix: this * v_in.
Definition: mat4f.hpp:871
constexpr Matrix4 & mul(const Matrix4 &a, const Matrix4 &b) noexcept
Multiply matrix: [this] = [a] x [b].
Definition: mat4f.hpp:773
Matrix4 & transpose(const Matrix4 &src) noexcept
Transpose the given src matrix into this matrix.
Definition: mat4f.hpp:488
constexpr Matrix4 & scale(const value_type s, Matrix4 &tmp) noexcept
Scale this matrix, i.e.
Definition: mat4f.hpp:1447
constexpr Vec4 & getColumn(const jau::nsize_t column, Vec4 &v_out) const noexcept
Get the named column of the given column-major matrix to v_out w/o boundary check.
Definition: mat4f.hpp:320
static bool mapObjToWin(const Vec3 &obj, const Matrix4 &mMv, const Matrix4 &mP, const Recti &viewport, Vec3 &winPos) noexcept
Map object coordinates to window coordinates.
Definition: mat4f.hpp:1468
constexpr Vec4 operator*(const Vec4 &rhs) const noexcept
Returns new Vec4, with this * v_in.
Definition: mat4f.hpp:820
constexpr Matrix4 & setToLookAt(const Vec3 &eye, const Vec3 &center, const Vec3 &up, Matrix4 &tmp) noexcept
Set this matrix to the look-at matrix based on given parameters.
Definition: mat4f.hpp:1299
constexpr Vec4 & getRow(const jau::nsize_t row, Vec4 &v_out) const noexcept
Get the named row of the given column-major matrix to v_out w/ boundary check.
Definition: mat4f.hpp:359
constexpr Matrix4 & operator*=(const Matrix4 &rhs) noexcept
Multiply matrix: [this] = [this] x [b].
Definition: mat4f.hpp:762
constexpr Matrix4 & loadIdentity() noexcept
Set this matrix to identity.
Definition: mat4f.hpp:250
constexpr_cxx26 Matrix4 & rotate(const value_type ang_rad, const value_type x, const value_type y, const value_type z, Matrix4 &tmp) noexcept
Rotate this matrix about give axis and angle in radians, i.e.
Definition: mat4f.hpp:1391
constexpr Matrix4 & operator=(const Matrix4 &o) noexcept
Copy assignment using the the values of the given src matrix.
Definition: mat4f.hpp:192
constexpr Vec4 & mulVec4(Vec4 &v_inout) const noexcept
Definition: mat4f.hpp:833
std::string toString(const std::string &rowPrefix) const noexcept
Returns a formatted string representation of this matrix.
Definition: mat4f.hpp:1852
constexpr Matrix4 & setToOrtho(const value_type left, const value_type right, const value_type bottom, const value_type top, const value_type zNear, const value_type zFar) noexcept
Set this matrix to orthogonal projection.
Definition: mat4f.hpp:1154
constexpr_cxx26 Matrix4 & setToRotationEuler(const value_type bankX, const value_type headingY, const value_type attitudeZ) noexcept
Set this matrix to rotation from the given Euler rotation angles in radians.
Definition: mat4f.hpp:1080
constexpr_cxx26 Matrix4 & rotate(const value_type ang_rad, const Vec3 &axis, Matrix4 &tmp) noexcept
Rotate this matrix about give axis and angle in radians, i.e.
Definition: mat4f.hpp:1403
constexpr std::vector< value_type > & get(std::vector< value_type > &dst, size_t dst_off) const noexcept
Get this matrix into the given FloatBuffer in column major order.
Definition: mat4f.hpp:424
constexpr Matrix4 & setToScale(const Vec3 &s) noexcept
Set this matrix to scale.
Definition: mat4f.hpp:982
Vector3F< value_type, std::is_floating_point_v< Value_type > > Vec3
Definition: mat4f.hpp:122
static bool mapWinToObj(const value_type winx, const value_type winy, const value_type winz, const Matrix4 &mMv, const Matrix4 &mP, const Recti &viewport, Vec3 &objPos, Matrix4 &mat4Tmp) noexcept
Map window coordinates to object coordinates.
Definition: mat4f.hpp:1547
constexpr Matrix4 & operator*=(const value_type s) noexcept
Multiply matrix with scalar: [this] = [this] x [s].
Definition: mat4f.hpp:703
constexpr Matrix4(const_iterator m) noexcept
Creates a new matrix based on given value_type[4*4] column major order.
Definition: mat4f.hpp:165
bool invert(const Matrix4 &src) noexcept
Invert the src matrix values into this matrix.
Definition: mat4f.hpp:597
constexpr_cxx26 Matrix4 & setToRotationEuler(const Vec3 &angradXYZ) noexcept
Set this matrix to rotation from the given Euler rotation angles in radians.
Definition: mat4f.hpp:1133
constexpr Matrix4 & setToTranslation(const Vec3 &t) noexcept
Set this matrix to translation.
Definition: mat4f.hpp:940
constexpr Matrix4() noexcept
Creates a new identity matrix.
Definition: mat4f.hpp:154
value_type * iterator
Definition: mat4f.hpp:119
const value_type * const_iterator
Definition: mat4f.hpp:120
constexpr Matrix4 & mul(const Matrix4 &b) noexcept
Multiply matrix: [this] = [this] x [b].
Definition: mat4f.hpp:717
constexpr Matrix4 & translate(const Vec3 &t, Matrix4 &tmp) noexcept
Translate this matrix, i.e.
Definition: mat4f.hpp:1425
constexpr Matrix4 & load(const_iterator src) noexcept
Load the values of the given matrix src to this matrix w/o boundary check.
Definition: mat4f.hpp:264
static bool mapWinToObj4(const value_type winx, const value_type winy, const value_type winz, const value_type clipw, const Matrix4 &invPMv, const Recti &viewport, const value_type near, const value_type far, Vec4 &objPos) noexcept
Map window coordinates to object coordinates.
Definition: mat4f.hpp:1735
value_type & reference
Definition: mat4f.hpp:117
constexpr Matrix4(const Matrix4 &o) noexcept
Creates a new matrix copying the values of the given src matrix.
Definition: mat4f.hpp:185
constexpr Vec3 & getRow(const jau::nsize_t row, Vec3 &v_out) const noexcept
Get the named row of the given column-major matrix to v_out w/o boundary check.
Definition: mat4f.hpp:383
constexpr Matrix4 & load(const Matrix4 &src) noexcept
Load the values of the given matrix src to this matrix w/o boundary check.
Definition: mat4f.hpp:289
constexpr Matrix4 & setToTranslation(const value_type x, const value_type y, const value_type z) noexcept
Set this matrix to translation.
Definition: mat4f.hpp:916
const value_type * const_pointer
Definition: mat4f.hpp:116
constexpr Matrix4(std::initializer_list< value_type > m) noexcept
Creates a new matrix based on given value_type initializer list in column major order.
Definition: mat4f.hpp:176
Quaternion implementation supporting Gimbal-Lock free rotations.
Definition: quaternion.hpp:62
Rectangle with x, y, width and height integer components.
Definition: recti.hpp:43
3D vector using three value_type components.
Definition: vec3f.hpp:51
constexpr Vector3F & set(const Vec2f &o, const value_type z_) noexcept
TODO constexpr bool operator<=>(const vec3f_t& rhs ) const noexcept { return ... }.
Definition: vec3f.hpp:143
value_type x
Definition: vec3f.hpp:72
constexpr Vector3F & normalize() noexcept
Normalize this vector in place.
Definition: vec3f.hpp:237
constexpr Vector3F cross(const Vector3F &b) const noexcept
cross product this x b
Definition: vec3f.hpp:285
value_type y
Definition: vec3f.hpp:73
value_type z
Definition: vec3f.hpp:74
4D vector using four value_type components.
Definition: vec4f.hpp:51
value_type w
Definition: vec4f.hpp:75
constexpr Vector4F & mul(const value_type sx, const value_type sy, const value_type sz, const value_type sw) noexcept
this = this * {sx, sy, sz, sw}, returns this.
Definition: vec4f.hpp:163
constexpr Vector4F & set(const Vec3f &o, const value_type w_) noexcept
TODO constexpr bool operator<=>(const vec4f_t& rhs ) const noexcept { return ... }.
Definition: vec4f.hpp:148
constexpr Vec3 & getVec3(Vec3 &out) const noexcept
out = { this.x, this.y, this.z } dropping w, returns out.
Definition: vec4f.hpp:128
constexpr Vector4F & add(const value_type dx, const value_type dy, const value_type dz, const value_type dw) noexcept
this = this + {dx, dy, dz, dw}, returns this.
Definition: vec4f.hpp:159
value_type x
Definition: vec4f.hpp:72
value_type y
Definition: vec4f.hpp:73
value_type z
Definition: vec4f.hpp:74
constexpr Vector4F & scale(const value_type s) noexcept
this = this * s, returns this.
Definition: vec4f.hpp:167
Providing frustum planes derived by different inputs (P*MV, ..) used to classify objects.
Definition: frustum.hpp:94
std::string to_string(const alphabet &v) noexcept
Definition: base_codec.hpp:97
#define constexpr_cxx26
Entry * get(const EUI48 &addr, const std::string &name, AddressNameEntryMatchFunc m) noexcept
Returns a matching BTSecurityRegistry::Entry with the given addr and/or name.
std::string & mat_to_string(std::string &sb, const std::string &rowPrefix, const std::string &f, const T a[], const jau::nsize_t rows, const jau::nsize_t columns, const bool rowMajorOrder) noexcept
Appends a matrix of floating points to the given string sb
Definition: float_math.hpp:537
std::enable_if< std::is_floating_point_v< T >, bool >::type constexpr equals(const T &a, const T &b, const T &epsilon=std::numeric_limits< T >::epsilon()) noexcept
Returns true if both values are equal, i.e.
Definition: float_math.hpp:394
uint_fast32_t nsize_t
Natural 'size_t' alternative using uint_fast32_t as its natural sized type.
Definition: int_types.hpp:53
constexpr T max(const T x, const T y) noexcept
Returns the maximum of two integrals (w/ branching) in O(1)
Definition: base_math.hpp:191
constexpr T abs(const T x) noexcept
Returns the absolute value of an arithmetic number (w/ branching) in O(1)
Definition: base_math.hpp:155
Matrix4< float > Mat4f
Definition: mat4f.hpp:1881
constexpr Matrix4< T > operator*(const Matrix4< T > &lhs, const Matrix4< T > &rhs) noexcept
Definition: mat4f.hpp:1859
std::ostream & operator<<(std::ostream &out, const Matrix4< T > &v) noexcept
Definition: mat4f.hpp:1877
constexpr const jau::fraction_i64 zero(0l, 1lu)
zero is 0/1
constexpr const jau::fraction_i64 one(1l, 1lu)
one is 10^0 or 1/1
Simple compound denoting a ray.
Definition: vec3f.hpp:399
uint8_t Value_type