| 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 |  | 
|---|