1 #ifndef _TCUMATRIX_HPP 2 #define _TCUMATRIX_HPP 3 /*------------------------------------------------------------------------- 4 * drawElements Quality Program Tester Core 5 * ---------------------------------------- 6 * 7 * Copyright 2014 The Android Open Source Project 8 * 9 * Licensed under the Apache License, Version 2.0 (the "License"); 10 * you may not use this file except in compliance with the License. 11 * You may obtain a copy of the License at 12 * 13 * http://www.apache.org/licenses/LICENSE-2.0 14 * 15 * Unless required by applicable law or agreed to in writing, software 16 * distributed under the License is distributed on an "AS IS" BASIS, 17 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 18 * See the License for the specific language governing permissions and 19 * limitations under the License. 20 * 21 *//*! 22 * \file 23 * \brief Templatized matrix class. 24 *//*--------------------------------------------------------------------*/ 25 26 #include "tcuDefs.hpp" 27 #include "tcuVector.hpp" 28 #include "tcuArray.hpp" 29 30 namespace tcu 31 { 32 33 // Templated matrix class. 34 template <typename T, int Rows, int Cols> 35 class Matrix 36 { 37 public: 38 typedef Vector<T, Rows> Element; 39 typedef T Scalar; 40 41 enum 42 { 43 SIZE = Cols, 44 ROWS = Rows, 45 COLS = Cols, 46 }; 47 48 Matrix (void); 49 explicit Matrix (const T& src); 50 explicit Matrix (const T src[Rows*Cols]); 51 Matrix (const Vector<T, Rows>& src); 52 Matrix (const Matrix<T, Rows, Cols>& src); 53 ~Matrix (void); 54 55 Matrix<T, Rows, Cols>& operator= (const Matrix<T, Rows, Cols>& src); 56 Matrix<T, Rows, Cols>& operator*= (const Matrix<T, Rows, Cols>& src); 57 58 void setRow (int rowNdx, const Vector<T, Cols>& vec); 59 void setColumn (int colNdx, const Vector<T, Rows>& vec); 60 61 Vector<T, Cols> getRow (int ndx) const; 62 Vector<T, Rows>& getColumn (int ndx); 63 const Vector<T, Rows>& getColumn (int ndx) const; 64 65 Vector<T, Rows>& operator[] (int ndx) { return getColumn(ndx); } 66 const Vector<T, Rows>& operator[] (int ndx) const { return getColumn(ndx); } 67 68 inline const T& operator() (int row, int col) const { return m_data[col][row]; } 69 inline T& operator() (int row, int col) { return m_data[col][row]; } 70 71 Array<T, Rows*Cols> getRowMajorData (void) const; 72 Array<T, Rows*Cols> getColumnMajorData (void) const; 73 74 private: 75 Vector<Vector<T, Rows>, Cols> m_data; 76 } DE_WARN_UNUSED_TYPE; 77 78 // Operators. 79 80 // Mat * Mat. 81 template <typename T, int Rows0, int Cols0, int Rows1, int Cols1> 82 Matrix<T, Rows0, Cols1> operator* (const Matrix<T, Rows0, Cols0>& a, const Matrix<T, Rows1, Cols1>& b); 83 84 // Mat * Vec (column vector). 85 template <typename T, int Rows, int Cols> 86 Vector<T, Rows> operator* (const Matrix<T, Rows, Cols>& mtx, const Vector<T, Cols>& vec); 87 88 // Vec * Mat (row vector). 89 template <typename T, int Rows, int Cols> 90 Vector<T, Cols> operator* (const Vector<T, Rows>& vec, const Matrix<T, Rows, Cols>& mtx); 91 92 template <typename T, int Rows, int Cols> 93 bool operator== (const Matrix<T, Rows, Cols>& lhs, const Matrix<T, Rows, Cols>& rhs); 94 95 template <typename T, int Rows, int Cols> 96 bool operator!= (const Matrix<T, Rows, Cols>& lhs, const Matrix<T, Rows, Cols>& rhs); 97 98 // Further operations 99 100 template <typename T, int Size> 101 struct SquareMatrixOps 102 { 103 static T doDeterminant (const Matrix<T, Size, Size>& mat); 104 static Matrix<T, Size, Size> doInverse (const Matrix<T, Size, Size>& mat); 105 }; 106 107 template <typename T> 108 struct SquareMatrixOps<T, 2> 109 { 110 static T doDeterminant (const Matrix<T, 2, 2>& mat); 111 static Matrix<T, 2, 2> doInverse (const Matrix<T, 2, 2>& mat); 112 }; 113 114 template <typename T> 115 struct SquareMatrixOps<T, 3> 116 { 117 static T doDeterminant (const Matrix<T, 3, 3>& mat); 118 static Matrix<T, 3, 3> doInverse (const Matrix<T, 3, 3>& mat); 119 }; 120 121 template <typename T> 122 struct SquareMatrixOps<T, 4> 123 { 124 static T doDeterminant (const Matrix<T, 4, 4>& mat); 125 static Matrix<T, 4, 4> doInverse (const Matrix<T, 4, 4>& mat); 126 }; 127 128 namespace matrix 129 { 130 131 template <typename T, int Size> 132 T determinant (const Matrix<T, Size, Size>& mat) 133 { 134 return SquareMatrixOps<T, Size>::doDeterminant(mat); 135 } 136 137 template <typename T, int Size> 138 Matrix<T, Size, Size> inverse (const Matrix<T, Size, Size>& mat) 139 { 140 return SquareMatrixOps<T, Size>::doInverse(mat); 141 } 142 143 } // matrix 144 145 // Template implementations. 146 147 template <typename T> 148 T SquareMatrixOps<T, 2>::doDeterminant (const Matrix<T, 2, 2>& mat) 149 { 150 return mat(0,0) * mat(1,1) - mat(1,0) * mat(0,1); 151 } 152 153 template <typename T> 154 T SquareMatrixOps<T, 3>::doDeterminant (const Matrix<T, 3, 3>& mat) 155 { 156 return + mat(0,0) * mat(1,1) * mat(2,2) 157 + mat(0,1) * mat(1,2) * mat(2,0) 158 + mat(0,2) * mat(1,0) * mat(2,1) 159 - mat(0,0) * mat(1,2) * mat(2,1) 160 - mat(0,1) * mat(1,0) * mat(2,2) 161 - mat(0,2) * mat(1,1) * mat(2,0); 162 } 163 164 template <typename T> 165 T SquareMatrixOps<T, 4>::doDeterminant (const Matrix<T, 4, 4>& mat) 166 { 167 using matrix::determinant; 168 169 const T minorMatrices[4][3*3] = 170 { 171 { 172 mat(1,1), mat(2,1), mat(3,1), 173 mat(1,2), mat(2,2), mat(3,2), 174 mat(1,3), mat(2,3), mat(3,3), 175 }, 176 { 177 mat(1,0), mat(2,0), mat(3,0), 178 mat(1,2), mat(2,2), mat(3,2), 179 mat(1,3), mat(2,3), mat(3,3), 180 }, 181 { 182 mat(1,0), mat(2,0), mat(3,0), 183 mat(1,1), mat(2,1), mat(3,1), 184 mat(1,3), mat(2,3), mat(3,3), 185 }, 186 { 187 mat(1,0), mat(2,0), mat(3,0), 188 mat(1,1), mat(2,1), mat(3,1), 189 mat(1,2), mat(2,2), mat(3,2), 190 } 191 }; 192 193 return + mat(0,0) * determinant(Matrix<T, 3, 3>(minorMatrices[0])) 194 - mat(0,1) * determinant(Matrix<T, 3, 3>(minorMatrices[1])) 195 + mat(0,2) * determinant(Matrix<T, 3, 3>(minorMatrices[2])) 196 - mat(0,3) * determinant(Matrix<T, 3, 3>(minorMatrices[3])); 197 } 198 199 template <typename T> 200 Matrix<T, 2, 2> SquareMatrixOps<T, 2>::doInverse (const Matrix<T, 2, 2>& mat) 201 { 202 using matrix::determinant; 203 204 const T det = determinant(mat); 205 Matrix<T, 2, 2> retVal; 206 207 retVal(0, 0) = mat(1, 1) / det; 208 retVal(0, 1) = -mat(0, 1) / det; 209 retVal(1, 0) = -mat(1, 0) / det; 210 retVal(1, 1) = mat(0, 0) / det; 211 212 return retVal; 213 } 214 215 template <typename T> 216 Matrix<T, 3, 3> SquareMatrixOps<T, 3>::doInverse (const Matrix<T, 3, 3>& mat) 217 { 218 // Blockwise inversion 219 using matrix::inverse; 220 221 const T areaA[2*2] = 222 { 223 mat(0,0), mat(0,1), 224 mat(1,0), mat(1,1) 225 }; 226 const T areaB[2] = 227 { 228 mat(0,2), 229 mat(1,2), 230 }; 231 const T areaC[2] = 232 { 233 mat(2,0), mat(2,1), 234 }; 235 const T areaD[1] = 236 { 237 mat(2,2) 238 }; 239 const T nullField[4] = { T(0.0f) }; 240 241 const Matrix<T, 2, 2> invA = inverse(Matrix<T, 2, 2>(areaA)); 242 const Matrix<T, 2, 1> matB = Matrix<T, 2, 1>(areaB); 243 const Matrix<T, 1, 2> matC = Matrix<T, 1, 2>(areaC); 244 const Matrix<T, 1, 1> matD = Matrix<T, 1, 1>(areaD); 245 246 const T schurComplement = T(1.0f) / (matD - matC*invA*matB)(0,0); 247 const Matrix<T, 2, 2> zeroMat = Matrix<T, 2, 2>(nullField); 248 249 const Matrix<T, 2, 2> blockA = invA + invA*matB*schurComplement*matC*invA; 250 const Matrix<T, 2, 1> blockB = (zeroMat-invA)*matB*schurComplement; 251 const Matrix<T, 1, 2> blockC = matC*invA*(-schurComplement); 252 const T blockD = schurComplement; 253 254 const T result[3*3] = 255 { 256 blockA(0,0), blockA(0,1), blockB(0,0), 257 blockA(1,0), blockA(1,1), blockB(1,0), 258 blockC(0,0), blockC(0,1), blockD, 259 }; 260 261 return Matrix<T, 3, 3>(result); 262 } 263 264 template <typename T> 265 Matrix<T, 4, 4> SquareMatrixOps<T, 4>::doInverse (const Matrix<T, 4, 4>& mat) 266 { 267 // Blockwise inversion 268 using matrix::inverse; 269 270 const T areaA[2*2] = 271 { 272 mat(0,0), mat(0,1), 273 mat(1,0), mat(1,1) 274 }; 275 const T areaB[2*2] = 276 { 277 mat(0,2), mat(0,3), 278 mat(1,2), mat(1,3) 279 }; 280 const T areaC[2*2] = 281 { 282 mat(2,0), mat(2,1), 283 mat(3,0), mat(3,1) 284 }; 285 const T areaD[2*2] = 286 { 287 mat(2,2), mat(2,3), 288 mat(3,2), mat(3,3) 289 }; 290 const T nullField[4] = { T(0.0f) }; 291 292 const Matrix<T, 2, 2> invA = inverse(Matrix<T, 2, 2>(areaA)); 293 const Matrix<T, 2, 2> matB = Matrix<T, 2, 2>(areaB); 294 const Matrix<T, 2, 2> matC = Matrix<T, 2, 2>(areaC); 295 const Matrix<T, 2, 2> matD = Matrix<T, 2, 2>(areaD); 296 297 const Matrix<T, 2, 2> schurComplement = inverse(matD - matC*invA*matB); 298 const Matrix<T, 2, 2> zeroMat = Matrix<T, 2, 2>(nullField); 299 300 const Matrix<T, 2, 2> blockA = invA + invA*matB*schurComplement*matC*invA; 301 const Matrix<T, 2, 2> blockB = (zeroMat-invA)*matB*schurComplement; 302 const Matrix<T, 2, 2> blockC = (zeroMat-schurComplement)*matC*invA; 303 const Matrix<T, 2, 2> blockD = schurComplement; 304 305 const T result[4*4] = 306 { 307 blockA(0,0), blockA(0,1), blockB(0,0), blockB(0,1), 308 blockA(1,0), blockA(1,1), blockB(1,0), blockB(1,1), 309 blockC(0,0), blockC(0,1), blockD(0,0), blockD(0,1), 310 blockC(1,0), blockC(1,1), blockD(1,0), blockD(1,1), 311 }; 312 313 return Matrix<T, 4, 4>(result); 314 } 315 316 // Initialize to identity. 317 template <typename T, int Rows, int Cols> 318 Matrix<T, Rows, Cols>::Matrix (void) 319 { 320 for (int row = 0; row < Rows; row++) 321 for (int col = 0; col < Cols; col++) 322 (*this)(row, col) = (row == col) ? T(1) : T(0); 323 } 324 325 // Initialize to diagonal matrix. 326 template <typename T, int Rows, int Cols> 327 Matrix<T, Rows, Cols>::Matrix (const T& src) 328 { 329 for (int row = 0; row < Rows; row++) 330 for (int col = 0; col < Cols; col++) 331 (*this)(row, col) = (row == col) ? src : T(0); 332 } 333 334 // Initialize from data array. 335 template <typename T, int Rows, int Cols> 336 Matrix<T, Rows, Cols>::Matrix (const T src[Rows*Cols]) 337 { 338 for (int row = 0; row < Rows; row++) 339 for (int col = 0; col < Cols; col++) 340 (*this)(row, col) = src[row*Cols + col]; 341 } 342 343 // Initialize to diagonal matrix. 344 template <typename T, int Rows, int Cols> 345 Matrix<T, Rows, Cols>::Matrix (const Vector<T, Rows>& src) 346 { 347 DE_STATIC_ASSERT(Rows == Cols); 348 for (int row = 0; row < Rows; row++) 349 for (int col = 0; col < Cols; col++) 350 (*this)(row, col) = (row == col) ? src.m_data[row] : T(0); 351 } 352 353 // Copy constructor. 354 template <typename T, int Rows, int Cols> 355 Matrix<T, Rows, Cols>::Matrix (const Matrix<T, Rows, Cols>& src) 356 { 357 *this = src; 358 } 359 360 // Destructor. 361 template <typename T, int Rows, int Cols> 362 Matrix<T, Rows, Cols>::~Matrix (void) 363 { 364 } 365 366 // Assignment operator. 367 template <typename T, int Rows, int Cols> 368 Matrix<T, Rows, Cols>& Matrix<T, Rows, Cols>::operator= (const Matrix<T, Rows, Cols>& src) 369 { 370 for (int row = 0; row < Rows; row++) 371 for (int col = 0; col < Cols; col++) 372 (*this)(row, col) = src(row, col); 373 return *this; 374 } 375 376 // Multipy and assign op 377 template <typename T, int Rows, int Cols> 378 Matrix<T, Rows, Cols>& Matrix<T, Rows, Cols>::operator*= (const Matrix<T, Rows, Cols>& src) 379 { 380 *this = *this * src; 381 return *this; 382 } 383 384 template <typename T, int Rows, int Cols> 385 void Matrix<T, Rows, Cols>::setRow (int rowNdx, const Vector<T, Cols>& vec) 386 { 387 for (int col = 0; col < Cols; col++) 388 (*this)(rowNdx, col) = vec.m_data[col]; 389 } 390 391 template <typename T, int Rows, int Cols> 392 void Matrix<T, Rows, Cols>::setColumn (int colNdx, const Vector<T, Rows>& vec) 393 { 394 m_data[colNdx] = vec; 395 } 396 397 template <typename T, int Rows, int Cols> 398 Vector<T, Cols> Matrix<T, Rows, Cols>::getRow (int rowNdx) const 399 { 400 Vector<T, Cols> res; 401 for (int col = 0; col < Cols; col++) 402 res[col] = (*this)(rowNdx, col); 403 return res; 404 } 405 406 template <typename T, int Rows, int Cols> 407 Vector<T, Rows>& Matrix<T, Rows, Cols>::getColumn (int colNdx) 408 { 409 return m_data[colNdx]; 410 } 411 412 template <typename T, int Rows, int Cols> 413 const Vector<T, Rows>& Matrix<T, Rows, Cols>::getColumn (int colNdx) const 414 { 415 return m_data[colNdx]; 416 } 417 418 template <typename T, int Rows, int Cols> 419 Array<T, Rows*Cols> Matrix<T, Rows, Cols>::getColumnMajorData (void) const 420 { 421 Array<T, Rows*Cols> a; 422 T* dst = a.getPtr(); 423 for (int col = 0; col < Cols; col++) 424 for (int row = 0; row < Rows; row++) 425 *dst++ = (*this)(row, col); 426 return a; 427 } 428 429 template <typename T, int Rows, int Cols> 430 Array<T, Rows*Cols> Matrix<T, Rows, Cols>::getRowMajorData (void) const 431 { 432 Array<T, Rows*Cols> a; 433 T* dst = a.getPtr(); 434 for (int row = 0; row < Rows; row++) 435 for (int col = 0; col < Cols; col++) 436 *dst++ = (*this)(row, col); 437 return a; 438 } 439 440 // Multiplication of two matrices. 441 template <typename T, int Rows0, int Cols0, int Rows1, int Cols1> 442 Matrix<T, Rows0, Cols1> operator* (const Matrix<T, Rows0, Cols0>& a, const Matrix<T, Rows1, Cols1>& b) 443 { 444 DE_STATIC_ASSERT(Cols0 == Rows1); 445 Matrix<T, Rows0, Cols1> res; 446 for (int row = 0; row < Rows0; row++) 447 { 448 for (int col = 0; col < Cols1; col++) 449 { 450 T v = T(0); 451 for (int ndx = 0; ndx < Cols0; ndx++) 452 v += a(row,ndx) * b(ndx,col); 453 res(row,col) = v; 454 } 455 } 456 return res; 457 } 458 459 // Multiply of matrix with column vector. 460 template <typename T, int Rows, int Cols> 461 Vector<T, Rows> operator* (const Matrix<T, Rows, Cols>& mtx, const Vector<T, Cols>& vec) 462 { 463 Vector<T, Rows> res; 464 for (int row = 0; row < Rows; row++) 465 { 466 T v = T(0); 467 for (int col = 0; col < Cols; col++) 468 v += mtx(row,col) * vec.m_data[col]; 469 res.m_data[row] = v; 470 } 471 return res; 472 } 473 474 // Multiply of matrix with row vector. 475 template <typename T, int Rows, int Cols> 476 Vector<T, Cols> operator* (const Vector<T, Rows>& vec, const Matrix<T, Rows, Cols>& mtx) 477 { 478 Vector<T, Cols> res; 479 for (int col = 0; col < Cols; col++) 480 { 481 T v = T(0); 482 for (int row = 0; row < Rows; row++) 483 v += mtx(row,col) * vec.m_data[row]; 484 res.m_data[col] = v; 485 } 486 return res; 487 } 488 489 // Common typedefs. 490 typedef Matrix<float, 2, 2> Matrix2f; 491 typedef Matrix<float, 3, 3> Matrix3f; 492 typedef Matrix<float, 4, 4> Matrix4f; 493 typedef Matrix<double, 2, 2> Matrix2d; 494 typedef Matrix<double, 3, 3> Matrix3d; 495 typedef Matrix<double, 4, 4> Matrix4d; 496 497 // GLSL-style naming \note CxR. 498 typedef Matrix2f Mat2; 499 typedef Matrix<float, 3, 2> Mat2x3; 500 typedef Matrix<float, 4, 2> Mat2x4; 501 typedef Matrix<float, 2, 3> Mat3x2; 502 typedef Matrix3f Mat3; 503 typedef Matrix<float, 4, 3> Mat3x4; 504 typedef Matrix<float, 2, 4> Mat4x2; 505 typedef Matrix<float, 3, 4> Mat4x3; 506 typedef Matrix4f Mat4; 507 508 // Matrix-scalar operators. 509 510 template <typename T, int Rows, int Cols> 511 Matrix<T, Rows, Cols> operator+ (const Matrix<T, Rows, Cols>& mtx, T scalar) 512 { 513 Matrix<T, Rows, Cols> res; 514 for (int col = 0; col < Cols; col++) 515 for (int row = 0; row < Rows; row++) 516 res(row, col) = mtx(row, col) + scalar; 517 return res; 518 } 519 520 template <typename T, int Rows, int Cols> 521 Matrix<T, Rows, Cols> operator- (const Matrix<T, Rows, Cols>& mtx, T scalar) 522 { 523 Matrix<T, Rows, Cols> res; 524 for (int col = 0; col < Cols; col++) 525 for (int row = 0; row < Rows; row++) 526 res(row, col) = mtx(row, col) - scalar; 527 return res; 528 } 529 530 template <typename T, int Rows, int Cols> 531 Matrix<T, Rows, Cols> operator* (const Matrix<T, Rows, Cols>& mtx, T scalar) 532 { 533 Matrix<T, Rows, Cols> res; 534 for (int col = 0; col < Cols; col++) 535 for (int row = 0; row < Rows; row++) 536 res(row, col) = mtx(row, col) * scalar; 537 return res; 538 } 539 540 template <typename T, int Rows, int Cols> 541 Matrix<T, Rows, Cols> operator/ (const Matrix<T, Rows, Cols>& mtx, T scalar) 542 { 543 Matrix<T, Rows, Cols> res; 544 for (int col = 0; col < Cols; col++) 545 for (int row = 0; row < Rows; row++) 546 res(row, col) = mtx(row, col) / scalar; 547 return res; 548 } 549 550 // Matrix-matrix component-wise operators. 551 552 template <typename T, int Rows, int Cols> 553 Matrix<T, Rows, Cols> operator+ (const Matrix<T, Rows, Cols>& a, const Matrix<T, Rows, Cols>& b) 554 { 555 Matrix<T, Rows, Cols> res; 556 for (int col = 0; col < Cols; col++) 557 for (int row = 0; row < Rows; row++) 558 res(row, col) = a(row, col) + b(row, col); 559 return res; 560 } 561 562 template <typename T, int Rows, int Cols> 563 Matrix<T, Rows, Cols> operator- (const Matrix<T, Rows, Cols>& a, const Matrix<T, Rows, Cols>& b) 564 { 565 Matrix<T, Rows, Cols> res; 566 for (int col = 0; col < Cols; col++) 567 for (int row = 0; row < Rows; row++) 568 res(row, col) = a(row, col) - b(row, col); 569 return res; 570 } 571 572 template <typename T, int Rows, int Cols> 573 Matrix<T, Rows, Cols> operator/ (const Matrix<T, Rows, Cols>& a, const Matrix<T, Rows, Cols>& b) 574 { 575 Matrix<T, Rows, Cols> res; 576 for (int col = 0; col < Cols; col++) 577 for (int row = 0; row < Rows; row++) 578 res(row, col) = a(row, col) / b(row, col); 579 return res; 580 } 581 582 template <typename T, int Rows, int Cols> 583 bool operator== (const Matrix<T, Rows, Cols>& lhs, const Matrix<T, Rows, Cols>& rhs) 584 { 585 for (int row = 0; row < Rows; row++) 586 for (int col = 0; col < Cols; col++) 587 if (lhs(row, col) != rhs(row, col)) 588 return false; 589 return true; 590 } 591 592 template <typename T, int Rows, int Cols> 593 bool operator!= (const Matrix<T, Rows, Cols>& lhs, const Matrix<T, Rows, Cols>& rhs) 594 { 595 return !(lhs == rhs); 596 } 597 598 } // tcu 599 600 #endif // _TCUMATRIX_HPP 601