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