matrix4.c 5.6 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234
  1. /******************************************************************************
  2. Copyright (C) 2013 by Hugh Bailey <[email protected]>
  3. This program is free software: you can redistribute it and/or modify
  4. it under the terms of the GNU General Public License as published by
  5. the Free Software Foundation, either version 2 of the License, or
  6. (at your option) any later version.
  7. This program is distributed in the hope that it will be useful,
  8. but WITHOUT ANY WARRANTY; without even the implied warranty of
  9. MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
  10. GNU General Public License for more details.
  11. You should have received a copy of the GNU General Public License
  12. along with this program. If not, see <http://www.gnu.org/licenses/>.
  13. ******************************************************************************/
  14. #include "math-defs.h"
  15. #include "matrix4.h"
  16. #include "matrix3.h"
  17. #include "quat.h"
  18. void matrix4_from_matrix3(struct matrix4 *dst, const struct matrix3 *m)
  19. {
  20. dst->x.m = m->x.m;
  21. dst->y.m = m->y.m;
  22. dst->z.m = m->z.m;
  23. dst->t.m = m->t.m;
  24. dst->t.w = 1.0f;
  25. }
  26. void matrix4_from_quat(struct matrix4 *dst, const struct quat *q)
  27. {
  28. float norm = quat_dot(q, q);
  29. float s = (norm > 0.0f) ? (2.0f/norm) : 0.0f;
  30. float xx = q->x * q->x * s;
  31. float yy = q->y * q->y * s;
  32. float zz = q->z * q->z * s;
  33. float xy = q->x * q->y * s;
  34. float xz = q->x * q->z * s;
  35. float yz = q->y * q->z * s;
  36. float wx = q->w * q->x * s;
  37. float wy = q->w * q->y * s;
  38. float wz = q->w * q->z * s;
  39. vec4_set(&dst->x, 1.0f - (yy + zz), xy + wz, xz - wy, 0.0f);
  40. vec4_set(&dst->y, xy - wz, 1.0f - (xx + zz), yz + wx, 0.0f);
  41. vec4_set(&dst->z, xz + wy, yz - wx, 1.0f - (xx + yy), 0.0f);
  42. vec4_set(&dst->t, 0.0f, 0.0f, 0.0f, 1.0f);
  43. }
  44. void matrix4_from_axisang(struct matrix4 *dst, const struct axisang *aa)
  45. {
  46. struct quat q;
  47. quat_from_axisang(&q, aa);
  48. matrix4_from_quat(dst, &q);
  49. }
  50. void matrix4_mul(struct matrix4 *dst, const struct matrix4 *m1,
  51. const struct matrix4 *m2)
  52. {
  53. const struct vec4 *m1v = (const struct vec4*)m1;
  54. const float *m2f = (const float*)m2;
  55. struct vec4 out[4];
  56. int i, j;
  57. for (i = 0; i < 4; i++) {
  58. for (j=0; j<4; j++) {
  59. struct vec4 temp;
  60. vec4_set(&temp, m2f[j], m2f[j+4], m2f[j+8], m2f[j+12]);
  61. out[i].ptr[j] = vec4_dot(&m1v[i], &temp);
  62. }
  63. }
  64. matrix4_copy(dst, (struct matrix4*)out);
  65. }
  66. static inline void get_3x3_submatrix(float *dst, const struct matrix4 *m,
  67. int i, int j)
  68. {
  69. const float *mf = (const float *)m;
  70. int ti, tj, idst, jdst;
  71. for (ti = 0; ti < 4; ti++) {
  72. if (ti < i)
  73. idst = ti;
  74. else if (ti > i)
  75. idst = ti-1;
  76. else
  77. continue;
  78. for (tj = 0; tj < 4; tj++) {
  79. if (tj < j)
  80. jdst = tj;
  81. else if (tj > j)
  82. jdst = tj-1;
  83. else
  84. continue;
  85. dst[(idst*3) + jdst] = mf[(ti*4) + tj];
  86. }
  87. }
  88. }
  89. static inline float get_3x3_determinant(const float *m)
  90. {
  91. return (m[0] * ((m[4]*m[8]) - (m[7]*m[5]))) -
  92. (m[1] * ((m[3]*m[8]) - (m[6]*m[5]))) +
  93. (m[2] * ((m[3]*m[7]) - (m[6]*m[4])));
  94. }
  95. float matrix4_determinant(const struct matrix4 *m)
  96. {
  97. const float *mf = (const float *)m;
  98. float det, result = 0.0f, i = 1.0f;
  99. float m3x3[9];
  100. int n;
  101. for (n = 0; n < 4; n++, i *= -1.0f) {
  102. get_3x3_submatrix(m3x3, m, 0, n);
  103. det = get_3x3_determinant(m3x3);
  104. result += mf[n] * det * i;
  105. }
  106. return result;
  107. }
  108. void matrix4_translate3v(struct matrix4 *dst, const struct matrix4 *m,
  109. const struct vec3 *v)
  110. {
  111. struct matrix4 temp;
  112. vec4_set(&temp.x, 1.0f, 0.0f, 0.0f, 0.0f);
  113. vec4_set(&temp.y, 0.0f, 1.0f, 0.0f, 0.0f);
  114. vec4_set(&temp.z, 0.0f, 0.0f, 1.0f, 0.0f);
  115. vec4_from_vec3(&temp.t, v);
  116. matrix4_mul(dst, m, &temp);
  117. }
  118. void matrix4_translate4v(struct matrix4 *dst, const struct matrix4 *m,
  119. const struct vec4 *v)
  120. {
  121. struct matrix4 temp;
  122. vec4_set(&temp.x, 1.0f, 0.0f, 0.0f, 0.0f);
  123. vec4_set(&temp.y, 0.0f, 1.0f, 0.0f, 0.0f);
  124. vec4_set(&temp.z, 0.0f, 0.0f, 1.0f, 0.0f);
  125. vec4_copy(&temp.t, v);
  126. matrix4_mul(dst, m, &temp);
  127. }
  128. void matrix4_rotate(struct matrix4 *dst, const struct matrix4 *m,
  129. const struct quat *q)
  130. {
  131. struct matrix4 temp;
  132. matrix4_from_quat(&temp, q);
  133. matrix4_mul(dst, m, &temp);
  134. }
  135. void matrix4_rotate_aa(struct matrix4 *dst, const struct matrix4 *m,
  136. const struct axisang *aa)
  137. {
  138. struct matrix4 temp;
  139. matrix4_from_axisang(&temp, aa);
  140. matrix4_mul(dst, m, &temp);
  141. }
  142. void matrix4_scale(struct matrix4 *dst, const struct matrix4 *m,
  143. const struct vec3 *v)
  144. {
  145. struct matrix4 temp;
  146. vec4_set(&temp.x, v->x, 0.0f, 0.0f, 0.0f);
  147. vec4_set(&temp.y, 0.0f, v->y, 0.0f, 0.0f);
  148. vec4_set(&temp.z, 0.0f, 0.0f, v->z, 0.0f);
  149. vec4_set(&temp.t, 0.0f, 0.0f, 0.0f, 1.0f);
  150. matrix4_mul(dst, m, &temp);
  151. }
  152. bool matrix4_inv(struct matrix4 *dst, const struct matrix4 *m)
  153. {
  154. struct vec4 *dstv = (struct vec4 *)dst;
  155. float det = matrix4_determinant(m);
  156. float m3x3[9];
  157. int i, j, sign;
  158. if (dst == m) {
  159. struct matrix4 temp = *m;
  160. return matrix4_inv(dst, &temp);
  161. }
  162. dstv = (struct vec4 *)dst;
  163. det = matrix4_determinant(m);
  164. if (fabs(det) < 0.0005f)
  165. return false;
  166. for (i = 0; i < 4; i++) {
  167. for (j = 0; j < 4; j++) {
  168. sign = 1 - ((i+j) % 2) * 2;
  169. get_3x3_submatrix(m3x3, m, i, j);
  170. dstv[j].ptr[i] = get_3x3_determinant(m3x3) *
  171. (float)sign / det;
  172. }
  173. }
  174. return true;
  175. }
  176. void matrix4_transpose(struct matrix4 *dst, const struct matrix4 *m)
  177. {
  178. struct matrix4 temp;
  179. /* TODO: Add SSE */
  180. temp.x.x = m->x.x;
  181. temp.x.y = m->y.x;
  182. temp.x.z = m->z.x;
  183. temp.x.w = m->t.x;
  184. temp.y.x = m->x.y;
  185. temp.y.y = m->y.y;
  186. temp.y.z = m->z.y;
  187. temp.y.w = m->t.y;
  188. temp.z.x = m->x.z;
  189. temp.z.y = m->y.z;
  190. temp.z.z = m->z.z;
  191. temp.z.w = m->t.z;
  192. temp.t.x = m->x.w;
  193. temp.t.y = m->y.w;
  194. temp.t.z = m->z.w;
  195. temp.t.w = m->t.w;
  196. matrix4_copy(dst, &temp);
  197. }