# HG changeset patch # User John Tsiombikas # Date 1315459842 -10800 # Node ID c0ae8e66844774c6cf806e060442cd4365a8be06 # Parent fd39c01989350558de1408dc52e7df47ca8cc5a6 added vmath library diff -r fd39c0198935 -r c0ae8e668447 libs/vmath/basis.h --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/libs/vmath/basis.h Thu Sep 08 08:30:42 2011 +0300 @@ -0,0 +1,39 @@ +#ifndef VMATH_BASIS_H_ +#define VMATH_BASIS_H_ + +#include "vector.h" +#include "matrix.h" + +enum { + LEFT_HANDED, + RIGHT_HANDED +}; + +#ifdef __cplusplus +extern "C" { +#endif /* __cplusplus */ + +void basis_matrix(mat4_t res, vec3_t i, vec3_t j, vec3_t k); +void basis_matrix_dir(mat4_t res, vec3_t dir); + +#ifdef __cplusplus +} /* extern "C" */ + +class Basis { +public: + Vector3 i, j, k; + + Basis(); + Basis(const Vector3 &i, const Vector3 &j, const Vector3 &k); + Basis(const Vector3 &dir, bool left_handed = true); + + void rotate(scalar_t x, scalar_t y, scalar_t z); + void rotate(const Vector3 &axis, scalar_t angle); + void rotate(const Matrix4x4 &mat); + void rotate(const Quaternion &quat); + + Matrix3x3 create_rotation_matrix() const; +}; +#endif /* __cplusplus */ + +#endif /* VMATH_BASIS_H_ */ diff -r fd39c0198935 -r c0ae8e668447 libs/vmath/basis_c.c --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/libs/vmath/basis_c.c Thu Sep 08 08:30:42 2011 +0300 @@ -0,0 +1,19 @@ +#include "basis.h" +#include "matrix.h" + +void basis_matrix(mat4_t res, vec3_t i, vec3_t j, vec3_t k) +{ + m4_identity(res); + m4_set_column(res, v4_cons(i.x, i.y, i.z, 1.0), 0); + m4_set_column(res, v4_cons(j.x, j.y, j.z, 1.0), 1); + m4_set_column(res, v4_cons(k.x, k.y, k.z, 1.0), 2); +} + +void basis_matrix_dir(mat4_t res, vec3_t dir) +{ + vec3_t k = v3_normalize(dir); + vec3_t j = {0, 1, 0}; + vec3_t i = v3_cross(j, k); + j = v3_cross(k, i); + basis_matrix(res, i, j, k); +} diff -r fd39c0198935 -r c0ae8e668447 libs/vmath/geom.c --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/libs/vmath/geom.c Thu Sep 08 08:30:42 2011 +0300 @@ -0,0 +1,131 @@ +#include +#include "geom.h" +#include "vector.h" + +plane_t plane_cons(scalar_t nx, scalar_t ny, scalar_t nz, scalar_t d) +{ + plane_t p; + p.norm.x = nx; + p.norm.y = ny; + p.norm.z = nz; + p.d = d; + return p; +} + +plane_t plane_poly(vec3_t v0, vec3_t v1, vec3_t v2) +{ + vec3_t a, b, norm; + + a = v3_sub(v1, v0); + b = v3_sub(v2, v0); + norm = v3_cross(a, b); + norm = v3_normalize(norm); + + return plane_ptnorm(v0, norm); +} + +plane_t plane_ptnorm(vec3_t pt, vec3_t normal) +{ + plane_t plane; + + plane.norm = normal; + plane.d = v3_dot(pt, normal); + + return plane; +} + +plane_t plane_invert(plane_t p) +{ + p.norm = v3_neg(p.norm); + p.d = -p.d; + return p; +} + +scalar_t plane_signed_dist(plane_t plane, vec3_t pt) +{ + vec3_t pp = plane_point(plane); + vec3_t pptopt = v3_sub(pt, pp); + return v3_dot(pptopt, plane.norm); +} + +scalar_t plane_dist(plane_t plane, vec3_t pt) +{ + return fabs(plane_signed_dist(plane, pt)); +} + +vec3_t plane_point(plane_t plane) +{ + return v3_scale(plane.norm, plane.d); +} + +int plane_ray_intersect(ray_t ray, plane_t plane, scalar_t *pos) +{ + vec3_t pt, orig_to_pt; + scalar_t ndotdir; + + pt = plane_point(plane); + ndotdir = v3_dot(plane.norm, ray.dir); + + if(fabs(ndotdir) < 1e-7) { + return 0; + } + + if(pos) { + orig_to_pt = v3_sub(pt, ray.origin); + *pos = v3_dot(plane.norm, orig_to_pt) / ndotdir; + } + return 1; +} + +sphere_t sphere_cons(scalar_t x, scalar_t y, scalar_t z, scalar_t rad) +{ + sphere_t sph; + sph.pos.x = x; + sph.pos.y = y; + sph.pos.z = z; + sph.rad = rad; + return sph; +} + +int sphere_ray_intersect(ray_t ray, sphere_t sph, scalar_t *pos) +{ + scalar_t a, b, c, d, sqrt_d, t1, t2, t; + + a = v3_dot(ray.dir, ray.dir); + b = 2.0 * ray.dir.x * (ray.origin.x - sph.pos.x) + + 2.0 * ray.dir.y * (ray.origin.y - sph.pos.y) + + 2.0 * ray.dir.z * (ray.origin.z - sph.pos.z); + c = v3_dot(sph.pos, sph.pos) + v3_dot(ray.origin, ray.origin) + + 2.0 * v3_dot(v3_neg(sph.pos), ray.origin) - sph.rad * sph.rad; + + d = b * b - 4.0 * a * c; + if(d < 0.0) { + return 0; + } + + sqrt_d = sqrt(d); + t1 = (-b + sqrt_d) / (2.0 * a); + t2 = (-b - sqrt_d) / (2.0 * a); + + if(t1 < 1e-7 || t1 > 1.0) { + t1 = t2; + } + if(t2 < 1e-7 || t2 > 1.0) { + t2 = t1; + } + t = t1 < t2 ? t1 : t2; + + if(t < 1e-7 || t > 1.0) { + return 0; + } + + if(pos) { + *pos = t; + } + return 1; +} + +int sphere_sphere_intersect(sphere_t sph1, sphere_t sph2, scalar_t *pos, scalar_t *rad) +{ + return -1; +} diff -r fd39c0198935 -r c0ae8e668447 libs/vmath/geom.h --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/libs/vmath/geom.h Thu Sep 08 08:30:42 2011 +0300 @@ -0,0 +1,55 @@ +#ifndef GEOM_H_ +#define GEOM_H_ + +#include "vector.h" +#include "ray.h" + +typedef struct { + vec3_t norm; + scalar_t d; +} plane_t; + +typedef struct { + vec3_t pos; + scalar_t rad; +} sphere_t; + +typedef struct { + vec3_t min, max; +} aabox_t; + +#ifdef __cplusplus +extern "C" { +#endif + +/* planes are good... you need planes, yes you do */ +plane_t plane_cons(scalar_t nx, scalar_t ny, scalar_t nz, scalar_t d); +plane_t plane_poly(vec3_t v0, vec3_t v1, vec3_t v2); +plane_t plane_ptnorm(vec3_t pt, vec3_t normal); + +plane_t plane_invert(plane_t p); + +scalar_t plane_signed_dist(plane_t plane, vec3_t pt); +scalar_t plane_dist(plane_t plane, vec3_t pt); +vec3_t plane_point(plane_t plane); + +int plane_ray_intersect(ray_t ray, plane_t plane, scalar_t *pos); + +/* spheres always come in handy */ +sphere_t sphere_cons(scalar_t x, scalar_t y, scalar_t z, scalar_t rad); + +int sphere_ray_intersect(ray_t ray, sphere_t sph, scalar_t *pos); +int sphere_sphere_intersect(sphere_t sph1, sphere_t sph2, scalar_t *pos, scalar_t *rad); + +#ifdef __cplusplus +} + +/* TODO +class Plane : public plane_t { +public: +}; +*/ + +#endif + +#endif /* GEOM_H_ */ diff -r fd39c0198935 -r c0ae8e668447 libs/vmath/matrix.h --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/libs/vmath/matrix.h Thu Sep 08 08:30:42 2011 +0300 @@ -0,0 +1,187 @@ +#ifndef VMATH_MATRIX_H_ +#define VMATH_MATRIX_H_ + +#include +#include "vmath_types.h" + +#ifdef __cplusplus +extern "C" { +#endif /* __cplusplus */ + +/* C matrix 3x3 functions */ +static inline void m3_identity(mat3_t m); +static inline void m3_cons(mat3_t m, + scalar_t m11, scalar_t m12, scalar_t m13, + scalar_t m21, scalar_t m22, scalar_t m23, + scalar_t m31, scalar_t m32, scalar_t m33); +static inline void m3_copy(mat3_t dest, mat3_t src); +void m3_to_m4(mat4_t dest, mat3_t src); + +void m3_print(FILE *fp, mat3_t m); + +/* C matrix 4x4 functions */ +static inline void m4_identity(mat4_t m); +static inline void m4_cons(mat4_t m, + scalar_t m11, scalar_t m12, scalar_t m13, scalar_t m14, + scalar_t m21, scalar_t m22, scalar_t m23, scalar_t m24, + scalar_t m31, scalar_t m32, scalar_t m33, scalar_t m34, + scalar_t m41, scalar_t m42, scalar_t m43, scalar_t m44); +static inline void m4_copy(mat4_t dest, mat4_t src); +void m4_to_m3(mat3_t dest, mat4_t src); + +static inline void m4_mult(mat4_t res, mat4_t m1, mat4_t m2); + +void m4_translate(mat4_t m, scalar_t x, scalar_t y, scalar_t z); +void m4_rotate(mat4_t m, scalar_t x, scalar_t y, scalar_t z); +void m4_rotate_x(mat4_t m, scalar_t angle); +void m4_rotate_y(mat4_t m, scalar_t angle); +void m4_rotate_z(mat4_t m, scalar_t angle); +void m4_rotate_axis(mat4_t m, scalar_t angle, scalar_t x, scalar_t y, scalar_t z); +void m4_rotate_quat(mat4_t m, quat_t q); +void m4_scale(mat4_t m, scalar_t x, scalar_t y, scalar_t z); +static inline void m4_set_column(mat4_t m, vec4_t v, int idx); +static inline void m4_set_row(mat4_t m, vec4_t v, int idx); + +void m4_transpose(mat4_t res, mat4_t m); +scalar_t m4_determinant(mat4_t m); +void m4_adjoint(mat4_t res, mat4_t m); +void m4_inverse(mat4_t res, mat4_t m); + +void m4_print(FILE *fp, mat4_t m); + +#ifdef __cplusplus +} + +/* when included from C++ source files, also define the matrix classes */ +#include + +/** 3x3 matrix */ +class Matrix3x3 { +private: + scalar_t m[3][3]; + +public: + + static Matrix3x3 identity; + + Matrix3x3(); + Matrix3x3( scalar_t m11, scalar_t m12, scalar_t m13, + scalar_t m21, scalar_t m22, scalar_t m23, + scalar_t m31, scalar_t m32, scalar_t m33); + Matrix3x3(const mat3_t cmat); + + Matrix3x3(const Matrix4x4 &mat4x4); + + /* binary operations matrix (op) matrix */ + friend Matrix3x3 operator +(const Matrix3x3 &m1, const Matrix3x3 &m2); + friend Matrix3x3 operator -(const Matrix3x3 &m1, const Matrix3x3 &m2); + friend Matrix3x3 operator *(const Matrix3x3 &m1, const Matrix3x3 &m2); + + friend void operator +=(Matrix3x3 &m1, const Matrix3x3 &m2); + friend void operator -=(Matrix3x3 &m1, const Matrix3x3 &m2); + friend void operator *=(Matrix3x3 &m1, const Matrix3x3 &m2); + + /* binary operations matrix (op) scalar and scalar (op) matrix */ + friend Matrix3x3 operator *(const Matrix3x3 &mat, scalar_t scalar); + friend Matrix3x3 operator *(scalar_t scalar, const Matrix3x3 &mat); + + friend void operator *=(Matrix3x3 &mat, scalar_t scalar); + + inline scalar_t *operator [](int index); + inline const scalar_t *operator [](int index) const; + + inline void reset_identity(); + + void translate(const Vector2 &trans); + void set_translation(const Vector2 &trans); + + void rotate(scalar_t angle); /* 2d rotation */ + void rotate(const Vector3 &euler_angles); /* 3d rotation with euler angles */ + void rotate(const Vector3 &axis, scalar_t angle); /* 3d axis/angle rotation */ + void set_rotation(scalar_t angle); + void set_rotation(const Vector3 &euler_angles); + void set_rotation(const Vector3 &axis, scalar_t angle); + + void scale(const Vector3 &scale_vec); + void set_scaling(const Vector3 &scale_vec); + + void set_column_vector(const Vector3 &vec, unsigned int col_index); + void set_row_vector(const Vector3 &vec, unsigned int row_index); + Vector3 get_column_vector(unsigned int col_index) const; + Vector3 get_row_vector(unsigned int row_index) const; + + void transpose(); + Matrix3x3 transposed() const; + scalar_t determinant() const; + Matrix3x3 inverse() const; + + friend std::ostream &operator <<(std::ostream &out, const Matrix3x3 &mat); +}; + +/** 4x4 matrix */ +class Matrix4x4 { +private: + scalar_t m[4][4]; + +public: + + static Matrix4x4 identity; + + Matrix4x4(); + Matrix4x4( scalar_t m11, scalar_t m12, scalar_t m13, scalar_t m14, + scalar_t m21, scalar_t m22, scalar_t m23, scalar_t m24, + scalar_t m31, scalar_t m32, scalar_t m33, scalar_t m34, + scalar_t m41, scalar_t m42, scalar_t m43, scalar_t m44); + Matrix4x4(const mat4_t cmat); + + Matrix4x4(const Matrix3x3 &mat3x3); + + /* binary operations matrix (op) matrix */ + friend Matrix4x4 operator +(const Matrix4x4 &m1, const Matrix4x4 &m2); + friend Matrix4x4 operator -(const Matrix4x4 &m1, const Matrix4x4 &m2); + friend Matrix4x4 operator *(const Matrix4x4 &m1, const Matrix4x4 &m2); + + friend void operator +=(Matrix4x4 &m1, const Matrix4x4 &m2); + friend void operator -=(Matrix4x4 &m1, const Matrix4x4 &m2); + friend inline void operator *=(Matrix4x4 &m1, const Matrix4x4 &m2); + + /* binary operations matrix (op) scalar and scalar (op) matrix */ + friend Matrix4x4 operator *(const Matrix4x4 &mat, scalar_t scalar); + friend Matrix4x4 operator *(scalar_t scalar, const Matrix4x4 &mat); + + friend void operator *=(Matrix4x4 &mat, scalar_t scalar); + + inline scalar_t *operator [](int index); + inline const scalar_t *operator [](int index) const; + + inline void reset_identity(); + + void translate(const Vector3 &trans); + void set_translation(const Vector3 &trans); + + void rotate(const Vector3 &euler_angles); /* 3d rotation with euler angles */ + void rotate(const Vector3 &axis, scalar_t angle); /* 3d axis/angle rotation */ + void set_rotation(const Vector3 &euler_angles); + void set_rotation(const Vector3 &axis, scalar_t angle); + + void scale(const Vector4 &scale_vec); + void set_scaling(const Vector4 &scale_vec); + + void set_column_vector(const Vector4 &vec, unsigned int col_index); + void set_row_vector(const Vector4 &vec, unsigned int row_index); + Vector4 get_column_vector(unsigned int col_index) const; + Vector4 get_row_vector(unsigned int row_index) const; + + void transpose(); + Matrix4x4 transposed() const; + scalar_t determinant() const; + Matrix4x4 adjoint() const; + Matrix4x4 inverse() const; + + friend std::ostream &operator <<(std::ostream &out, const Matrix4x4 &mat); +}; +#endif /* __cplusplus */ + +#include "matrix.inl" + +#endif /* VMATH_MATRIX_H_ */ diff -r fd39c0198935 -r c0ae8e668447 libs/vmath/matrix_c.c --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/libs/vmath/matrix_c.c Thu Sep 08 08:30:42 2011 +0300 @@ -0,0 +1,246 @@ +#include +#include "matrix.h" +#include "vector.h" +#include "quat.h" + +void m3_to_m4(mat4_t dest, mat3_t src) +{ + int i, j; + + memset(dest, 0, sizeof(mat4_t)); + for(i=0; i<3; i++) { + for(j=0; j<3; j++) { + dest[i][j] = src[i][j]; + } + } + dest[3][3] = 1.0; +} + +void m3_print(FILE *fp, mat3_t m) +{ + int i; + for(i=0; i<3; i++) { + fprintf(fp, "[ %12.5f %12.5f %12.5f ]\n", (float)m[i][0], (float)m[i][1], (float)m[i][2]); + } +} + +/* C matrix 4x4 functions */ +void m4_to_m3(mat3_t dest, mat4_t src) +{ + int i, j; + for(i=0; i<3; i++) { + for(j=0; j<3; j++) { + dest[i][j] = src[i][j]; + } + } +} + +void m4_translate(mat4_t m, scalar_t x, scalar_t y, scalar_t z) +{ + mat4_t tm; + m4_identity(tm); + tm[0][3] = x; + tm[1][3] = y; + tm[2][3] = z; + m4_mult(m, m, tm); +} + +void m4_rotate(mat4_t m, scalar_t x, scalar_t y, scalar_t z) +{ + m4_rotate_x(m, x); + m4_rotate_y(m, y); + m4_rotate_z(m, z); +} + +void m4_rotate_x(mat4_t m, scalar_t angle) +{ + mat4_t rm; + m4_identity(rm); + rm[1][1] = cos(angle); rm[1][2] = -sin(angle); + rm[2][1] = sin(angle); rm[2][2] = cos(angle); + m4_mult(m, m, rm); +} + +void m4_rotate_y(mat4_t m, scalar_t angle) +{ + mat4_t rm; + m4_identity(rm); + rm[0][0] = cos(angle); rm[0][2] = sin(angle); + rm[2][0] = -sin(angle); rm[2][2] = cos(angle); + m4_mult(m, m, rm); +} + +void m4_rotate_z(mat4_t m, scalar_t angle) +{ + mat4_t rm; + m4_identity(rm); + rm[0][0] = cos(angle); rm[0][1] = -sin(angle); + rm[1][0] = sin(angle); rm[1][1] = cos(angle); + m4_mult(m, m, rm); +} + +void m4_rotate_axis(mat4_t m, scalar_t angle, scalar_t x, scalar_t y, scalar_t z) +{ + mat4_t xform; + scalar_t sina = sin(angle); + scalar_t cosa = cos(angle); + scalar_t one_minus_cosa = 1.0 - cosa; + scalar_t nxsq = x * x; + scalar_t nysq = y * y; + scalar_t nzsq = z * z; + + m4_identity(xform); + xform[0][0] = nxsq + (1.0 - nxsq) * cosa; + xform[0][1] = x * y * one_minus_cosa - z * sina; + xform[0][2] = x * z * one_minus_cosa + y * sina; + xform[1][0] = x * y * one_minus_cosa + z * sina; + xform[1][1] = nysq + (1.0 - nysq) * cosa; + xform[1][2] = y * z * one_minus_cosa - x * sina; + xform[2][0] = x * z * one_minus_cosa - y * sina; + xform[2][1] = y * z * one_minus_cosa + x * sina; + xform[2][2] = nzsq + (1.0 - nzsq) * cosa; + + m4_mult(m, m, xform); +} + +void m4_rotate_quat(mat4_t m, quat_t q) +{ + mat4_t rm; + quat_to_mat4(rm, q); + m4_mult(m, m, rm); +} + +void m4_scale(mat4_t m, scalar_t x, scalar_t y, scalar_t z) +{ + mat4_t sm; + m4_identity(sm); + sm[0][0] = x; + sm[1][1] = y; + sm[2][2] = z; + m4_mult(m, m, sm); +} + +void m4_transpose(mat4_t res, mat4_t m) +{ + int i, j; + mat4_t tmp; + m4_copy(tmp, m); + + for(i=0; i<4; i++) { + for(j=0; j<4; j++) { + res[i][j] = tmp[j][i]; + } + } +} + +scalar_t m4_determinant(mat4_t m) +{ + scalar_t det11 = (m[1][1] * (m[2][2] * m[3][3] - m[3][2] * m[2][3])) - + (m[1][2] * (m[2][1] * m[3][3] - m[3][1] * m[2][3])) + + (m[1][3] * (m[2][1] * m[3][2] - m[3][1] * m[2][2])); + + scalar_t det12 = (m[1][0] * (m[2][2] * m[3][3] - m[3][2] * m[2][3])) - + (m[1][2] * (m[2][0] * m[3][3] - m[3][0] * m[2][3])) + + (m[1][3] * (m[2][0] * m[3][2] - m[3][0] * m[2][2])); + + scalar_t det13 = (m[1][0] * (m[2][1] * m[3][3] - m[3][1] * m[2][3])) - + (m[1][1] * (m[2][0] * m[3][3] - m[3][0] * m[2][3])) + + (m[1][3] * (m[2][0] * m[3][1] - m[3][0] * m[2][1])); + + scalar_t det14 = (m[1][0] * (m[2][1] * m[3][2] - m[3][1] * m[2][2])) - + (m[1][1] * (m[2][0] * m[3][2] - m[3][0] * m[2][2])) + + (m[1][2] * (m[2][0] * m[3][1] - m[3][0] * m[2][1])); + + return m[0][0] * det11 - m[0][1] * det12 + m[0][2] * det13 - m[0][3] * det14; +} + +void m4_adjoint(mat4_t res, mat4_t m) +{ + int i, j; + mat4_t coef; + + coef[0][0] = (m[1][1] * (m[2][2] * m[3][3] - m[3][2] * m[2][3])) - + (m[1][2] * (m[2][1] * m[3][3] - m[3][1] * m[2][3])) + + (m[1][3] * (m[2][1] * m[3][2] - m[3][1] * m[2][2])); + coef[0][1] = (m[1][0] * (m[2][2] * m[3][3] - m[3][2] * m[2][3])) - + (m[1][2] * (m[2][0] * m[3][3] - m[3][0] * m[2][3])) + + (m[1][3] * (m[2][0] * m[3][2] - m[3][0] * m[2][2])); + coef[0][2] = (m[1][0] * (m[2][1] * m[3][3] - m[3][1] * m[2][3])) - + (m[1][1] * (m[2][0] * m[3][3] - m[3][0] * m[2][3])) + + (m[1][3] * (m[2][0] * m[3][1] - m[3][0] * m[2][1])); + coef[0][3] = (m[1][0] * (m[2][1] * m[3][2] - m[3][1] * m[2][2])) - + (m[1][1] * (m[2][0] * m[3][2] - m[3][0] * m[2][2])) + + (m[1][2] * (m[2][0] * m[3][1] - m[3][0] * m[2][1])); + + coef[1][0] = (m[0][1] * (m[2][2] * m[3][3] - m[3][2] * m[2][3])) - + (m[0][2] * (m[2][1] * m[3][3] - m[3][1] * m[2][3])) + + (m[0][3] * (m[2][1] * m[3][2] - m[3][1] * m[2][2])); + coef[1][1] = (m[0][0] * (m[2][2] * m[3][3] - m[3][2] * m[2][3])) - + (m[0][2] * (m[2][0] * m[3][3] - m[3][0] * m[2][3])) + + (m[0][3] * (m[2][0] * m[3][2] - m[3][0] * m[2][2])); + coef[1][2] = (m[0][0] * (m[2][1] * m[3][3] - m[3][1] * m[2][3])) - + (m[0][1] * (m[2][0] * m[3][3] - m[3][0] * m[2][3])) + + (m[0][3] * (m[2][0] * m[3][1] - m[3][0] * m[2][1])); + coef[1][3] = (m[0][0] * (m[2][1] * m[3][2] - m[3][1] * m[2][2])) - + (m[0][1] * (m[2][0] * m[3][2] - m[3][0] * m[2][2])) + + (m[0][2] * (m[2][0] * m[3][1] - m[3][0] * m[2][1])); + + coef[2][0] = (m[0][1] * (m[1][2] * m[3][3] - m[3][2] * m[1][3])) - + (m[0][2] * (m[1][1] * m[3][3] - m[3][1] * m[1][3])) + + (m[0][3] * (m[1][1] * m[3][2] - m[3][1] * m[1][2])); + coef[2][1] = (m[0][0] * (m[1][2] * m[3][3] - m[3][2] * m[1][3])) - + (m[0][2] * (m[1][0] * m[3][3] - m[3][0] * m[1][3])) + + (m[0][3] * (m[1][0] * m[3][2] - m[3][0] * m[1][2])); + coef[2][2] = (m[0][0] * (m[1][1] * m[3][3] - m[3][1] * m[1][3])) - + (m[0][1] * (m[1][0] * m[3][3] - m[3][0] * m[1][3])) + + (m[0][3] * (m[1][0] * m[3][1] - m[3][0] * m[1][1])); + coef[2][3] = (m[0][0] * (m[1][1] * m[3][2] - m[3][1] * m[1][2])) - + (m[0][1] * (m[1][0] * m[3][2] - m[3][0] * m[1][2])) + + (m[0][2] * (m[1][0] * m[3][1] - m[3][0] * m[1][1])); + + coef[3][0] = (m[0][1] * (m[1][2] * m[2][3] - m[2][2] * m[1][3])) - + (m[0][2] * (m[1][1] * m[2][3] - m[2][1] * m[1][3])) + + (m[0][3] * (m[1][1] * m[2][2] - m[2][1] * m[1][2])); + coef[3][1] = (m[0][0] * (m[1][2] * m[2][3] - m[2][2] * m[1][3])) - + (m[0][2] * (m[1][0] * m[2][3] - m[2][0] * m[1][3])) + + (m[0][3] * (m[1][0] * m[2][2] - m[2][0] * m[1][2])); + coef[3][2] = (m[0][0] * (m[1][1] * m[2][3] - m[2][1] * m[1][3])) - + (m[0][1] * (m[1][0] * m[2][3] - m[2][0] * m[1][3])) + + (m[0][3] * (m[1][0] * m[2][1] - m[2][0] * m[1][1])); + coef[3][3] = (m[0][0] * (m[1][1] * m[2][2] - m[2][1] * m[1][2])) - + (m[0][1] * (m[1][0] * m[2][2] - m[2][0] * m[1][2])) + + (m[0][2] * (m[1][0] * m[2][1] - m[2][0] * m[1][1])); + + m4_transpose(res, coef); + + for(i=0; i<4; i++) { + for(j=0; j<4; j++) { + res[i][j] = j % 2 ? -res[i][j] : res[i][j]; + if(i % 2) res[i][j] = -res[i][j]; + } + } +} + +void m4_inverse(mat4_t res, mat4_t m) +{ + int i, j; + mat4_t adj; + scalar_t det; + + m4_adjoint(adj, m); + det = m4_determinant(m); + + for(i=0; i<4; i++) { + for(j=0; j<4; j++) { + res[i][j] = adj[i][j] / det; + } + } +} + +void m4_print(FILE *fp, mat4_t m) +{ + int i; + for(i=0; i<4; i++) { + fprintf(fp, "[ %12.5f %12.5f %12.5f %12.5f ]\n", (float)m[i][0], (float)m[i][1], (float)m[i][2], (float)m[i][3]); + } +} diff -r fd39c0198935 -r c0ae8e668447 libs/vmath/quat.h --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/libs/vmath/quat.h Thu Sep 08 08:30:42 2011 +0300 @@ -0,0 +1,101 @@ +#ifndef VMATH_QUATERNION_H_ +#define VMATH_QUATERNION_H_ + +#include +#include "vmath_types.h" +#include "vector.h" + +#ifdef __cplusplus +extern "C" { +#endif /* __cplusplus */ + +#define quat_cons(s, x, y, z) v4_cons(x, y, z, s) +#define quat_vec(q) v3_cons((q).x, (q).y, (q).z) +#define quat_s(q) ((q).w) +#define quat_identity() quat_cons(1.0, 0.0, 0.0, 0.0) +void quat_print(FILE *fp, quat_t q); + +#define quat_add v4_add +#define quat_sub v4_sub +#define quat_neg v4_neg + +static inline quat_t quat_mul(quat_t q1, quat_t q2); + +static inline quat_t quat_conjugate(quat_t q); + +#define quat_length v4_length +#define quat_length_sq v4_length_sq + +#define quat_normalize v4_normalize +static inline quat_t quat_inverse(quat_t q); + +quat_t quat_rotate(quat_t q, scalar_t angle, scalar_t x, scalar_t y, scalar_t z); +quat_t quat_rotate_quat(quat_t q, quat_t rotq); + +static inline void quat_to_mat3(mat3_t res, quat_t q); +static inline void quat_to_mat4(mat4_t res, quat_t q); + +#define quat_lerp quat_slerp +quat_t quat_slerp(quat_t q1, quat_t q2, scalar_t t); + + +#ifdef __cplusplus +} /* extern "C" */ + +#include + +/* Quaternion */ +class Quaternion { +public: + scalar_t s; + Vector3 v; + + Quaternion(); + Quaternion(scalar_t s, const Vector3 &v); + Quaternion(scalar_t s, scalar_t x, scalar_t y, scalar_t z); + Quaternion(const Vector3 &axis, scalar_t angle); + Quaternion(const quat_t &quat); + + Quaternion operator +(const Quaternion &quat) const; + Quaternion operator -(const Quaternion &quat) const; + Quaternion operator -() const; + Quaternion operator *(const Quaternion &quat) const; + + void operator +=(const Quaternion &quat); + void operator -=(const Quaternion &quat); + void operator *=(const Quaternion &quat); + + void reset_identity(); + + Quaternion conjugate() const; + + scalar_t length() const; + scalar_t length_sq() const; + + void normalize(); + Quaternion normalized() const; + + Quaternion inverse() const; + + void set_rotation(const Vector3 &axis, scalar_t angle); + void rotate(const Vector3 &axis, scalar_t angle); + /* note: this is a totally different operation from the above + * this treats the quaternion as signifying direction and rotates + * it by a rotation quaternion by rot * q * rot' + */ + void rotate(const Quaternion &q); + + Matrix3x3 get_rotation_matrix() const; + + friend Quaternion slerp(const Quaternion &q1, const Quaternion &q2, scalar_t t); + + friend std::ostream &operator <<(std::ostream &out, const Quaternion &q); +}; + +Quaternion slerp(const Quaternion &q1, const Quaternion &q2, scalar_t t); +inline Quaternion lerp(const Quaternion &q1, const Quaternion &q2, scalar_t t); +#endif /* __cplusplus */ + +#include "quat.inl" + +#endif /* VMATH_QUATERNION_H_ */ diff -r fd39c0198935 -r c0ae8e668447 libs/vmath/quat_c.c --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/libs/vmath/quat_c.c Thu Sep 08 08:30:42 2011 +0300 @@ -0,0 +1,42 @@ +#include +#include +#include "quat.h" + +void quat_print(FILE *fp, quat_t q) +{ + fprintf(fp, "([ %.4f %.4f %.4f ] %.4f)", q.x, q.y, q.z, q.w); +} + +quat_t quat_rotate(quat_t q, scalar_t angle, scalar_t x, scalar_t y, scalar_t z) +{ + quat_t rq; + scalar_t half_angle = angle * 0.5; + scalar_t sin_half = sin(half_angle); + + rq.w = cos(half_angle); + rq.x = x * sin_half; + rq.y = y * sin_half; + rq.z = z * sin_half; + + return quat_mul(q, rq); +} + +quat_t quat_rotate_quat(quat_t q, quat_t rotq) +{ + return quat_mul(quat_mul(rotq, q), quat_conjugate(rotq)); +} + +quat_t quat_slerp(quat_t q1, quat_t q2, scalar_t t) +{ + quat_t res; + scalar_t angle = acos(q1.w * q2.w + q1.x * q2.x + q1.y * q2.y + q1.z * q2.z); + scalar_t a = sin((1.0f - t) * angle); + scalar_t b = sin(t * angle); + scalar_t c = sin(angle); + + res.x = (q1.x * a + q2.x * b) / c; + res.y = (q1.y * a + q2.y * b) / c; + res.z = (q1.z * a + q2.z * b) / c; + res.w = (q1.w * a + q2.w * b) / c; + return quat_normalize(res); +} diff -r fd39c0198935 -r c0ae8e668447 libs/vmath/ray.h --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/libs/vmath/ray.h Thu Sep 08 08:30:42 2011 +0300 @@ -0,0 +1,55 @@ +#ifndef VMATH_RAY_H_ +#define VMATH_RAY_H_ + +#include "matrix.h" +#include "vector.h" + +typedef struct { + vec3_t origin, dir; +} ray_t; + +#ifdef __cplusplus +extern "C" { +#endif /* __cplusplus */ + +static inline ray_t ray_cons(vec3_t origin, vec3_t dir); +ray_t ray_transform(ray_t r, mat4_t m); + +#ifdef __cplusplus +} /* __cplusplus */ + +#include + +class Ray { +private: + std::stack ior_stack; + +public: + /* enviornmental index of refraction, normally 1.0 */ + static scalar_t env_ior; + + Vector3 origin, dir; + scalar_t energy; + int iter; + scalar_t ior; + long time; + + Ray(); + Ray(const Vector3 &origin, const Vector3 &dir); + + void transform(const Matrix4x4 &xform); + Ray transformed(const Matrix4x4 &xform) const; + + void enter(scalar_t new_ior); + void leave(); + + scalar_t calc_ior(bool entering, scalar_t mat_ior = 1.0) const; +}; + +inline Ray reflect_ray(const Ray &inray, const Vector3 &norm); +inline Ray refract_ray(const Ray &inray, const Vector3 &norm, scalar_t ior, bool entering, scalar_t ray_mag = -1.0); +#endif /* __cplusplus */ + +#include "ray.inl" + +#endif /* VMATH_RAY_H_ */ diff -r fd39c0198935 -r c0ae8e668447 libs/vmath/ray_c.c --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/libs/vmath/ray_c.c Thu Sep 08 08:30:42 2011 +0300 @@ -0,0 +1,18 @@ +#include "ray.h" +#include "vector.h" + +ray_t ray_transform(ray_t r, mat4_t xform) +{ + mat4_t upper; + vec3_t dir; + + m4_copy(upper, xform); + upper[0][3] = upper[1][3] = upper[2][3] = upper[3][0] = upper[3][1] = upper[3][2] = 0.0; + upper[3][3] = 1.0; + + dir = v3_sub(r.dir, r.origin); + dir = v3_transform(dir, upper); + r.origin = v3_transform(r.origin, xform); + r.dir = v3_add(dir, r.origin); + return r; +} diff -r fd39c0198935 -r c0ae8e668447 libs/vmath/sphvec.h --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/libs/vmath/sphvec.h Thu Sep 08 08:30:42 2011 +0300 @@ -0,0 +1,18 @@ +#ifndef VMATH_SPHVEC_H_ +#define VMATH_SPHVEC_H_ + +#include "vmath_types.h" + +#ifdef __cplusplus +/* Vector in spherical coordinates */ +class SphVector { +public: + scalar_t theta, phi, r; + + SphVector(scalar_t theta = 0.0, scalar_t phi = 0.0, scalar_t r = 1.0); + SphVector(const Vector3 &cvec); + SphVector &operator =(const Vector3 &cvec); +}; +#endif /* __cplusplus */ + +#endif /* VMATH_SPHVEC_H_ */ diff -r fd39c0198935 -r c0ae8e668447 libs/vmath/vector.h --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/libs/vmath/vector.h Thu Sep 08 08:30:42 2011 +0300 @@ -0,0 +1,274 @@ +#ifndef VMATH_VECTOR_H_ +#define VMATH_VECTOR_H_ + +#include +#include "vmath_types.h" + +#ifdef __cplusplus +extern "C" { +#endif /* __cplusplus */ + +/* C 2D vector functions */ +static inline vec2_t v2_cons(scalar_t x, scalar_t y); +static inline void v2_print(FILE *fp, vec2_t v); + +static inline vec2_t v2_add(vec2_t v1, vec2_t v2); +static inline vec2_t v2_sub(vec2_t v1, vec2_t v2); +static inline vec2_t v2_scale(vec2_t v, scalar_t s); +static inline scalar_t v2_dot(vec2_t v1, vec2_t v2); +static inline scalar_t v2_length(vec2_t v); +static inline scalar_t v2_length_sq(vec2_t v); +static inline vec2_t v2_normalize(vec2_t v); + +static inline vec2_t v2_lerp(vec2_t v1, vec2_t v2, scalar_t t); + +/* C 3D vector functions */ +static inline vec3_t v3_cons(scalar_t x, scalar_t y, scalar_t z); +static inline void v3_print(FILE *fp, vec3_t v); + +static inline vec3_t v3_add(vec3_t v1, vec3_t v2); +static inline vec3_t v3_sub(vec3_t v1, vec3_t v2); +static inline vec3_t v3_neg(vec3_t v); +static inline vec3_t v3_mul(vec3_t v1, vec3_t v2); +static inline vec3_t v3_scale(vec3_t v1, scalar_t s); +static inline scalar_t v3_dot(vec3_t v1, vec3_t v2); +static inline vec3_t v3_cross(vec3_t v1, vec3_t v2); +static inline scalar_t v3_length(vec3_t v); +static inline scalar_t v3_length_sq(vec3_t v); +static inline vec3_t v3_normalize(vec3_t v); +static inline vec3_t v3_transform(vec3_t v, mat4_t m); + +static inline vec3_t v3_rotate(vec3_t v, scalar_t x, scalar_t y, scalar_t z); +static inline vec3_t v3_rotate_axis(vec3_t v, scalar_t angle, scalar_t x, scalar_t y, scalar_t z); +static inline vec3_t v3_rotate_quat(vec3_t v, quat_t q); + +static inline vec3_t v3_reflect(vec3_t v, vec3_t n); + +static inline vec3_t v3_lerp(vec3_t v1, vec3_t v2, scalar_t t); + +/* C 4D vector functions */ +static inline vec4_t v4_cons(scalar_t x, scalar_t y, scalar_t z, scalar_t w); +static inline void v4_print(FILE *fp, vec4_t v); + +static inline vec4_t v4_add(vec4_t v1, vec4_t v2); +static inline vec4_t v4_sub(vec4_t v1, vec4_t v2); +static inline vec4_t v4_neg(vec4_t v); +static inline vec4_t v4_mul(vec4_t v1, vec4_t v2); +static inline vec4_t v4_scale(vec4_t v, scalar_t s); +static inline scalar_t v4_dot(vec4_t v1, vec4_t v2); +static inline scalar_t v4_length(vec4_t v); +static inline scalar_t v4_length_sq(vec4_t v); +static inline vec4_t v4_normalize(vec4_t v); +static inline vec4_t v4_transform(vec4_t v, mat4_t m); + +#ifdef __cplusplus +} /* extern "C" */ + +/* when included from C++ source files, also define the vector classes */ +#include + +/** 2D Vector */ +class Vector2 { +public: + scalar_t x, y; + + explicit Vector2(scalar_t x = 0.0, scalar_t y = 0.0); + Vector2(const vec2_t &vec); + Vector2(const Vector3 &vec); + Vector2(const Vector4 &vec); + + inline scalar_t &operator [](int elem); + inline const scalar_t &operator [](int elem) const; + + inline scalar_t length() const; + inline scalar_t length_sq() const; + void normalize(); + Vector2 normalized() const; + + void transform(const Matrix3x3 &mat); + Vector2 transformed(const Matrix3x3 &mat) const; + + void rotate(scalar_t angle); + Vector2 rotated(scalar_t angle) const; + + Vector2 reflection(const Vector2 &normal) const; + Vector2 refraction(const Vector2 &normal, scalar_t src_ior, scalar_t dst_ior) const; +}; + +/* unary operations */ +inline Vector2 operator -(const Vector2 &vec); + +/* binary vector (op) vector operations */ +inline scalar_t dot_product(const Vector2 &v1, const Vector2 &v2); + +inline Vector2 operator +(const Vector2 &v1, const Vector2 &v2); +inline Vector2 operator -(const Vector2 &v1, const Vector2 &v2); +inline Vector2 operator *(const Vector2 &v1, const Vector2 &v2); +inline Vector2 operator /(const Vector2 &v1, const Vector2 &v2); +inline bool operator ==(const Vector2 &v1, const Vector2 &v2); + +inline void operator +=(Vector2 &v1, const Vector2 &v2); +inline void operator -=(Vector2 &v1, const Vector2 &v2); +inline void operator *=(Vector2 &v1, const Vector2 &v2); +inline void operator /=(Vector2 &v1, const Vector2 &v2); + +/* binary vector (op) scalar and scalar (op) vector operations */ +inline Vector2 operator +(const Vector2 &vec, scalar_t scalar); +inline Vector2 operator +(scalar_t scalar, const Vector2 &vec); +inline Vector2 operator -(const Vector2 &vec, scalar_t scalar); +inline Vector2 operator *(const Vector2 &vec, scalar_t scalar); +inline Vector2 operator *(scalar_t scalar, const Vector2 &vec); +inline Vector2 operator /(const Vector2 &vec, scalar_t scalar); + +inline void operator +=(Vector2 &vec, scalar_t scalar); +inline void operator -=(Vector2 &vec, scalar_t scalar); +inline void operator *=(Vector2 &vec, scalar_t scalar); +inline void operator /=(Vector2 &vec, scalar_t scalar); + +std::ostream &operator <<(std::ostream &out, const Vector2 &vec); + +inline Vector2 lerp(const Vector2 &a, const Vector2 &b, scalar_t t); +inline Vector2 catmull_rom_spline(const Vector2 &v0, const Vector2 &v1, + const Vector2 &v2, const Vector2 &v3, scalar_t t); + +/* 3D Vector */ +class Vector3 { +public: + scalar_t x, y, z; + + explicit Vector3(scalar_t x = 0.0, scalar_t y = 0.0, scalar_t z = 0.0); + Vector3(const vec3_t &vec); + Vector3(const Vector2 &vec); + Vector3(const Vector4 &vec); + Vector3(const SphVector &sph); + + Vector3 &operator =(const SphVector &sph); + + inline scalar_t &operator [](int elem); + inline const scalar_t &operator [](int elem) const; + + inline scalar_t length() const; + inline scalar_t length_sq() const; + void normalize(); + Vector3 normalized() const; + + void transform(const Matrix3x3 &mat); + Vector3 transformed(const Matrix3x3 &mat) const; + void transform(const Matrix4x4 &mat); + Vector3 transformed(const Matrix4x4 &mat) const; + void transform(const Quaternion &quat); + Vector3 transformed(const Quaternion &quat) const; + + void rotate(const Vector3 &euler); + Vector3 rotated(const Vector3 &euler) const; + + Vector3 reflection(const Vector3 &normal) const; + Vector3 refraction(const Vector3 &normal, scalar_t src_ior, scalar_t dst_ior) const; + Vector3 refraction(const Vector3 &normal, scalar_t ior) const; +}; + +/* unary operations */ +inline Vector3 operator -(const Vector3 &vec); + +/* binary vector (op) vector operations */ +inline scalar_t dot_product(const Vector3 &v1, const Vector3 &v2); +inline Vector3 cross_product(const Vector3 &v1, const Vector3 &v2); + +inline Vector3 operator +(const Vector3 &v1, const Vector3 &v2); +inline Vector3 operator -(const Vector3 &v1, const Vector3 &v2); +inline Vector3 operator *(const Vector3 &v1, const Vector3 &v2); +inline Vector3 operator /(const Vector3 &v1, const Vector3 &v2); +inline bool operator ==(const Vector3 &v1, const Vector3 &v2); + +inline void operator +=(Vector3 &v1, const Vector3 &v2); +inline void operator -=(Vector3 &v1, const Vector3 &v2); +inline void operator *=(Vector3 &v1, const Vector3 &v2); +inline void operator /=(Vector3 &v1, const Vector3 &v2); + +/* binary vector (op) scalar and scalar (op) vector operations */ +inline Vector3 operator +(const Vector3 &vec, scalar_t scalar); +inline Vector3 operator +(scalar_t scalar, const Vector3 &vec); +inline Vector3 operator -(const Vector3 &vec, scalar_t scalar); +inline Vector3 operator *(const Vector3 &vec, scalar_t scalar); +inline Vector3 operator *(scalar_t scalar, const Vector3 &vec); +inline Vector3 operator /(const Vector3 &vec, scalar_t scalar); + +inline void operator +=(Vector3 &vec, scalar_t scalar); +inline void operator -=(Vector3 &vec, scalar_t scalar); +inline void operator *=(Vector3 &vec, scalar_t scalar); +inline void operator /=(Vector3 &vec, scalar_t scalar); + +std::ostream &operator <<(std::ostream &out, const Vector3 &vec); + +inline Vector3 lerp(const Vector3 &a, const Vector3 &b, scalar_t t); +inline Vector3 catmull_rom_spline(const Vector3 &v0, const Vector3 &v1, + const Vector3 &v2, const Vector3 &v3, scalar_t t); + +/* 4D Vector */ +class Vector4 { +public: + scalar_t x, y, z, w; + + explicit Vector4(scalar_t x = 0.0, scalar_t y = 0.0, scalar_t z = 0.0, scalar_t w = 0.0); + Vector4(const vec4_t &vec); + Vector4(const Vector2 &vec); + Vector4(const Vector3 &vec); + + inline scalar_t &operator [](int elem); + inline const scalar_t &operator [](int elem) const; + + inline scalar_t length() const; + inline scalar_t length_sq() const; + void normalize(); + Vector4 normalized() const; + + void transform(const Matrix4x4 &mat); + Vector4 transformed(const Matrix4x4 &mat) const; + + Vector4 reflection(const Vector4 &normal) const; + Vector4 refraction(const Vector4 &normal, scalar_t src_ior, scalar_t dst_ior) const; +}; + + +/* unary operations */ +inline Vector4 operator -(const Vector4 &vec); + +/* binary vector (op) vector operations */ +inline scalar_t dot_product(const Vector4 &v1, const Vector4 &v2); +inline Vector4 cross_product(const Vector4 &v1, const Vector4 &v2, const Vector4 &v3); + +inline Vector4 operator +(const Vector4 &v1, const Vector4 &v2); +inline Vector4 operator -(const Vector4 &v1, const Vector4 &v2); +inline Vector4 operator *(const Vector4 &v1, const Vector4 &v2); +inline Vector4 operator /(const Vector4 &v1, const Vector4 &v2); +inline bool operator ==(const Vector4 &v1, const Vector4 &v2); + +inline void operator +=(Vector4 &v1, const Vector4 &v2); +inline void operator -=(Vector4 &v1, const Vector4 &v2); +inline void operator *=(Vector4 &v1, const Vector4 &v2); +inline void operator /=(Vector4 &v1, const Vector4 &v2); + +/* binary vector (op) scalar and scalar (op) vector operations */ +inline Vector4 operator +(const Vector4 &vec, scalar_t scalar); +inline Vector4 operator +(scalar_t scalar, const Vector4 &vec); +inline Vector4 operator -(const Vector4 &vec, scalar_t scalar); +inline Vector4 operator *(const Vector4 &vec, scalar_t scalar); +inline Vector4 operator *(scalar_t scalar, const Vector4 &vec); +inline Vector4 operator /(const Vector4 &vec, scalar_t scalar); + +inline void operator +=(Vector4 &vec, scalar_t scalar); +inline void operator -=(Vector4 &vec, scalar_t scalar); +inline void operator *=(Vector4 &vec, scalar_t scalar); +inline void operator /=(Vector4 &vec, scalar_t scalar); + +std::ostream &operator <<(std::ostream &out, const Vector4 &vec); + +inline Vector4 lerp(const Vector4 &v0, const Vector4 &v1, scalar_t t); +inline Vector4 catmull_rom_spline(const Vector4 &v0, const Vector4 &v1, + const Vector4 &v2, const Vector4 &v3, scalar_t t); + +#endif /* __cplusplus */ + +#include "vector.inl" + +#endif /* VMATH_VECTOR_H_ */ diff -r fd39c0198935 -r c0ae8e668447 libs/vmath/vmath.c --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/libs/vmath/vmath.c Thu Sep 08 08:30:42 2011 +0300 @@ -0,0 +1,317 @@ +#include +#include +#include "vmath.h" + +/** Numerical calculation of integrals using simpson's rule */ +scalar_t integral(scalar_t (*f)(scalar_t), scalar_t low, scalar_t high, int samples) +{ + int i; + scalar_t h = (high - low) / (scalar_t)samples; + scalar_t sum = 0.0; + + for(i=0; i +#include "vmath_types.h" + +#ifndef M_PI +#define M_PI PI +#endif + +#ifndef M_E +#define M_E 2.718281828459045 +#endif + +#define PI 3.141592653589793 +#define HALF_PI 1.570796326794897 +#define QUARTER_PI 0.785398163397448 +#define TWO_PI 6.283185307179586 + + +#define RAD_TO_DEG(a) ((((scalar_t)a) * 360.0) / TWO_PI) +#define DEG_TO_RAD(a) (((scalar_t)a) * (PI / 180.0)) + +#define SQ(x) ((x) * (x)) + +#define MIN(a, b) ((a) < (b) ? (a) : (b)) +#define MAX(a, b) ((a) > (b) ? (a) : (b)) + +#ifndef __GNUC__ +#define round(x) ((x) >= 0 ? (x) + 0.5 : (x) - 0.5) +#endif + +#ifdef __cplusplus +extern "C" { +#endif /* __cplusplus */ + +static inline scalar_t frand(scalar_t range); +static inline vec3_t sphrand(scalar_t rad); + +scalar_t integral(scalar_t (*f)(scalar_t), scalar_t low, scalar_t high, int samples); +scalar_t gaussian(scalar_t x, scalar_t mean, scalar_t sdev); + +static inline scalar_t lerp(scalar_t a, scalar_t b, scalar_t t); + +scalar_t bspline(scalar_t a, scalar_t b, scalar_t c, scalar_t d, scalar_t t); +scalar_t spline(scalar_t a, scalar_t b, scalar_t c, scalar_t d, scalar_t t); +scalar_t bezier(scalar_t a, scalar_t b, scalar_t c, scalar_t d, scalar_t t); + +scalar_t noise1(scalar_t x); +scalar_t noise2(scalar_t x, scalar_t y); +scalar_t noise3(scalar_t x, scalar_t y, scalar_t z); + +scalar_t fbm1(scalar_t x, int octaves); +scalar_t fbm2(scalar_t x, scalar_t y, int octaves); +scalar_t fbm3(scalar_t x, scalar_t y, scalar_t z, int octaves); + +scalar_t turbulence1(scalar_t x, int octaves); +scalar_t turbulence2(scalar_t x, scalar_t y, int octaves); +scalar_t turbulence3(scalar_t x, scalar_t y, scalar_t z, int octaves); + +#ifdef __cplusplus +} +#endif /* __cplusplus */ + +#include "vmath.inl" + +#include "vector.h" +#include "matrix.h" +#include "quat.h" +#include "sphvec.h" +#include "ray.h" +#include "geom.h" + +#endif /* VMATH_H_ */ diff -r fd39c0198935 -r c0ae8e668447 libs/vmath/vmath_config.h --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/libs/vmath/vmath_config.h Thu Sep 08 08:30:42 2011 +0300 @@ -0,0 +1,19 @@ +#ifndef VMATH_CONFIG_H_ +#define VMATH_CONFIG_H_ + +#if (__STDC_VERSION__ < 199999) +#if defined(__GNUC__) || defined(_MSC_VER) +#define inline __inline +#else +#define inline + +#ifdef VECTOR_H_ +#warning "compiling vector operations without inline, performance might suffer" +#endif /* VECTOR_H_ */ + +#endif /* gcc/msvc */ +#endif /* not C99 */ + +#define SINGLE_PRECISION_MATH + +#endif /* VMATH_CONFIG_H_ */ diff -r fd39c0198935 -r c0ae8e668447 libs/vmath/vmath_types.h --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/libs/vmath/vmath_types.h Thu Sep 08 08:30:42 2011 +0300 @@ -0,0 +1,40 @@ +#ifndef VMATH_TYPES_H_ +#define VMATH_TYPES_H_ + +#include "vmath_config.h" + +#define SMALL_NUMBER 1.e-4 +#define XSMALL_NUMBER 1.e-8 +#define ERROR_MARGIN 1.e-6 + + +#ifdef SINGLE_PRECISION_MATH +typedef float scalar_t; +#else +typedef double scalar_t; +#endif /* floating point precision */ + +/* vectors */ +typedef struct { scalar_t x, y; } vec2_t; +typedef struct { scalar_t x, y, z; } vec3_t; +typedef struct { scalar_t x, y, z, w; } vec4_t; + +/* quaternions */ +typedef vec4_t quat_t; + +/* matrices */ +typedef scalar_t mat3_t[3][3]; +typedef scalar_t mat4_t[4][4]; + + +#ifdef __cplusplus +class Vector2; +class Vector3; +class Vector4; +class Quaternion; +class Matrix3x3; +class Matrix4x4; +class SphVector; +#endif /* __cplusplus */ + +#endif /* VMATH_TYPES_H_ */