| 1 | // Copyright 2009-2021 Intel Corporation |
| 2 | // SPDX-License-Identifier: Apache-2.0 |
| 3 | |
| 4 | #pragma once |
| 5 | |
| 6 | #include "vec2.h" |
| 7 | |
| 8 | namespace embree |
| 9 | { |
| 10 | //////////////////////////////////////////////////////////////////////////////// |
| 11 | /// 2D Linear Transform (2x2 Matrix) |
| 12 | //////////////////////////////////////////////////////////////////////////////// |
| 13 | |
| 14 | template<typename T> struct LinearSpace2 |
| 15 | { |
| 16 | typedef T Vector; |
| 17 | typedef typename T::Scalar Scalar; |
| 18 | |
| 19 | /*! default matrix constructor */ |
| 20 | __forceinline LinearSpace2 ( ) {} |
| 21 | __forceinline LinearSpace2 ( const LinearSpace2& other ) { vx = other.vx; vy = other.vy; } |
| 22 | __forceinline LinearSpace2& operator=( const LinearSpace2& other ) { vx = other.vx; vy = other.vy; return *this; } |
| 23 | |
| 24 | template<typename L1> __forceinline LinearSpace2( const LinearSpace2<L1>& s ) : vx(s.vx), vy(s.vy) {} |
| 25 | |
| 26 | /*! matrix construction from column vectors */ |
| 27 | __forceinline LinearSpace2(const Vector& vx, const Vector& vy) |
| 28 | : vx(vx), vy(vy) {} |
| 29 | |
| 30 | /*! matrix construction from row mayor data */ |
| 31 | __forceinline LinearSpace2(const Scalar& m00, const Scalar& m01, |
| 32 | const Scalar& m10, const Scalar& m11) |
| 33 | : vx(m00,m10), vy(m01,m11) {} |
| 34 | |
| 35 | /*! compute the determinant of the matrix */ |
| 36 | __forceinline const Scalar det() const { return vx.x*vy.y - vx.y*vy.x; } |
| 37 | |
| 38 | /*! compute adjoint matrix */ |
| 39 | __forceinline const LinearSpace2 adjoint() const { return LinearSpace2(vy.y,-vy.x,-vx.y,vx.x); } |
| 40 | |
| 41 | /*! compute inverse matrix */ |
| 42 | __forceinline const LinearSpace2 inverse() const { return adjoint()/det(); } |
| 43 | |
| 44 | /*! compute transposed matrix */ |
| 45 | __forceinline const LinearSpace2 transposed() const { return LinearSpace2(vx.x,vx.y,vy.x,vy.y); } |
| 46 | |
| 47 | /*! returns first row of matrix */ |
| 48 | __forceinline Vector row0() const { return Vector(vx.x,vy.x); } |
| 49 | |
| 50 | /*! returns second row of matrix */ |
| 51 | __forceinline Vector row1() const { return Vector(vx.y,vy.y); } |
| 52 | |
| 53 | //////////////////////////////////////////////////////////////////////////////// |
| 54 | /// Constants |
| 55 | //////////////////////////////////////////////////////////////////////////////// |
| 56 | |
| 57 | __forceinline LinearSpace2( ZeroTy ) : vx(zero), vy(zero) {} |
| 58 | __forceinline LinearSpace2( OneTy ) : vx(one, zero), vy(zero, one) {} |
| 59 | |
| 60 | /*! return matrix for scaling */ |
| 61 | static __forceinline LinearSpace2 scale(const Vector& s) { |
| 62 | return LinearSpace2(s.x, 0, |
| 63 | 0 , s.y); |
| 64 | } |
| 65 | |
| 66 | /*! return matrix for rotation */ |
| 67 | static __forceinline LinearSpace2 rotate(const Scalar& r) { |
| 68 | Scalar s = sin(r), c = cos(r); |
| 69 | return LinearSpace2(c, -s, |
| 70 | s, c); |
| 71 | } |
| 72 | |
| 73 | /*! return closest orthogonal matrix (i.e. a general rotation including reflection) */ |
| 74 | LinearSpace2 orthogonal() const |
| 75 | { |
| 76 | LinearSpace2 m = *this; |
| 77 | |
| 78 | // mirrored? |
| 79 | Scalar mirror(one); |
| 80 | if (m.det() < Scalar(zero)) { |
| 81 | m.vx = -m.vx; |
| 82 | mirror = -mirror; |
| 83 | } |
| 84 | |
| 85 | // rotation |
| 86 | for (int i = 0; i < 99; i++) { |
| 87 | const LinearSpace2 m_next = 0.5 * (m + m.transposed().inverse()); |
| 88 | const LinearSpace2 d = m_next - m; |
| 89 | m = m_next; |
| 90 | // norm^2 of difference small enough? |
| 91 | if (max(dot(d.vx, d.vx), dot(d.vy, d.vy)) < 1e-8) |
| 92 | break; |
| 93 | } |
| 94 | |
| 95 | // rotation * mirror_x |
| 96 | return LinearSpace2(mirror*m.vx, m.vy); |
| 97 | } |
| 98 | |
| 99 | public: |
| 100 | |
| 101 | /*! the column vectors of the matrix */ |
| 102 | Vector vx,vy; |
| 103 | }; |
| 104 | |
| 105 | //////////////////////////////////////////////////////////////////////////////// |
| 106 | // Unary Operators |
| 107 | //////////////////////////////////////////////////////////////////////////////// |
| 108 | |
| 109 | template<typename T> __forceinline LinearSpace2<T> operator -( const LinearSpace2<T>& a ) { return LinearSpace2<T>(-a.vx,-a.vy); } |
| 110 | template<typename T> __forceinline LinearSpace2<T> operator +( const LinearSpace2<T>& a ) { return LinearSpace2<T>(+a.vx,+a.vy); } |
| 111 | template<typename T> __forceinline LinearSpace2<T> rcp ( const LinearSpace2<T>& a ) { return a.inverse(); } |
| 112 | |
| 113 | //////////////////////////////////////////////////////////////////////////////// |
| 114 | // Binary Operators |
| 115 | //////////////////////////////////////////////////////////////////////////////// |
| 116 | |
| 117 | template<typename T> __forceinline LinearSpace2<T> operator +( const LinearSpace2<T>& a, const LinearSpace2<T>& b ) { return LinearSpace2<T>(a.vx+b.vx,a.vy+b.vy); } |
| 118 | template<typename T> __forceinline LinearSpace2<T> operator -( const LinearSpace2<T>& a, const LinearSpace2<T>& b ) { return LinearSpace2<T>(a.vx-b.vx,a.vy-b.vy); } |
| 119 | |
| 120 | template<typename T> __forceinline LinearSpace2<T> operator*(const typename T::Scalar & a, const LinearSpace2<T>& b) { return LinearSpace2<T>(a*b.vx, a*b.vy); } |
| 121 | template<typename T> __forceinline T operator*(const LinearSpace2<T>& a, const T & b) { return b.x*a.vx + b.y*a.vy; } |
| 122 | template<typename T> __forceinline LinearSpace2<T> operator*(const LinearSpace2<T>& a, const LinearSpace2<T>& b) { return LinearSpace2<T>(a*b.vx, a*b.vy); } |
| 123 | |
| 124 | template<typename T> __forceinline LinearSpace2<T> operator/(const LinearSpace2<T>& a, const typename T::Scalar & b) { return LinearSpace2<T>(a.vx/b, a.vy/b); } |
| 125 | template<typename T> __forceinline LinearSpace2<T> operator/(const LinearSpace2<T>& a, const LinearSpace2<T>& b) { return a * rcp(b); } |
| 126 | |
| 127 | template<typename T> __forceinline LinearSpace2<T>& operator *=( LinearSpace2<T>& a, const LinearSpace2<T>& b ) { return a = a * b; } |
| 128 | template<typename T> __forceinline LinearSpace2<T>& operator /=( LinearSpace2<T>& a, const LinearSpace2<T>& b ) { return a = a / b; } |
| 129 | |
| 130 | //////////////////////////////////////////////////////////////////////////////// |
| 131 | /// Comparison Operators |
| 132 | //////////////////////////////////////////////////////////////////////////////// |
| 133 | |
| 134 | template<typename T> __forceinline bool operator ==( const LinearSpace2<T>& a, const LinearSpace2<T>& b ) { return a.vx == b.vx && a.vy == b.vy; } |
| 135 | template<typename T> __forceinline bool operator !=( const LinearSpace2<T>& a, const LinearSpace2<T>& b ) { return a.vx != b.vx || a.vy != b.vy; } |
| 136 | |
| 137 | //////////////////////////////////////////////////////////////////////////////// |
| 138 | /// Output Operators |
| 139 | //////////////////////////////////////////////////////////////////////////////// |
| 140 | |
| 141 | template<typename T> static embree_ostream operator<<(embree_ostream cout, const LinearSpace2<T>& m) { |
| 142 | return cout << "{ vx = " << m.vx << ", vy = " << m.vy << "}" ; |
| 143 | } |
| 144 | |
| 145 | /*! Shortcuts for common linear spaces. */ |
| 146 | typedef LinearSpace2<Vec2f> LinearSpace2f; |
| 147 | typedef LinearSpace2<Vec2fa> LinearSpace2fa; |
| 148 | } |
| 149 | |