matrix.cpp 4.1 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207
  1. /*
  2. * DoRayMe - a quick and dirty Raytracer
  3. * Matrix implementation
  4. *
  5. * Created by Manoël Trapier
  6. * Copyright (c) 2020 986-Studio.
  7. *
  8. */
  9. #include <stdio.h>
  10. #include <matrix.h>
  11. #include <tuple.h>
  12. #include <math_helper.h>
  13. Matrix::Matrix(int width)
  14. {
  15. int i;
  16. this->size = width;
  17. for(i = 0; i < width*width; i++)
  18. {
  19. this->data[i] = 0;
  20. }
  21. };
  22. Matrix::Matrix(double values[], int width)
  23. {
  24. int x, y;
  25. this->size = width;
  26. for(y = 0; y < this->size; y++)
  27. {
  28. for (x = 0 ; x < this->size ; x++)
  29. {
  30. this->data[this->size * x + y] = values[this->size * x + y];
  31. }
  32. }
  33. };
  34. bool Matrix::operator==(const Matrix &b) const
  35. {
  36. int i;
  37. if (this->size != b.size)
  38. {
  39. /* If they are not the same size don't even bother */
  40. return false;
  41. }
  42. for(i = 0; i < this->size*this->size; i++)
  43. {
  44. if (!double_equal(this->data[i], b.data[i]))
  45. {
  46. return false;
  47. }
  48. }
  49. return true;
  50. }
  51. bool Matrix::operator!=(const Matrix &b) const
  52. {
  53. int i;
  54. if (this->size != b.size)
  55. {
  56. /* If they are not the same size don't even bother */
  57. return true;
  58. }
  59. for(i = 0; i < this->size*this->size; i++)
  60. {
  61. if (!double_equal(this->data[i], b.data[i]))
  62. {
  63. return true;
  64. }
  65. }
  66. return false;
  67. }
  68. Matrix Matrix::operator*(const Matrix &b) const
  69. {
  70. int x, y, k;
  71. Matrix ret = Matrix(this->size);
  72. if (this->size == b.size)
  73. {
  74. for (y = 0 ; y < this->size ; y++)
  75. {
  76. for (x = 0 ; x < this->size ; x++)
  77. {
  78. double v = 0;
  79. for (k = 0 ; k < this->size ; k++)
  80. {
  81. v += this->get(x, k) * b.get(k, y);
  82. }
  83. ret.set(x, y, v);
  84. }
  85. }
  86. }
  87. return ret;
  88. }
  89. //#define FastGet4(_x, _y) (this->data[4 * (_x) + (_y)])
  90. /* TODO: Check if we can optimise this function. It is called a lot */
  91. /*
  92. Tuple Matrix::operator*(const Tuple &b) const
  93. {
  94. return Tuple(b.x * FastGet4(0, 0) + b.y * FastGet4(0, 1) + b.z * FastGet4(0, 2) + b.w * FastGet4(0, 3),
  95. b.x * FastGet4(1, 0) + b.y * FastGet4(1, 1) + b.z * FastGet4(1, 2) + b.w * FastGet4(1, 3),
  96. b.x * FastGet4(2, 0) + b.y * FastGet4(2, 1) + b.z * FastGet4(2, 2) + b.w * FastGet4(2, 3),
  97. b.x * FastGet4(3, 0) + b.y * FastGet4(3, 1) + b.z * FastGet4(3, 2) + b.w * FastGet4(3, 3));
  98. }
  99. */
  100. Matrix Matrix::identity()
  101. {
  102. int i;
  103. for(i = 0; i < this->size; i++)
  104. {
  105. this->set(i, i, 1);
  106. }
  107. return *this;
  108. }
  109. Matrix Matrix::transpose()
  110. {
  111. int x, y;
  112. Matrix ret = Matrix(this->size);
  113. for (y = 0 ; y < this->size ; y++)
  114. {
  115. for (x = 0 ; x < this->size ; x++)
  116. {
  117. ret.set(y, x, this->get(x, y));
  118. }
  119. }
  120. return ret;
  121. }
  122. Matrix Matrix::submatrix(int row, int column)
  123. {
  124. int i, j;
  125. int x = 0, y = 0;
  126. Matrix ret = Matrix(this->size - 1);
  127. for (i = 0 ; i < this->size ; i++)
  128. {
  129. if (i == row)
  130. {
  131. /* Skip that row */
  132. continue;
  133. }
  134. y = 0;
  135. for (j = 0 ; j < this->size ; j++)
  136. {
  137. if (j == column)
  138. {
  139. /* skip that column */
  140. continue;
  141. }
  142. ret.set(x, y, this->get(i, j));
  143. y++;
  144. }
  145. x++;
  146. }
  147. return ret;
  148. }
  149. double Matrix::determinant()
  150. {
  151. double det = 0;
  152. if (this->size == 2)
  153. {
  154. det = this->data[0] * this->data[3] - this->data[1] * this->data[2];
  155. }
  156. else
  157. {
  158. int col;
  159. for(col = 0; col < this->size; col++)
  160. {
  161. det = det + this->get(0, col) * this->cofactor(0, col);
  162. }
  163. }
  164. return det;
  165. }
  166. Matrix Matrix::inverse()
  167. {
  168. Matrix ret = Matrix(this->size);
  169. if (this->isInvertible())
  170. {
  171. int row, col;
  172. double c;
  173. for (row = 0; row < this->size; row++)
  174. {
  175. for (col = 0; col < this->size; col++)
  176. {
  177. c = this->cofactor(row, col);
  178. ret.set(col, row, c / this->determinant());
  179. }
  180. }
  181. }
  182. return ret;
  183. }