its ot'
[vg.git] / src / vg / vg_m.h
index c71f3abd5af081f077c6547855d55de8a1841208..8ecadbdbdf1a19eddaa336f5c6f82efc2a4fa1d9 100644 (file)
@@ -28,11 +28,21 @@ static inline float vg_randf(void)
    return (float)rand()/(float)(RAND_MAX);
 }
 
-static inline float vg_randint(int max)
+static inline int vg_randint(int max)
 {
    return rand()%max;
 }
 
+static float stable_force( float current, float diff )
+{
+   float fnew = current + diff;
+
+   if( fnew * current < 0.0f )
+      return 0.0f;
+
+   return fnew;
+}
+
 #define VG_MIN( A, B ) ((A)<(B)?(A):(B))
 #define VG_MAX( A, B ) ((A)>(B)?(A):(B))
 
@@ -177,6 +187,12 @@ static inline void v2_normalize( v2f a )
        v2_muls( a, 1.f / v2_length( a ), a );
 }
 
+static inline void v2_floor( v2f a, v2f b )
+{
+   b[0] = floorf( a[0] );
+   b[1] = floorf( a[1] );
+}
+
 /*
  * Vector 3
  */
@@ -368,7 +384,15 @@ static inline void v4_zero( v4f a )
        a[0] = 0.f; a[1] = 0.f; a[2] = 0.f; a[3] = 0.f;
 }
 
-static inline void v4_muladds( v3f a, v3f b, float s, v3f d )
+static inline void v4_muls( v4f a, float s, v4f d )
+{
+       d[0] = a[0]*s; 
+   d[1] = a[1]*s;
+   d[2] = a[2]*s;
+   d[3] = a[3]*s;
+}
+
+static inline void v4_muladds( v4f a, v4f b, float s, v4f d )
 {
        d[0] = a[0]+b[0]*s; 
    d[1] = a[1]+b[1]*s;
@@ -658,6 +682,20 @@ static inline void m4x3_mulv( m4x3f m, v3f v, v3f d )
        v3_copy( res, d );
 }
 
+/* 
+ * Transform plane ( xyz, distance )
+ */
+static inline void m4x3_mulp( m4x3f m, v4f p, v4f d )
+{
+   v3f o;
+
+   v3_muls( p, p[3], o );
+   m4x3_mulv( m, o, o );
+   m3x3_mulv( m, p, d );
+   
+   d[3] = v3_dot( o, d );
+}
+
 /*
  * Affine transforms
  */
@@ -676,6 +714,13 @@ static inline void m4x3_scale( m4x3f m, float s )
        v3_muls( m[2], s, m[2] );
 }
 
+static inline void m4x3_scalev( m4x3f m, v3f v )
+{
+   v3_muls(m[0], v[0], m[0]);
+   v3_muls(m[1], v[1], m[1]);
+   v3_muls(m[2], v[2], m[2]);
+}
+
 static inline void m4x3_rotate_x( m4x3f m, float angle )
 {
        m4x3f t = M4X3_IDENTITY;
@@ -763,6 +808,15 @@ static inline void box_copy( boxf a, boxf b )
        v3_copy( a[1], b[1] );
 }
 
+static inline int box_overlap( boxf a, boxf b )
+{
+   return
+       ( a[0][0] <= b[1][0] && a[1][0] >= b[0][0] ) &&
+       ( a[0][1] <= b[1][1] && a[1][1] >= b[0][1] ) &&
+       ( a[0][2] <= b[1][2] && a[1][2] >= b[0][2] )
+       ;
+}
+
 static inline void box_init_inf( boxf box )
 {
        v3_fill( box[0],  INFINITY );
@@ -778,14 +832,35 @@ static inline void m4x3_transform_aabb( m4x3f m, boxf box )
        v3_fill( box[0],  INFINITY );
        v3_fill( box[1], -INFINITY );
 
-       m4x3_expand_aabb_point( m, box, a );
+       m4x3_expand_aabb_point( m, box, (v3f){ a[0], a[1], a[2] } );
        m4x3_expand_aabb_point( m, box, (v3f){ a[0], b[1], a[2] } );
-       m4x3_expand_aabb_point( m, box, (v3f){ b[0], a[1], a[2] } );
        m4x3_expand_aabb_point( m, box, (v3f){ b[0], b[1], a[2] } );
-       m4x3_expand_aabb_point( m, box, b );
+       m4x3_expand_aabb_point( m, box, (v3f){ b[0], a[1], a[2] } );
+
+       m4x3_expand_aabb_point( m, box, (v3f){ a[0], a[1], b[2] } );
        m4x3_expand_aabb_point( m, box, (v3f){ a[0], b[1], b[2] } );
-       m4x3_expand_aabb_point( m, box, (v3f){ b[0], a[1], b[2] } );
        m4x3_expand_aabb_point( m, box, (v3f){ b[0], b[1], b[2] } );
+       m4x3_expand_aabb_point( m, box, (v3f){ b[0], a[1], b[2] } );
+}
+
+int ray_aabb( boxf box, v3f co, v3f dir, float dist )
+{
+   v3f v0, v1;
+   float tmin, tmax;
+
+   v3_sub( box[0], co, v0 );
+   v3_sub( box[1], co, v1 );
+   v3_div( v0, dir, v0 );
+   v3_div( v1, dir, v1 );
+   
+   tmin = vg_minf( v0[0], v1[0] );
+   tmax = vg_maxf( v0[0], v1[0] );
+   tmin = vg_maxf( tmin, vg_minf( v0[1], v1[1] ));
+   tmax = vg_minf( tmax, vg_maxf( v0[1], v1[1] ));
+   tmin = vg_maxf( tmin, vg_minf( v0[2], v1[2] ));
+   tmax = vg_minf( tmax, vg_maxf( v0[2], v1[2] ));
+
+   return tmax >= tmin && tmin < dist && tmax > 0;
 }
 
 static inline void m4x3_lookat( m4x3f m, v3f pos, v3f target, v3f up )
@@ -811,6 +886,10 @@ static inline void m4x3_lookat( m4x3f m, v3f pos, v3f target, v3f up )
                                                                { 0.0f, 1.0f, 0.0f, 0.0f },\
                                                                { 0.0f, 0.0f, 1.0f, 0.0f },\
                                                                { 0.0f, 0.0f, 0.0f, 1.0f }}
+#define M4X4_ZERO                      {{0.0f, 0.0f, 0.0f, 0.0f },\
+                                                               { 0.0f, 0.0f, 0.0f, 0.0f },\
+                                                               { 0.0f, 0.0f, 0.0f, 0.0f },\
+                                                               { 0.0f, 0.0f, 0.0f, 0.0f }}
 
 static void m4x4_projection( m4x4f m, float angle,
       float ratio, float near, float far )
@@ -860,6 +939,12 @@ static inline void m4x4_identity( m4x4f a )
        m4x4_copy( id, a );
 }
 
+static inline void m4x4_zero( m4x4f a )
+{
+       m4x4f zero = M4X4_ZERO;
+       m4x4_copy( zero, a );
+}
+
 static inline void m4x4_mul( m4x4f a, m4x4f b, m4x4f d )
 {
    float a00 = a[0][0], a01 = a[0][1], a02 = a[0][2], a03 = a[0][3],
@@ -890,6 +975,75 @@ static inline void m4x4_mul( m4x4f a, m4x4f b, m4x4f d )
   d[3][3] = a03*b30 + a13*b31 + a23*b32 + a33*b33;
 }
 
+static inline void m4x4_mulv( m4x4f m, v4f v, v4f d ) 
+{
+       v4f res;
+  
+       res[0] = m[0][0]*v[0] + m[1][0]*v[1] + m[2][0]*v[2] + m[3][0]*v[3];
+       res[1] = m[0][1]*v[0] + m[1][1]*v[1] + m[2][1]*v[2] + m[3][1]*v[3];
+       res[2] = m[0][2]*v[0] + m[1][2]*v[1] + m[2][2]*v[2] + m[3][2]*v[3];
+       res[3] = m[0][3]*v[0] + m[1][3]*v[1] + m[2][3]*v[2] + m[3][3]*v[3];
+
+       v4_copy( res, d );
+}
+
+static inline void m4x4_inv( m4x4f a, m4x4f d )
+{
+   float a00 = a[0][0], a01 = a[0][1], a02 = a[0][2], a03 = a[0][3],
+         a10 = a[1][0], a11 = a[1][1], a12 = a[1][2], a13 = a[1][3],
+         a20 = a[2][0], a21 = a[2][1], a22 = a[2][2], a23 = a[2][3],
+         a30 = a[3][0], a31 = a[3][1], a32 = a[3][2], a33 = a[3][3],
+         det,
+         t[6];
+
+   t[0] = a22*a33 - a32*a23;
+   t[1] = a21*a33 - a31*a23;
+   t[2] = a21*a32 - a31*a22;
+   t[3] = a20*a33 - a30*a23;
+   t[4] = a20*a32 - a30*a22;
+   t[5] = a20*a31 - a30*a21;
+
+   d[0][0] =  a11*t[0] - a12*t[1] + a13*t[2];
+   d[1][0] =-(a10*t[0] - a12*t[3] + a13*t[4]);
+   d[2][0] =  a10*t[1] - a11*t[3] + a13*t[5];
+   d[3][0] =-(a10*t[2] - a11*t[4] + a12*t[5]);
+
+   d[0][1] =-(a01*t[0] - a02*t[1] + a03*t[2]);
+   d[1][1] =  a00*t[0] - a02*t[3] + a03*t[4];
+   d[2][1] =-(a00*t[1] - a01*t[3] + a03*t[5]);
+   d[3][1] =  a00*t[2] - a01*t[4] + a02*t[5];
+
+   t[0] = a12*a33 - a32*a13;
+   t[1] = a11*a33 - a31*a13;
+   t[2] = a11*a32 - a31*a12;
+   t[3] = a10*a33 - a30*a13;
+   t[4] = a10*a32 - a30*a12;
+   t[5] = a10*a31 - a30*a11;
+
+   d[0][2] =  a01*t[0] - a02*t[1] + a03*t[2];
+   d[1][2] =-(a00*t[0] - a02*t[3] + a03*t[4]);
+   d[2][2] =  a00*t[1] - a01*t[3] + a03*t[5];
+   d[3][2] =-(a00*t[2] - a01*t[4] + a02*t[5]);
+
+   t[0] = a12*a23 - a22*a13;
+   t[1] = a11*a23 - a21*a13;
+   t[2] = a11*a22 - a21*a12;
+   t[3] = a10*a23 - a20*a13;
+   t[4] = a10*a22 - a20*a12;
+   t[5] = a10*a21 - a20*a11;
+
+   d[0][3] =-(a01*t[0] - a02*t[1] + a03*t[2]);
+   d[1][3] =  a00*t[0] - a02*t[3] + a03*t[4];
+   d[2][3] =-(a00*t[1] - a01*t[3] + a03*t[5]);
+   d[3][3] =  a00*t[2] - a01*t[4] + a02*t[5];
+   
+   det = 1.0f / (a00*d[0][0] + a01*d[1][0] + a02*d[2][0] + a03*d[3][0]);
+   v4_muls( d[0], det, d[0] );
+   v4_muls( d[1], det, d[1] );
+   v4_muls( d[2], det, d[2] );
+   v4_muls( d[3], det, d[3] );
+}
+
 /*
  * Planes (double precision)
  */
@@ -1073,3 +1227,48 @@ static void m3x3_q( m3x3f m, v4f q )
       q[3] = rinv * (m[0][1] - m[1][0]);
    }
 }
+
+static int ray_tri( v3f tri[3], v3f co, v3f dir, float *dist )
+{
+   float const kEpsilon = 0.00001f;
+
+   v3f v0, v1, h, s, q, n;
+   float a,f,u,v,t;
+
+   float *pa = tri[0],
+         *pb = tri[1],
+         *pc = tri[2];
+
+   v3_sub( pb, pa, v0 );
+   v3_sub( pc, pa, v1 );
+   v3_cross( dir, v1, h );
+   v3_cross( v0, v1, n );
+
+   if( v3_dot( n, dir ) > 0.0f ) /* Backface culling */
+      return 0;
+   
+   /* Parralel */
+   a = v3_dot( v0, h );
+   if( a > -kEpsilon && a < kEpsilon )
+      return 0;
+
+   f = 1.0f/a;
+   v3_sub( co, pa, s );
+
+   u = f * v3_dot(s, h);
+   if( u < 0.0f || u > 1.0f )
+      return 0;
+   
+   v3_cross( s, v0, q );
+   v = f * v3_dot( dir, q );
+   if( v < 0.0f || u+v > 1.0f )
+      return 0;
+
+   t = f * v3_dot(v1, q);
+   if( t > kEpsilon )
+   {
+      *dist = t;
+      return 1;
+   }
+   else return 0;
+}