1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
use num::Zero;
use alga::general::{Real, Id};
use na::{self, Vector3};
use query::algorithms::johnson_simplex::JohnsonSimplex;
use query::{Ray, RayCast, RayIntersection};
use query::ray_internal;
use shape::Triangle;
use math::{Point, Isometry};
use utils;
impl<P: Point, M: Isometry<P>> RayCast<P, M> for Triangle<P> {
#[inline]
fn toi_and_normal_with_ray(&self, m: &M, ray: &Ray<P>, solid: bool) -> Option<RayIntersection<P::Vector>> {
let ls_ray = ray.inverse_transform_by(m);
let res = if na::dimension::<P::Vector>() == 3 {
triangle_ray_intersection(self.a(), self.b(), self.c(), &ls_ray).map(|(r, _)| r)
}
else {
ray_internal::implicit_toi_and_normal_with_ray(&Id::new(), self,
&mut JohnsonSimplex::<P>::new_w_tls(),
&ls_ray,
solid)
};
res.map(|mut r| { r.normal = m.rotate_vector(&r.normal); r })
}
}
pub fn triangle_ray_intersection<P: Point>(a: &P, b: &P, c: &P, ray: &Ray<P>)
-> Option<(RayIntersection<P::Vector>, Vector3<P::Real>)> {
let ab = *b - *a;
let ac = *c - *a;
let n = utils::cross3(&ab, &ac);
let d = na::dot(&n, &ray.dir);
if d.is_zero() {
return None;
}
let ap = ray.origin - *a;
let t = na::dot(&ap, &n);
if (t < na::zero() && d < na::zero()) ||
(t > na::zero() && d > na::zero()) {
return None;
}
let d = d.abs();
let e = -utils::cross3(&ray.dir, &ap);
let mut v;
let mut w;
let toi;
let normal;
if t < na::zero() {
v = -na::dot(&ac, &e);
if v < na::zero() || v > d {
return None;
}
w = na::dot(&ab, &e);
if w < na::zero() || v + w > d {
return None;
}
let invd = na::one::<P::Real>() / d;
toi = -t * invd;
normal = -na::normalize(&n);
v = v * invd;
w = w * invd;
}
else {
v = na::dot(&ac, &e);
if v < na::zero() || v > d {
return None;
}
w = -na::dot(&ab, &e);
if w < na::zero() || v + w > d {
return None;
}
let invd = na::one::<P::Real>() / d;
toi = t * invd;
normal = na::normalize(&n);
v = v * invd;
w = w * invd;
}
Some((RayIntersection::new(toi, normal), Vector3::new(-v - w + na::one(), v, w)))
}