1// Copyright 2009-2021 Intel Corporation 2// SPDX-License-Identifier: Apache-2.0 3 4#pragma once 5 6#include "vec.isph" 7#include "quaternion.isph" 8 9struct LinearSpace3f 10{ 11 Vec3f vx; 12 Vec3f vy; 13 Vec3f vz; 14}; 15 16struct LinearSpace3fa 17{ 18 Vec3fa vx; 19 Vec3fa vy; 20 Vec3fa vz; 21}; 22 23//////////////////////////////////////////////////////////////////////////////// 24/// Constructors 25//////////////////////////////////////////////////////////////////////////////// 26 27inline uniform LinearSpace3f make_LinearSpace3f(const uniform Vec3f &x, const uniform Vec3f &y, const uniform Vec3f &z) { 28 uniform LinearSpace3f l; l.vx = x; l.vy = y; l.vz = z; return l; 29} 30 31inline varying LinearSpace3f make_LinearSpace3f(const varying Vec3f &x, const varying Vec3f &y, const varying Vec3f &z) { 32 varying LinearSpace3f l; l.vx = x; l.vy = y; l.vz = z; return l; 33} 34 35inline uniform LinearSpace3f make_LinearSpace3f(const uniform LinearSpace3fa &s) { 36 uniform LinearSpace3f l; l.vx = make_Vec3f(s.vx); l.vy = make_Vec3f(s.vy); l.vz = make_Vec3f(s.vz); return l; 37} 38 39inline varying LinearSpace3f make_LinearSpace3f(const varying LinearSpace3fa &s) { 40 varying LinearSpace3f l; l.vx = make_Vec3f(s.vx); l.vy = make_Vec3f(s.vy); l.vz = make_Vec3f(s.vz); return l; 41} 42 43/*! construction from quaternion */ 44inline uniform LinearSpace3f make_LinearSpace3f(const uniform Quaternion3f& q) 45{ 46 uniform LinearSpace3f l; 47 l.vx = make_Vec3f(q.r*q.r + q.i*q.i - q.j*q.j - q.k*q.k, 2.0f*(q.i*q.j + q.r*q.k), 2.0f*(q.i*q.k - q.r*q.j)); 48 l.vy = make_Vec3f(2.0f*(q.i*q.j - q.r*q.k), (q.r*q.r - q.i*q.i + q.j*q.j - q.k*q.k), 2.0f*(q.j*q.k + q.r*q.i)); 49 l.vz = make_Vec3f(2.0f*(q.i*q.k + q.r*q.j), 2.0f*(q.j*q.k - q.r*q.i), q.r*q.r - q.i*q.i - q.j*q.j + q.k*q.k); 50 return l; 51} 52 53inline varying LinearSpace3f make_LinearSpace3f(const varying Quaternion3f& q) 54{ 55 LinearSpace3f l; 56 l.vx = make_Vec3f(q.r*q.r + q.i*q.i - q.j*q.j - q.k*q.k, 2.0f*(q.i*q.j + q.r*q.k), 2.0f*(q.i*q.k - q.r*q.j)); 57 l.vy = make_Vec3f(2.0f*(q.i*q.j - q.r*q.k), (q.r*q.r - q.i*q.i + q.j*q.j - q.k*q.k), 2.0f*(q.j*q.k + q.r*q.i)); 58 l.vz = make_Vec3f(2.0f*(q.i*q.k + q.r*q.j), 2.0f*(q.j*q.k - q.r*q.i), q.r*q.r - q.i*q.i - q.j*q.j + q.k*q.k); 59 return l; 60} 61 62inline uniform LinearSpace3f make_LinearSpace3f_identity() 63{ 64 return make_LinearSpace3f(make_Vec3f(1.f,0.f,0.f), 65 make_Vec3f(0.f,1.f,0.f), 66 make_Vec3f(0.f,0.f,1.f)); 67} 68 69/*! return scale matrix */ 70inline uniform LinearSpace3f make_LinearSpace3f_scale(const uniform Vec3f& s) 71{ 72 return make_LinearSpace3f( make_Vec3f( s.x, 0.0f, 0.0f), 73 make_Vec3f(0.0f, s.y, 0.0f), 74 make_Vec3f(0.0f, 0.0f, s.z)); 75} 76 77/*! return matrix for rotation around arbitrary axis */ 78inline uniform LinearSpace3f make_LinearSpace3f_rotate(const uniform Vec3f& _u, const uniform float& r) 79{ 80 uniform Vec3f u = normalize(_u); 81 uniform float s = sin(r), c = cos(r); 82 return make_LinearSpace3f(make_Vec3f(u.x*u.x+(1-u.x*u.x)*c, u.x*u.y*(1-c)+u.z*s, u.x*u.z*(1-c)-u.y*s), 83 make_Vec3f(u.x*u.y*(1-c)-u.z*s, u.y*u.y+(1-u.y*u.y)*c, u.y*u.z*(1-c)+u.x*s), 84 make_Vec3f(u.x*u.z*(1-c)+u.y*s, u.y*u.z*(1-c)-u.x*s, u.z*u.z+(1-u.z*u.z)*c)); 85} 86 87/*! return matrix for rotation around arbitrary axis */ 88inline varying LinearSpace3f make_LinearSpace3f_rotate(const varying Vec3f& _u, const varying float& r) 89{ 90 varying Vec3f u = normalize(_u); 91 varying float s = sin(r), c = cos(r); 92 return make_LinearSpace3f(make_Vec3f(u.x*u.x+(1-u.x*u.x)*c, u.x*u.y*(1-c)+u.z*s, u.x*u.z*(1-c)-u.y*s), 93 make_Vec3f(u.x*u.y*(1-c)-u.z*s, u.y*u.y+(1-u.y*u.y)*c, u.y*u.z*(1-c)+u.x*s), 94 make_Vec3f(u.x*u.z*(1-c)+u.y*s, u.y*u.z*(1-c)-u.x*s, u.z*u.z+(1-u.z*u.z)*c)); 95} 96 97/*! return matrix for rotation around arbitrary axis */ 98inline varying LinearSpace3f make_LinearSpace3f_rotate(const uniform Vec3f& _u, const varying float& r) 99{ 100 varying Vec3f u = normalize(_u); 101 varying float s = sin(r), c = cos(r); 102 return make_LinearSpace3f(make_Vec3f(u.x*u.x+(1-u.x*u.x)*c, u.x*u.y*(1-c)+u.z*s, u.x*u.z*(1-c)-u.y*s), 103 make_Vec3f(u.x*u.y*(1-c)-u.z*s, u.y*u.y+(1-u.y*u.y)*c, u.y*u.z*(1-c)+u.x*s), 104 make_Vec3f(u.x*u.z*(1-c)+u.y*s, u.y*u.z*(1-c)-u.x*s, u.z*u.z+(1-u.z*u.z)*c)); 105} 106 107//////////////////////////////////////////////////////////////////////////////// 108// Unary Operators 109//////////////////////////////////////////////////////////////////////////////// 110 111inline uniform LinearSpace3f neg(const uniform LinearSpace3f &l) { return make_LinearSpace3f(neg(l.vx),neg(l.vy),neg(l.vz)); } 112inline varying LinearSpace3f neg(const varying LinearSpace3f &l) { return make_LinearSpace3f(neg(l.vx),neg(l.vy),neg(l.vz)); } 113 114/*! compute the determinant of the matrix */ 115inline uniform float det(const uniform LinearSpace3f &l) { return dot(l.vx,cross(l.vy,l.vz)); } 116inline varying float det(const varying LinearSpace3f &l) { return dot(l.vx,cross(l.vy,l.vz)); } 117 118/*! compute transposed matrix */ 119inline uniform LinearSpace3f transposed(const uniform LinearSpace3f &l) { 120 return make_LinearSpace3f(make_Vec3f(l.vx.x,l.vy.x,l.vz.x), 121 make_Vec3f(l.vx.y,l.vy.y,l.vz.y), 122 make_Vec3f(l.vx.z,l.vy.z,l.vz.z)); 123} 124inline varying LinearSpace3f transposed(const varying LinearSpace3f &l) { 125 return make_LinearSpace3f(make_Vec3f(l.vx.x,l.vy.x,l.vz.x), 126 make_Vec3f(l.vx.y,l.vy.y,l.vz.y), 127 make_Vec3f(l.vx.z,l.vy.z,l.vz.z)); 128} 129 130/*! compute adjoint matrix */ 131inline uniform LinearSpace3f adjoint(const uniform LinearSpace3f &l) { 132 return transposed(make_LinearSpace3f(cross(l.vy,l.vz),cross(l.vz,l.vx),cross(l.vx,l.vy))); 133} 134inline varying LinearSpace3f adjoint(const varying LinearSpace3f &l) { 135 return transposed(make_LinearSpace3f(cross(l.vy,l.vz),cross(l.vz,l.vx),cross(l.vx,l.vy))); 136} 137 138/*! calculates orthogonal coordinate frame with z-Vector pointing towards N */ 139inline uniform LinearSpace3f frame(const uniform Vec3f &N) 140{ 141 const uniform Vec3f dx0 = make_Vec3f(0.0f,N.z,-N.y); 142 const uniform Vec3f dx1 = make_Vec3f(-N.z,0.0f,N.x); 143 const uniform Vec3f dx = normalize(dot(dx0,dx0) > dot(dx1,dx1) ? dx0 : dx1); 144 const uniform Vec3f dy = normalize(cross(N,dx)); 145 return make_LinearSpace3f(dx,dy,N); 146} 147 148inline varying LinearSpace3f frame(const varying Vec3f &N) 149{ 150 const varying Vec3f dx0 = make_Vec3f(0.0f,N.z,-N.y); 151 const varying Vec3f dx1 = make_Vec3f(-N.z,0.0f,N.x); 152 const varying Vec3f dx = normalize(dot(dx0,dx0) > dot(dx1,dx1) ? dx0 : dx1); 153 const varying Vec3f dy = normalize(cross(N,dx)); 154 return make_LinearSpace3f(dx,dy,N); 155} 156 157//////////////////////////////////////////////////////////////////////////////// 158/// Binary Operators 159//////////////////////////////////////////////////////////////////////////////// 160 161inline uniform LinearSpace3f operator+(const uniform LinearSpace3f &a, const uniform LinearSpace3f &b) { return make_LinearSpace3f(a.vx+b.vx, a.vy+b.vy, a.vz+b.vz); } 162inline varying LinearSpace3f operator+(const varying LinearSpace3f &a, const varying LinearSpace3f &b) { return make_LinearSpace3f(a.vx+b.vx, a.vy+b.vy, a.vz+b.vz); } 163 164inline uniform LinearSpace3f operator-(const uniform LinearSpace3f &a, const uniform LinearSpace3f &b) { return make_LinearSpace3f(a.vx-b.vx, a.vy-b.vy, a.vz-b.vz); } 165inline varying LinearSpace3f operator-(const varying LinearSpace3f &a, const varying LinearSpace3f &b) { return make_LinearSpace3f(a.vx-b.vx, a.vy-b.vy, a.vz-b.vz); } 166 167inline uniform Vec3f operator*(const uniform LinearSpace3f &l, const uniform Vec3f &v) { return v.x*l.vx + v.y*l.vy + v.z*l.vz; } 168inline varying Vec3f operator*(const uniform LinearSpace3f &l, const varying Vec3f &v) { return v.x*l.vx + v.y*l.vy + v.z*l.vz; } 169inline varying Vec3f operator*(const varying LinearSpace3f &l, const varying Vec3f &v) { return v.x*l.vx + v.y*l.vy + v.z*l.vz; } 170 171inline uniform LinearSpace3f operator*(const uniform float &a, const uniform LinearSpace3f &b) { return make_LinearSpace3f(a*b.vx, a*b.vy, a*b.vz); } 172inline uniform LinearSpace3f operator*(const uniform LinearSpace3f &a, const uniform float &b) { return make_LinearSpace3f(a.vx*b, a.vy*b, a.vz*b); } 173inline uniform LinearSpace3f operator*(const uniform LinearSpace3f &a, const uniform LinearSpace3f &b) { return make_LinearSpace3f(a*b.vx, a*b.vy, a*b.vz); } 174 175inline varying LinearSpace3f operator*(const varying float &a, const varying LinearSpace3f &b) { return make_LinearSpace3f(a*b.vx, a*b.vy, a*b.vz); } 176inline varying LinearSpace3f operator*(const varying LinearSpace3f &a, const varying float &b) { return make_LinearSpace3f(a.vx*b, a.vy*b, a.vz*b); } 177inline varying LinearSpace3f operator*(const varying LinearSpace3f &a, const varying LinearSpace3f &b) { return make_LinearSpace3f(a*b.vx, a*b.vy, a*b.vz); } 178 179inline varying LinearSpace3f operator*(const varying float &a, const uniform LinearSpace3f &b) { return make_LinearSpace3f(a*b.vx, a*b.vy, a*b.vz); } 180inline varying LinearSpace3f operator*(const uniform LinearSpace3f &a, const varying float &b) { return make_LinearSpace3f(a.vx*b, a.vy*b, a.vz*b); } 181 182inline uniform Vec3f xfmVector(const uniform LinearSpace3f l, const uniform Vec3f v) { return v.x*l.vx + v.y*l.vy + v.z*l.vz; } 183inline varying Vec3f xfmVector(const uniform LinearSpace3f l, const varying Vec3f v) { return v.x*l.vx + v.y*l.vy + v.z*l.vz; } 184inline varying Vec3f xfmVector(const varying LinearSpace3f l, const varying Vec3f v) { return v.x*l.vx + v.y*l.vy + v.z*l.vz; } 185 186//////////////////////////////////////////////////////////////////////////////// 187/// Comparison Operators 188//////////////////////////////////////////////////////////////////////////////// 189 190inline uniform bool eq(const uniform LinearSpace3f &a, const uniform LinearSpace3f &b) { return eq(a.vx,b.vx) & eq(a.vy,b.vy) & eq(a.vz,b.vz); } 191inline varying bool eq(const varying LinearSpace3f &a, const varying LinearSpace3f &b) { return eq(a.vx,b.vx) & eq(a.vy,b.vy) & eq(a.vz,b.vz); } 192 193inline uniform bool ne(const uniform LinearSpace3f &a, const uniform LinearSpace3f &b) { return ne(a.vx,b.vx) | ne(a.vy,b.vy) | ne(a.vz,b.vz); } 194inline varying bool ne(const varying LinearSpace3f &a, const varying LinearSpace3f &b) { return ne(a.vx,b.vx) | ne(a.vy,b.vy) | ne(a.vz,b.vz); } 195 196//////////////////////////////////////////////////////////////////////////////// 197// Unary Operators 198//////////////////////////////////////////////////////////////////////////////// 199 200/*! compute inverse matrix */ 201inline uniform LinearSpace3f rcp(const uniform LinearSpace3f &l) { return rcp(det(l))*adjoint(l); } 202inline varying LinearSpace3f rcp(const varying LinearSpace3f &l) { return rcp(det(l))*adjoint(l); } 203 204 205//////////////////////////////////////////////////////////////////////////////// 206// Interpolation 207//////////////////////////////////////////////////////////////////////////////// 208 209inline uniform LinearSpace3f lerp(uniform float factor, const uniform LinearSpace3f& a, const uniform LinearSpace3f& b) { 210 return make_LinearSpace3f(lerp(factor, a.vx, b.vx), lerp(factor, a.vy, b.vy), lerp(factor, a.vz, b.vz)); 211} 212 213inline varying LinearSpace3f lerp(varying float factor, const uniform LinearSpace3f& a, const uniform LinearSpace3f& b) { 214 return make_LinearSpace3f(lerp(factor, a.vx, b.vx), lerp(factor, a.vy, b.vy), lerp(factor, a.vz, b.vz)); 215} 216 217inline varying LinearSpace3f lerp(varying float factor, const varying LinearSpace3f& a, const varying LinearSpace3f& b) { 218 return make_LinearSpace3f(lerp(factor, a.vx, b.vx), lerp(factor, a.vy, b.vy), lerp(factor, a.vz, b.vz)); 219} 220 221//////////////////////////////////////////////////////////////////////////////// 222// Norms 223//////////////////////////////////////////////////////////////////////////////// 224 225inline uniform float norm_1(const uniform LinearSpace3f& l) 226{ 227 uniform float sum = 0; 228 sum += reduce_add(abs(l.vx)); 229 sum += reduce_add(abs(l.vy)); 230 sum += reduce_add(abs(l.vz)); 231 return sum / 9.f; 232} 233 234inline uniform Quaternion3f make_Quaternion3f(const uniform LinearSpace3f & l) 235{ 236 uniform Quaternion3f q; 237 const uniform float tr = l.vx.x + l.vy.y + l.vz.z; 238 if (tr > 0) 239 { 240 const uniform float S = sqrt(tr+1.0) * 2; 241 q.r = 0.25 * S; 242 q.i = (l.vy.z - l.vz.y) / S; 243 q.j = (l.vz.x - l.vx.z) / S; 244 q.k = (l.vx.y - l.vy.x) / S; 245 } 246 else if ((l.vx.x > l.vy.y)&(l.vx.x > l.vz.z)) 247 { 248 const uniform float S = sqrt(1.0 + l.vx.x - l.vy.y - l.vz.z) * 2; 249 q.r = (l.vy.z - l.vz.y) / S; 250 q.i = 0.25 * S; 251 q.j = (l.vy.x + l.vx.y) / S; 252 q.k = (l.vz.x + l.vx.z) / S; 253 } 254 else if (l.vy.y > l.vz.z) 255 { 256 const uniform float S = sqrt(1.0 + l.vy.y - l.vx.x - l.vz.z) * 2; 257 q.r = (l.vz.x - l.vx.z) / S; 258 q.i = (l.vy.x + l.vx.y) / S; 259 q.j = 0.25 * S; 260 q.k = (l.vz.y + l.vy.z) / S; 261 } 262 else 263 { 264 const uniform float S = sqrt(1.0 + l.vz.z - l.vx.x - l.vy.y) * 2; 265 q.r = (l.vx.y - l.vy.x) / S; 266 q.i = (l.vz.x + l.vx.z) / S; 267 q.j = (l.vz.y + l.vy.z) / S; 268 q.k = 0.25 * S; 269 } 270 return q; 271} 272 273inline varying Quaternion3f make_Quaternion3f(const varying LinearSpace3f & l) 274{ 275 Quaternion3f q; 276 const float tr = l.vx.x + l.vy.y + l.vz.z; 277 if (tr > 0) 278 { 279 const float S = sqrt(tr+1.0) * 2; 280 q.r = 0.25 * S; 281 q.i = (l.vy.z - l.vz.y) / S; 282 q.j = (l.vz.x - l.vx.z) / S; 283 q.k = (l.vx.y - l.vy.x) / S; 284 } 285 else if ((l.vx.x > l.vy.y)&(l.vx.x > l.vz.z)) 286 { 287 const float S = sqrt(1.0 + l.vx.x - l.vy.y - l.vz.z) * 2; 288 q.r = (l.vy.z - l.vz.y) / S; 289 q.i = 0.25 * S; 290 q.j = (l.vy.x + l.vx.y) / S; 291 q.k = (l.vz.x + l.vx.z) / S; 292 } 293 else if (l.vy.y > l.vz.z) 294 { 295 const float S = sqrt(1.0 + l.vy.y - l.vx.x - l.vz.z) * 2; 296 q.r = (l.vz.x - l.vx.z) / S; 297 q.i = (l.vy.x + l.vx.y) / S; 298 q.j = 0.25 * S; 299 q.k = (l.vz.y + l.vy.z) / S; 300 } 301 else 302 { 303 const float S = sqrt(1.0 + l.vz.z - l.vx.x - l.vy.y) * 2; 304 q.r = (l.vx.y - l.vy.x) / S; 305 q.i = (l.vz.x + l.vx.z) / S; 306 q.j = (l.vz.y + l.vy.z) / S; 307 q.k = 0.25 * S; 308 } 309 return q; 310}