Coverage for trimesh/ray/ray_triangle.py: 91%
100 statements
« prev ^ index » next coverage.py v7.14.1, created at 2026-06-24 04:40 +0000
1"""
2A basic slow implementation of ray- triangle queries.
3"""
5import numpy as np
7from .. import caching, grouping, intersections, util
8from .. import triangles as triangles_mod
9from ..constants import tol
10from .ray_util import contains_points
13class RayMeshIntersector:
14 """
15 An object to query a mesh for ray intersections.
16 Precomputes an r-tree for each triangle on the mesh.
17 """
19 def __init__(self, mesh):
20 self.mesh = mesh
21 self._cache = caching.Cache(self.mesh.__hash__)
23 def intersects_id(
24 self,
25 ray_origins,
26 ray_directions,
27 return_locations=False,
28 multiple_hits=True,
29 **kwargs,
30 ):
31 """
32 Find the intersections between the current mesh and an
33 array of rays.
35 Parameters
36 ------------
37 ray_origins : (m, 3) float
38 Ray origin points
39 ray_directions : (m, 3) float
40 Ray direction vectors
41 multiple_hits : bool
42 Consider multiple hits of each ray or not
43 return_locations : bool
44 Return hit locations or not
46 Returns
47 -----------
48 index_triangle : (h,) int
49 Index of triangles hit
50 index_ray : (h,) int
51 Index of ray that hit triangle
52 locations : (h, 3) float
53 [optional] Position of intersection in space
54 """
55 (index_tri, index_ray, locations) = ray_triangle_id(
56 triangles=self.mesh.triangles,
57 ray_origins=ray_origins,
58 ray_directions=ray_directions,
59 tree=self.mesh.triangles_tree,
60 multiple_hits=multiple_hits,
61 triangles_normal=self.mesh.face_normals,
62 )
63 if return_locations:
64 if len(index_tri) == 0:
65 return index_tri, index_ray, locations
66 unique = grouping.unique_rows(np.column_stack((locations, index_ray)))[0]
67 return index_tri[unique], index_ray[unique], locations[unique]
68 return index_tri, index_ray
70 def intersects_location(self, ray_origins, ray_directions, **kwargs):
71 """
72 Return unique cartesian locations where rays hit the mesh.
73 If you are counting the number of hits a ray had, this method
74 should be used as if only the triangle index is used on- edge hits
75 will be counted twice.
77 Parameters
78 ------------
79 ray_origins : (m, 3) float
80 Ray origin points
81 ray_directions : (m, 3) float
82 Ray direction vectors
84 Returns
85 ---------
86 locations : (n) sequence of (m,3) float
87 Intersection points
88 index_ray : (n,) int
89 Array of ray indexes
90 index_tri: (n,) int
91 Array of triangle (face) indexes
92 """
93 (index_tri, index_ray, locations) = self.intersects_id(
94 ray_origins=ray_origins,
95 ray_directions=ray_directions,
96 return_locations=True,
97 **kwargs,
98 )
99 return locations, index_ray, index_tri
101 def intersects_first(self, ray_origins, ray_directions, **kwargs):
102 """
103 Find the index of the first triangle a ray hits.
106 Parameters
107 ----------
108 ray_origins : (n, 3) float
109 Origins of rays
110 ray_directions : (n, 3) float
111 Direction (vector) of rays
113 Returns
114 ----------
115 triangle_index : (n,) int
116 Index of triangle ray hit, or -1 if not hit
117 """
119 (index_tri, index_ray) = self.intersects_id(
120 ray_origins=ray_origins,
121 ray_directions=ray_directions,
122 return_locations=False,
123 multiple_hits=False,
124 **kwargs,
125 )
127 # put the result into the form of "one triangle index per ray"
128 result = np.ones(len(ray_origins), dtype=np.int64) * -1
129 result[index_ray] = index_tri
131 return result
133 def intersects_any(self, ray_origins, ray_directions, **kwargs):
134 """
135 Find out if each ray hit any triangle on the mesh.
137 Parameters
138 ------------
139 ray_origins : (m, 3) float
140 Ray origin points
141 ray_directions : (m, 3) float
142 Ray direction vectors
144 Returns
145 ---------
146 hit : (m,) bool
147 Whether any ray hit any triangle on the mesh
148 """
149 _index_tri, index_ray = self.intersects_id(ray_origins, ray_directions)
150 hit_any = np.zeros(len(ray_origins), dtype=bool)
151 hit_idx = np.unique(index_ray)
152 if len(hit_idx) > 0:
153 hit_any[hit_idx] = True
154 return hit_any
156 def contains_points(self, points):
157 """
158 Check if a mesh contains a list of points, using ray tests.
160 If the point is on the surface of the mesh the behavior
161 is undefined.
163 Parameters
164 ------------
165 points : (n, 3) float
166 Points in space
168 Returns
169 ---------
170 contains : (n,) bool
171 Whether point is inside mesh or not
172 """
174 return contains_points(self, points)
177def ray_triangle_id(
178 triangles,
179 ray_origins,
180 ray_directions,
181 triangles_normal=None,
182 tree=None,
183 multiple_hits=True,
184):
185 """
186 Find the intersections between a group of triangles and rays
188 Parameters
189 -------------
190 triangles : (n, 3, 3) float
191 Triangles in space
192 ray_origins : (m, 3) float
193 Ray origin points
194 ray_directions : (m, 3) float
195 Ray direction vectors
196 triangles_normal : (n, 3) float
197 Normal vector of triangles, optional
198 tree : rtree.Index
199 Rtree object holding triangle bounds
201 Returns
202 -----------
203 index_triangle : (h,) int
204 Index of triangles hit
205 index_ray : (h,) int
206 Index of ray that hit triangle
207 locations : (h, 3) float
208 Position of intersection in space
209 """
210 triangles = np.asanyarray(triangles, dtype=np.float64)
211 ray_origins = np.asanyarray(ray_origins, dtype=np.float64)
212 ray_directions = np.asanyarray(ray_directions, dtype=np.float64)
214 # if we didn't get passed an r-tree for the bounds of each
215 # triangle create one here
216 if tree is None:
217 tree = triangles_mod.bounds_tree(triangles)
219 # find the list of likely triangles and which ray they
220 # correspond with, via rtree queries
221 ray_candidates, ray_id = ray_triangle_candidates(
222 ray_origins=ray_origins, ray_directions=ray_directions, tree=tree
223 )
225 # get subsets which are corresponding rays and triangles
226 # (c,3,3) triangle candidates
227 triangle_candidates = triangles[ray_candidates]
228 # (c,3) origins and vectors for the rays
229 line_origins = ray_origins[ray_id]
230 line_directions = ray_directions[ray_id]
232 # get the plane origins and normals from the triangle candidates
233 plane_origins = triangle_candidates[:, 0, :]
234 if triangles_normal is None:
235 plane_normals, triangle_ok = triangles_mod.normals(triangle_candidates)
236 if not triangle_ok.all():
237 raise ValueError("Invalid triangles!")
238 else:
239 plane_normals = triangles_normal[ray_candidates]
241 # find the intersection location of the rays with the planes
242 location, valid = intersections.planes_lines(
243 plane_origins=plane_origins,
244 plane_normals=plane_normals,
245 line_origins=line_origins,
246 line_directions=line_directions,
247 )
249 if len(triangle_candidates) == 0 or not valid.any():
250 # we got no hits so return early with empty array
251 return (
252 np.array([], dtype=np.int64),
253 np.array([], dtype=np.int64),
254 np.array([], dtype=np.float64),
255 )
257 # find the barycentric coordinates of each plane intersection on the
258 # triangle candidates
259 barycentric = triangles_mod.points_to_barycentric(
260 triangle_candidates[valid], location
261 )
263 # the plane intersection is inside the triangle if all barycentric
264 # coordinates are between 0.0 and 1.0
265 hit = np.logical_and(
266 (barycentric > -tol.zero).all(axis=1), (barycentric < (1 + tol.zero)).all(axis=1)
267 )
269 # the result index of the triangle is a candidate with a valid
270 # plane intersection and a triangle which contains the plane
271 # intersection point
272 index_tri = ray_candidates[valid][hit]
273 # the ray index is a subset with a valid plane intersection and
274 # contained by a triangle
275 index_ray = ray_id[valid][hit]
276 # locations are already valid plane intersections, just mask by hits
277 location = location[hit]
279 # only return points that are forward from the origin
280 vector = location - ray_origins[index_ray]
281 distance = util.diagonal_dot(vector, ray_directions[index_ray])
282 forward = distance > -1e-6
284 index_tri = index_tri[forward]
285 index_ray = index_ray[forward]
286 location = location[forward]
287 distance = distance[forward]
289 if multiple_hits:
290 return index_tri, index_ray, location
292 # since we are not returning multiple hits, we need to
293 # figure out which hit is first
294 if len(index_ray) == 0:
295 return index_tri, index_ray, location
297 # find the first hit
298 first = np.array([g[distance[g].argmin()] for g in grouping.group(index_ray)])
300 return index_tri[first], index_ray[first], location[first]
303def ray_triangle_candidates(ray_origins, ray_directions, tree):
304 """
305 Do broad- phase search for triangles that the rays
306 may intersect.
308 Does this by creating a bounding box for the ray as it
309 passes through the volume occupied by the tree
311 Parameters
312 ------------
313 ray_origins : (m, 3) float
314 Ray origin points.
315 ray_directions : (m, 3) float
316 Ray direction vectors
317 tree : rtree object
318 Ccontains AABB of each triangle
320 Returns
321 ----------
322 ray_candidates : (n,) int
323 Triangle indexes
324 ray_id : (n,) int
325 Corresponding ray index for a triangle candidate
326 """
327 bounding = ray_bounds(
328 ray_origins=ray_origins, ray_directions=ray_directions, bounds=tree.bounds
329 )
331 index = []
332 candidates = []
333 for i, bounds in enumerate(bounding):
334 cand = list(tree.intersection(bounds))
335 candidates.extend(cand)
336 index.extend([i] * len(cand))
337 return np.array(candidates, dtype=np.int64), np.array(index, dtype=np.int64)
340def ray_bounds(ray_origins, ray_directions, bounds, buffer_dist=1e-5):
341 """
342 Given a set of rays and a bounding box for the volume of interest
343 where the rays will be passing through, find the bounding boxes
344 of the rays as they pass through the volume.
346 Parameters
347 ------------
348 ray_origins: (m,3) float, ray origin points
349 ray_directions: (m,3) float, ray direction vectors
350 bounds: (2,3) bounding box (min, max)
351 buffer_dist: float, distance to pad zero width bounding boxes
353 Returns
354 ---------
355 ray_bounding: (n) set of AABB of rays passing through volume
356 """
358 ray_origins = np.asanyarray(ray_origins, dtype=np.float64)
359 ray_directions = np.asanyarray(ray_directions, dtype=np.float64)
361 # bounding box we are testing against
362 bounds = np.asanyarray(bounds)
364 # find the primary axis of the vector
365 axis = np.abs(ray_directions).argmax(axis=1)
366 axis_bound = bounds.reshape((2, -1)).T[axis]
367 axis_ori = np.array([ray_origins[i][a] for i, a in enumerate(axis)]).reshape((-1, 1))
368 axis_dir = np.array([ray_directions[i][a] for i, a in enumerate(axis)]).reshape(
369 (-1, 1)
370 )
372 # parametric equation of a line
373 # point = direction*t + origin
374 # p = dt + o
375 # t = (p-o)/d
376 nonzero = (axis_dir != 0.0).reshape(-1)
377 t = np.zeros_like(axis_bound)
378 t[nonzero] = (axis_bound[nonzero] - axis_ori[nonzero]) / axis_dir[nonzero]
380 # prevent the bounding box from including triangles
381 # behind the ray origin
382 t[t < buffer_dist] = buffer_dist
384 # the value of t for both the upper and lower bounds
385 t_a = t[:, 0].reshape((-1, 1))
386 t_b = t[:, 1].reshape((-1, 1))
388 # the cartesian point for where the line hits the plane defined by
389 # axis
390 on_a = (ray_directions * t_a) + ray_origins
391 on_b = (ray_directions * t_b) + ray_origins
393 on_plane = np.column_stack((on_a, on_b)).reshape((-1, 2, ray_directions.shape[1]))
395 ray_bounding = np.hstack((on_plane.min(axis=1), on_plane.max(axis=1)))
396 # pad the bounding box by TOL_BUFFER
397 # not sure if this is necessary, but if the ray is axis aligned
398 # this function will otherwise return zero volume bounding boxes
399 # which may or may not screw up the r-tree intersection queries
400 ray_bounding += np.array([-1, -1, -1, 1, 1, 1]) * buffer_dist
402 return ray_bounding