libs/mathlib/ray.c

   1 /*\r
   2 Copyright (C) 1999-2007 id Software, Inc. and contributors.\r
   3 For a list of contributors, see the accompanying CONTRIBUTORS file.\r
   4 \r
   5 This file is part of GtkRadiant.\r
   6 \r
   7 GtkRadiant is free software; you can redistribute it and/or modify\r
   8 it under the terms of the GNU General Public License as published by\r
   9 the Free Software Foundation; either version 2 of the License, or\r
  10 (at your option) any later version.\r
  11 \r
  12 GtkRadiant is distributed in the hope that it will be useful,\r
  13 but WITHOUT ANY WARRANTY; without even the implied warranty of\r
  14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the\r
  15 GNU General Public License for more details.\r
  16 \r
  17 You should have received a copy of the GNU General Public License\r
  18 along with GtkRadiant; if not, write to the Free Software\r
  19 Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301  USA\r
  20 */\r
  21 \r
  22 #include "mathlib.h"\r
  23 /*! for memcpy */\r
  24 #include <memory.h>\r
  25 \r
  26 vec3_t identity = { 0,0,0 };\r
  27 \r
  28 void ray_construct_for_vec3(ray_t *ray, const vec3_t origin, const vec3_t direction)\r
  29 {\r
  30   VectorCopy(origin, ray->origin);\r
  31   VectorCopy(direction, ray->direction);\r
  32 }\r
  33   \r
  34 void ray_transform(ray_t *ray, const m4x4_t matrix)\r
  35 {\r
  36   m4x4_transform_point(matrix, ray->origin);\r
  37   m4x4_transform_normal(matrix, ray->direction);\r
  38 }\r
  39 \r
  40 vec_t ray_intersect_point(const ray_t *ray, const vec3_t point, vec_t epsilon, vec_t divergence)\r
  41 {\r
  42   vec3_t displacement;\r
  43   vec_t depth;\r
  44 \r
  45   // calc displacement of test point from ray origin\r
  46   VectorSubtract(point, ray->origin, displacement);\r
  47   // calc length of displacement vector along ray direction\r
  48         depth = DotProduct(displacement, ray->direction);\r
  49   if(depth < 0.0f) return (vec_t)VEC_MAX;\r
  50   // calc position of closest point on ray to test point\r
  51   VectorMA (ray->origin, depth, ray->direction, displacement);\r
  52   // calc displacement of test point from closest point\r
  53         VectorSubtract(point, displacement, displacement);\r
  54   // calc length of displacement, subtract depth-dependant epsilon\r
  55   if (VectorLength(displacement) - (epsilon + (depth * divergence)) > 0.0f) return (vec_t)VEC_MAX;\r
  56   return depth;\r
  57 }\r
  58 \r
  59 // Tomas Moller and Ben Trumbore. Fast, minimum storage ray-triangle intersection. Journal of graphics tools, 2(1):21-28, 1997\r
  60 \r
  61 #define EPSILON 0.000001\r
  62 \r
  63 vec_t ray_intersect_triangle(const ray_t *ray, qboolean bCullBack, const vec3_t vert0, const vec3_t vert1, const vec3_t vert2)\r
  64 {\r
  65   float edge1[3], edge2[3], tvec[3], pvec[3], qvec[3];\r
  66   float det,inv_det;\r
  67   float u, v;\r
  68   vec_t depth = (vec_t)VEC_MAX;\r
  69   \r
  70   /* find vectors for two edges sharing vert0 */\r
  71   VectorSubtract(vert1, vert0, edge1);\r
  72   VectorSubtract(vert2, vert0, edge2);\r
  73   \r
  74   /* begin calculating determinant - also used to calculate U parameter */\r
  75   CrossProduct(ray->direction, edge2, pvec);\r
  76   \r
  77   /* if determinant is near zero, ray lies in plane of triangle */\r
  78   det = DotProduct(edge1, pvec);\r
  79   \r
  80   if (bCullBack == qtrue)\r
  81   {\r
  82     if (det < EPSILON)\r
  83       return depth;\r
  84     \r
  85     // calculate distance from vert0 to ray origin\r
  86     VectorSubtract(ray->origin, vert0, tvec);\r
  87     \r
  88     // calculate U parameter and test bounds\r
  89     u = DotProduct(tvec, pvec);\r
  90     if (u < 0.0 || u > det)\r
  91       return depth;\r
  92     \r
  93     // prepare to test V parameter\r
  94     CrossProduct(tvec, edge1, qvec);\r
  95     \r
  96     // calculate V parameter and test bounds\r
  97     v = DotProduct(ray->direction, qvec);\r
  98     if (v < 0.0 || u + v > det)\r
  99       return depth;\r
 100     \r
 101     // calculate t, scale parameters, ray intersects triangle\r
 102     depth = DotProduct(edge2, qvec);\r
 103     inv_det = 1.0f / det;\r
 104     depth *= inv_det;\r
 105     //u *= inv_det;\r
 106     //v *= inv_det;\r
 107   }\r
 108   else\r
 109   {\r
 110     /* the non-culling branch */\r
 111     if (det > -EPSILON && det < EPSILON)\r
 112       return depth;\r
 113     inv_det = 1.0f / det;\r
 114     \r
 115     /* calculate distance from vert0 to ray origin */\r
 116     VectorSubtract(ray->origin, vert0, tvec);\r
 117     \r
 118     /* calculate U parameter and test bounds */\r
 119     u = DotProduct(tvec, pvec) * inv_det;\r
 120     if (u < 0.0 || u > 1.0)\r
 121       return depth;\r
 122     \r
 123     /* prepare to test V parameter */\r
 124     CrossProduct(tvec, edge1, qvec);\r
 125     \r
 126     /* calculate V parameter and test bounds */\r
 127     v = DotProduct(ray->direction, qvec) * inv_det;\r
 128     if (v < 0.0 || u + v > 1.0)\r
 129       return depth;\r
 130     \r
 131     /* calculate t, ray intersects triangle */\r
 132     depth = DotProduct(edge2, qvec) * inv_det;\r
 133   }\r
 134   return depth;\r
 135 }\r