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 30namespace tcu 31{ 32 33// Templated matrix class. 34template <typename T, int Rows, int Cols> 35class Matrix 36{ 37public: 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 74private: 75 Vector<Vector<T, Rows>, Cols> m_data; 76} DE_WARN_UNUSED_TYPE; 77 78// Operators. 79 80// Mat * Mat. 81template <typename T, int Rows0, int Cols0, int Rows1, int Cols1> 82Matrix<T, Rows0, Cols1> operator* (const Matrix<T, Rows0, Cols0>& a, const Matrix<T, Rows1, Cols1>& b); 83 84// Mat * Vec (column vector). 85template <typename T, int Rows, int Cols> 86Vector<T, Rows> operator* (const Matrix<T, Rows, Cols>& mtx, const Vector<T, Cols>& vec); 87 88// Vec * Mat (row vector). 89template <typename T, int Rows, int Cols> 90Vector<T, Cols> operator* (const Vector<T, Rows>& vec, const Matrix<T, Rows, Cols>& mtx); 91 92// Further operations 93 94template <typename T, int Size> 95struct SquareMatrixOps 96{ 97 static T doDeterminant (const Matrix<T, Size, Size>& mat); 98 static Matrix<T, Size, Size> doInverse (const Matrix<T, Size, Size>& mat); 99}; 100 101template <typename T> 102struct SquareMatrixOps<T, 2> 103{ 104 static T doDeterminant (const Matrix<T, 2, 2>& mat); 105 static Matrix<T, 2, 2> doInverse (const Matrix<T, 2, 2>& mat); 106}; 107 108template <typename T> 109struct SquareMatrixOps<T, 3> 110{ 111 static T doDeterminant (const Matrix<T, 3, 3>& mat); 112 static Matrix<T, 3, 3> doInverse (const Matrix<T, 3, 3>& mat); 113}; 114 115template <typename T> 116struct SquareMatrixOps<T, 4> 117{ 118 static T doDeterminant (const Matrix<T, 4, 4>& mat); 119 static Matrix<T, 4, 4> doInverse (const Matrix<T, 4, 4>& mat); 120}; 121 122namespace matrix 123{ 124 125template <typename T, int Size> 126T determinant (const Matrix<T, Size, Size>& mat) 127{ 128 return SquareMatrixOps<T, Size>::doDeterminant(mat); 129} 130 131template <typename T, int Size> 132Matrix<T, Size, Size> inverse (const Matrix<T, Size, Size>& mat) 133{ 134 return SquareMatrixOps<T, Size>::doInverse(mat); 135} 136 137} // matrix 138 139// Template implementations. 140 141template <typename T> 142T SquareMatrixOps<T, 2>::doDeterminant (const Matrix<T, 2, 2>& mat) 143{ 144 return mat(0,0) * mat(1,1) - mat(1,0) * mat(0,1); 145} 146 147template <typename T> 148T SquareMatrixOps<T, 3>::doDeterminant (const Matrix<T, 3, 3>& mat) 149{ 150 return + mat(0,0) * mat(1,1) * mat(2,2) 151 + mat(0,1) * mat(1,2) * mat(2,0) 152 + mat(0,2) * mat(1,0) * mat(2,1) 153 - mat(0,0) * mat(1,2) * mat(2,1) 154 - mat(0,1) * mat(1,0) * mat(2,2) 155 - mat(0,2) * mat(1,1) * mat(2,0); 156} 157 158template <typename T> 159T SquareMatrixOps<T, 4>::doDeterminant (const Matrix<T, 4, 4>& mat) 160{ 161 using matrix::determinant; 162 163 const T minorMatrices[4][3*3] = 164 { 165 { 166 mat(1,1), mat(2,1), mat(3,1), 167 mat(1,2), mat(2,2), mat(3,2), 168 mat(1,3), mat(2,3), mat(3,3), 169 }, 170 { 171 mat(1,0), mat(2,0), mat(3,0), 172 mat(1,2), mat(2,2), mat(3,2), 173 mat(1,3), mat(2,3), mat(3,3), 174 }, 175 { 176 mat(1,0), mat(2,0), mat(3,0), 177 mat(1,1), mat(2,1), mat(3,1), 178 mat(1,3), mat(2,3), mat(3,3), 179 }, 180 { 181 mat(1,0), mat(2,0), mat(3,0), 182 mat(1,1), mat(2,1), mat(3,1), 183 mat(1,2), mat(2,2), mat(3,2), 184 } 185 }; 186 187 return + mat(0,0) * determinant(Matrix<T, 3, 3>(minorMatrices[0])) 188 - mat(0,1) * determinant(Matrix<T, 3, 3>(minorMatrices[1])) 189 + mat(0,2) * determinant(Matrix<T, 3, 3>(minorMatrices[2])) 190 - mat(0,3) * determinant(Matrix<T, 3, 3>(minorMatrices[3])); 191} 192 193template <typename T> 194Matrix<T, 2, 2> SquareMatrixOps<T, 2>::doInverse (const Matrix<T, 2, 2>& mat) 195{ 196 using matrix::determinant; 197 198 const T det = determinant(mat); 199 Matrix<T, 2, 2> retVal; 200 201 retVal(0, 0) = mat(1, 1) / det; 202 retVal(0, 1) = -mat(0, 1) / det; 203 retVal(1, 0) = -mat(1, 0) / det; 204 retVal(1, 1) = mat(0, 0) / det; 205 206 return retVal; 207} 208 209template <typename T> 210Matrix<T, 3, 3> SquareMatrixOps<T, 3>::doInverse (const Matrix<T, 3, 3>& mat) 211{ 212 // Blockwise inversion 213 using matrix::inverse; 214 215 const T areaA[2*2] = 216 { 217 mat(0,0), mat(0,1), 218 mat(1,0), mat(1,1) 219 }; 220 const T areaB[2] = 221 { 222 mat(0,2), 223 mat(1,2), 224 }; 225 const T areaC[2] = 226 { 227 mat(2,0), mat(2,1), 228 }; 229 const T areaD[1] = 230 { 231 mat(2,2) 232 }; 233 const T nullField[4] = { T(0.0f) }; 234 235 const Matrix<T, 2, 2> invA = inverse(Matrix<T, 2, 2>(areaA)); 236 const Matrix<T, 2, 1> matB = Matrix<T, 2, 1>(areaB); 237 const Matrix<T, 1, 2> matC = Matrix<T, 1, 2>(areaC); 238 const Matrix<T, 1, 1> matD = Matrix<T, 1, 1>(areaD); 239 240 const T schurComplement = T(1.0f) / (matD - matC*invA*matB)(0,0); 241 const Matrix<T, 2, 2> zeroMat = Matrix<T, 2, 2>(nullField); 242 243 const Matrix<T, 2, 2> blockA = invA + invA*matB*schurComplement*matC*invA; 244 const Matrix<T, 2, 1> blockB = (zeroMat-invA)*matB*schurComplement; 245 const Matrix<T, 1, 2> blockC = matC*invA*(-schurComplement); 246 const T blockD = schurComplement; 247 248 const T result[3*3] = 249 { 250 blockA(0,0), blockA(0,1), blockB(0,0), 251 blockA(1,0), blockA(1,1), blockB(1,0), 252 blockC(0,0), blockC(0,1), blockD, 253 }; 254 255 return Matrix<T, 3, 3>(result); 256} 257 258template <typename T> 259Matrix<T, 4, 4> SquareMatrixOps<T, 4>::doInverse (const Matrix<T, 4, 4>& mat) 260{ 261 // Blockwise inversion 262 using matrix::inverse; 263 264 const T areaA[2*2] = 265 { 266 mat(0,0), mat(0,1), 267 mat(1,0), mat(1,1) 268 }; 269 const T areaB[2*2] = 270 { 271 mat(0,2), mat(0,3), 272 mat(1,2), mat(1,3) 273 }; 274 const T areaC[2*2] = 275 { 276 mat(2,0), mat(2,1), 277 mat(3,0), mat(3,1) 278 }; 279 const T areaD[2*2] = 280 { 281 mat(2,2), mat(2,3), 282 mat(3,2), mat(3,3) 283 }; 284 const T nullField[4] = { T(0.0f) }; 285 286 const Matrix<T, 2, 2> invA = inverse(Matrix<T, 2, 2>(areaA)); 287 const Matrix<T, 2, 2> matB = Matrix<T, 2, 2>(areaB); 288 const Matrix<T, 2, 2> matC = Matrix<T, 2, 2>(areaC); 289 const Matrix<T, 2, 2> matD = Matrix<T, 2, 2>(areaD); 290 291 const Matrix<T, 2, 2> schurComplement = inverse(matD - matC*invA*matB); 292 const Matrix<T, 2, 2> zeroMat = Matrix<T, 2, 2>(nullField); 293 294 const Matrix<T, 2, 2> blockA = invA + invA*matB*schurComplement*matC*invA; 295 const Matrix<T, 2, 2> blockB = (zeroMat-invA)*matB*schurComplement; 296 const Matrix<T, 2, 2> blockC = (zeroMat-schurComplement)*matC*invA; 297 const Matrix<T, 2, 2> blockD = schurComplement; 298 299 const T result[4*4] = 300 { 301 blockA(0,0), blockA(0,1), blockB(0,0), blockB(0,1), 302 blockA(1,0), blockA(1,1), blockB(1,0), blockB(1,1), 303 blockC(0,0), blockC(0,1), blockD(0,0), blockD(0,1), 304 blockC(1,0), blockC(1,1), blockD(1,0), blockD(1,1), 305 }; 306 307 return Matrix<T, 4, 4>(result); 308} 309 310// Initialize to identity. 311template <typename T, int Rows, int Cols> 312Matrix<T, Rows, Cols>::Matrix (void) 313{ 314 for (int row = 0; row < Rows; row++) 315 for (int col = 0; col < Cols; col++) 316 (*this)(row, col) = (row == col) ? T(1) : T(0); 317} 318 319// Initialize to diagonal matrix. 320template <typename T, int Rows, int Cols> 321Matrix<T, Rows, Cols>::Matrix (const T& src) 322{ 323 for (int row = 0; row < Rows; row++) 324 for (int col = 0; col < Cols; col++) 325 (*this)(row, col) = (row == col) ? src : T(0); 326} 327 328// Initialize from data array. 329template <typename T, int Rows, int Cols> 330Matrix<T, Rows, Cols>::Matrix (const T src[Rows*Cols]) 331{ 332 for (int row = 0; row < Rows; row++) 333 for (int col = 0; col < Cols; col++) 334 (*this)(row, col) = src[row*Cols + col]; 335} 336 337// Initialize to diagonal matrix. 338template <typename T, int Rows, int Cols> 339Matrix<T, Rows, Cols>::Matrix (const Vector<T, Rows>& src) 340{ 341 DE_STATIC_ASSERT(Rows == Cols); 342 for (int row = 0; row < Rows; row++) 343 for (int col = 0; col < Cols; col++) 344 (*this)(row, col) = (row == col) ? src.m_data[row] : T(0); 345} 346 347// Copy constructor. 348template <typename T, int Rows, int Cols> 349Matrix<T, Rows, Cols>::Matrix (const Matrix<T, Rows, Cols>& src) 350{ 351 *this = src; 352} 353 354// Destructor. 355template <typename T, int Rows, int Cols> 356Matrix<T, Rows, Cols>::~Matrix (void) 357{ 358} 359 360// Assignment operator. 361template <typename T, int Rows, int Cols> 362Matrix<T, Rows, Cols>& Matrix<T, Rows, Cols>::operator= (const Matrix<T, Rows, Cols>& src) 363{ 364 for (int row = 0; row < Rows; row++) 365 for (int col = 0; col < Cols; col++) 366 (*this)(row, col) = src(row, col); 367 return *this; 368} 369 370// Multipy and assign op 371template <typename T, int Rows, int Cols> 372Matrix<T, Rows, Cols>& Matrix<T, Rows, Cols>::operator*= (const Matrix<T, Rows, Cols>& src) 373{ 374 *this = *this * src; 375 return *this; 376} 377 378template <typename T, int Rows, int Cols> 379void Matrix<T, Rows, Cols>::setRow (int rowNdx, const Vector<T, Cols>& vec) 380{ 381 for (int col = 0; col < Cols; col++) 382 (*this)(rowNdx, col) = vec.m_data[col]; 383} 384 385template <typename T, int Rows, int Cols> 386void Matrix<T, Rows, Cols>::setColumn (int colNdx, const Vector<T, Rows>& vec) 387{ 388 m_data[colNdx] = vec; 389} 390 391template <typename T, int Rows, int Cols> 392Vector<T, Cols> Matrix<T, Rows, Cols>::getRow (int rowNdx) const 393{ 394 Vector<T, Cols> res; 395 for (int col = 0; col < Cols; col++) 396 res[col] = (*this)(rowNdx, col); 397 return res; 398} 399 400template <typename T, int Rows, int Cols> 401Vector<T, Rows>& Matrix<T, Rows, Cols>::getColumn (int colNdx) 402{ 403 return m_data[colNdx]; 404} 405 406template <typename T, int Rows, int Cols> 407const Vector<T, Rows>& Matrix<T, Rows, Cols>::getColumn (int colNdx) const 408{ 409 return m_data[colNdx]; 410} 411 412template <typename T, int Rows, int Cols> 413Array<T, Rows*Cols> Matrix<T, Rows, Cols>::getColumnMajorData (void) const 414{ 415 Array<T, Rows*Cols> a; 416 T* dst = a.getPtr(); 417 for (int col = 0; col < Cols; col++) 418 for (int row = 0; row < Rows; row++) 419 *dst++ = (*this)(row, col); 420 return a; 421} 422 423template <typename T, int Rows, int Cols> 424Array<T, Rows*Cols> Matrix<T, Rows, Cols>::getRowMajorData (void) const 425{ 426 Array<T, Rows*Cols> a; 427 T* dst = a.getPtr(); 428 for (int row = 0; row < Rows; row++) 429 for (int col = 0; col < Cols; col++) 430 *dst++ = (*this)(row, col); 431 return a; 432} 433 434// Multiplication of two matrices. 435template <typename T, int Rows0, int Cols0, int Rows1, int Cols1> 436Matrix<T, Rows0, Cols1> operator* (const Matrix<T, Rows0, Cols0>& a, const Matrix<T, Rows1, Cols1>& b) 437{ 438 DE_STATIC_ASSERT(Cols0 == Rows1); 439 Matrix<T, Rows0, Cols1> res; 440 for (int row = 0; row < Rows0; row++) 441 { 442 for (int col = 0; col < Cols1; col++) 443 { 444 T v = T(0); 445 for (int ndx = 0; ndx < Cols0; ndx++) 446 v += a(row,ndx) * b(ndx,col); 447 res(row,col) = v; 448 } 449 } 450 return res; 451} 452 453// Multiply of matrix with column vector. 454template <typename T, int Rows, int Cols> 455Vector<T, Rows> operator* (const Matrix<T, Rows, Cols>& mtx, const Vector<T, Cols>& vec) 456{ 457 Vector<T, Rows> res; 458 for (int row = 0; row < Rows; row++) 459 { 460 T v = T(0); 461 for (int col = 0; col < Cols; col++) 462 v += mtx(row,col) * vec.m_data[col]; 463 res.m_data[row] = v; 464 } 465 return res; 466} 467 468// Multiply of matrix with row vector. 469template <typename T, int Rows, int Cols> 470Vector<T, Cols> operator* (const Vector<T, Rows>& vec, const Matrix<T, Rows, Cols>& mtx) 471{ 472 Vector<T, Cols> res; 473 for (int col = 0; col < Cols; col++) 474 { 475 T v = T(0); 476 for (int row = 0; row < Rows; row++) 477 v += mtx(row,col) * vec.m_data[row]; 478 res.m_data[col] = v; 479 } 480 return res; 481} 482 483// Common typedefs. 484typedef Matrix<float, 2, 2> Matrix2f; 485typedef Matrix<float, 3, 3> Matrix3f; 486typedef Matrix<float, 4, 4> Matrix4f; 487typedef Matrix<double, 2, 2> Matrix2d; 488typedef Matrix<double, 3, 3> Matrix3d; 489typedef Matrix<double, 4, 4> Matrix4d; 490 491// GLSL-style naming \note CxR. 492typedef Matrix2f Mat2; 493typedef Matrix<float, 3, 2> Mat2x3; 494typedef Matrix<float, 4, 2> Mat2x4; 495typedef Matrix<float, 2, 3> Mat3x2; 496typedef Matrix3f Mat3; 497typedef Matrix<float, 4, 3> Mat3x4; 498typedef Matrix<float, 2, 4> Mat4x2; 499typedef Matrix<float, 3, 4> Mat4x3; 500typedef Matrix4f Mat4; 501 502// Matrix-scalar operators. 503 504template <typename T, int Rows, int Cols> 505Matrix<T, Rows, Cols> operator+ (const Matrix<T, Rows, Cols>& mtx, T scalar) 506{ 507 Matrix<T, Rows, Cols> res; 508 for (int col = 0; col < Cols; col++) 509 for (int row = 0; row < Rows; row++) 510 res(row, col) = mtx(row, col) + scalar; 511 return res; 512} 513 514template <typename T, int Rows, int Cols> 515Matrix<T, Rows, Cols> operator- (const Matrix<T, Rows, Cols>& mtx, T scalar) 516{ 517 Matrix<T, Rows, Cols> res; 518 for (int col = 0; col < Cols; col++) 519 for (int row = 0; row < Rows; row++) 520 res(row, col) = mtx(row, col) - scalar; 521 return res; 522} 523 524template <typename T, int Rows, int Cols> 525Matrix<T, Rows, Cols> operator* (const Matrix<T, Rows, Cols>& mtx, T scalar) 526{ 527 Matrix<T, Rows, Cols> res; 528 for (int col = 0; col < Cols; col++) 529 for (int row = 0; row < Rows; row++) 530 res(row, col) = mtx(row, col) * scalar; 531 return res; 532} 533 534template <typename T, int Rows, int Cols> 535Matrix<T, Rows, Cols> operator/ (const Matrix<T, Rows, Cols>& mtx, T scalar) 536{ 537 Matrix<T, Rows, Cols> res; 538 for (int col = 0; col < Cols; col++) 539 for (int row = 0; row < Rows; row++) 540 res(row, col) = mtx(row, col) / scalar; 541 return res; 542} 543 544// Matrix-matrix component-wise operators. 545 546template <typename T, int Rows, int Cols> 547Matrix<T, Rows, Cols> operator+ (const Matrix<T, Rows, Cols>& a, const Matrix<T, Rows, Cols>& b) 548{ 549 Matrix<T, Rows, Cols> res; 550 for (int col = 0; col < Cols; col++) 551 for (int row = 0; row < Rows; row++) 552 res(row, col) = a(row, col) + b(row, col); 553 return res; 554} 555 556template <typename T, int Rows, int Cols> 557Matrix<T, Rows, Cols> operator- (const Matrix<T, Rows, Cols>& a, const Matrix<T, Rows, Cols>& b) 558{ 559 Matrix<T, Rows, Cols> res; 560 for (int col = 0; col < Cols; col++) 561 for (int row = 0; row < Rows; row++) 562 res(row, col) = a(row, col) - b(row, col); 563 return res; 564} 565 566template <typename T, int Rows, int Cols> 567Matrix<T, Rows, Cols> operator/ (const Matrix<T, Rows, Cols>& a, const Matrix<T, Rows, Cols>& b) 568{ 569 Matrix<T, Rows, Cols> res; 570 for (int col = 0; col < Cols; col++) 571 for (int row = 0; row < Rows; row++) 572 res(row, col) = a(row, col) / b(row, col); 573 return res; 574} 575 576} // tcu 577 578#endif // _TCUMATRIX_HPP 579