1 /*
2 * MOAB, a Mesh-Oriented datABase, is a software component for creating,
3 * storing and accessing finite element mesh data.
4 *
5 * Copyright 2004 Sandia Corporation. Under the terms of Contract
6 * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government
7 * retains certain rights in this software.
8 *
9 * This library is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Lesser General Public
11 * License as published by the Free Software Foundation; either
12 * version 2.1 of the License, or (at your option) any later version.
13 *
14 */
15
16 /**\file Geometry.hpp
17 *\author Jason Kraftcheck (kraftche@cae.wisc.edu)
18 *\date 2006-07-27
19 */
20
21 #ifndef MB_GEOM_UTIL_HPP
22 #define MB_GEOM_UTIL_HPP
23
24 #include "moab/CartVect.hpp"
25 #include <cmath>
26
27 namespace moab
28 {
29
30 /** \class GeomUtil
31 * \brief Functions for computational geometry on triangles, rays, and boxes
32 */
33 namespace GeomUtil
34 {
35
36 /** Given a line segment and an axis-aligned box,
37 * return the sub-segment of the line segment that
38 * itersects the box.
39 *
40 * Can be used to intersect ray with box by passing seg_end
41 * as HUGE_VAL or std::numeric_limits<double>::maximum().
42 *
43 *\param box_min Minimum corner of axis-aligned box
44 *\param box_max Maximum corner of axis-aligned box
45 *\param seg_pt A point in the line containing the segement
46 *\param seg_unit_dir A unit vector in the direction of the line
47 * containing the semgent.
48 *\param seg_start The distance from seg_pt in the direction of
49 * seg_unit_dir at which the segment begins.
50 * As input, the start of the original segment, as output, the
51 * start of the sub-segment intersecting the box.
52 * Note: seg_start must be less than seg_end
53 *\param seg_end The distance from seg_pt in the direction of
54 * seg_unit_dir at which the segment ends.
55 * As input, the end of the original segment, as output, the
56 * end of the sub-segment intersecting the box.
57 * Note: seg_start must be less than seg_end
58 *\return true if line semgent intersects box, false otherwise.
59 */
60 bool segment_box_intersect( CartVect box_min,
61 CartVect box_max,
62 const CartVect& seg_pt,
63 const CartVect& seg_unit_dir,
64 double& seg_start,
65 double& seg_end );
66
67 /**\brief Test for intersection between a ray and a triangle.
68 *\param ray_point The start point of the ray.
69 *\param ray_unit_direciton The direction of the ray. Must be a unit vector.
70 *\param t_out Output: The distance along the ray from ray_point in the
71 * direction of ray_unit_direction at which the ray
72 * itersected the triangle.
73 *\param ray_length Optional: If non-null, a pointer a maximum length
74 * for the ray, converting this function to a segment-tri-
75 * intersect test.
76 *\return true if intersection, false otherwise.
77 */
78 bool ray_tri_intersect( const CartVect vertices[3],
79 const CartVect& ray_point,
80 const CartVect& ray_unit_direction,
81 double& t_out,
82 const double* ray_length = 0 );
83
84 /**\brief Plücker test for intersection between a ray and a triangle.
85 *\param vertices Nodes of the triangle.
86 *\param ray_point The start point of the ray.
87 *\param ray_unit_direction The direction of the ray. Must be a unit vector.
88 *\param t_out Output: The distance along the ray from ray_point in the
89 * direction of ray_unit_direction at which the ray
90 * intersected the triangle.
91 *\param nonneg_ray_length Optional: If non-null, a maximum length for the ray,
92 * converting this function to a segment-tri-intersect
93 * test.
94 *\param neg_ray_length Optional: If non-null, a maximum length for the ray
95 * behind the origin, converting this function to a
96 * segment-tri-intersect test.
97 *\param orientation Optional: Reject intersections without the specified
98 * orientation of ray with respect to triangle normal
99 * vector. Indicate desired orientation by passing
100 * 1 (forward), -1 (reverse), or 0 (no preference).
101 *\param int_type Optional Output: The type of intersection; used to
102 * identify edge/node intersections.
103 *\return true if intersection, false otherwise.
104 */
105 enum intersection_type
106 {
107 NONE = 0,
108 INTERIOR,
109 NODE0,
110 NODE1,
111 NODE2,
112 EDGE0,
113 EDGE1,
114 EDGE2
115 };
116 /* intersection type is determined by which of the intermediate values are == 0. There
117 are three such values that can therefore be encoded in a 3 bit integer.
118 0 = none are == 0 -> interior type
119 1 = pip0 == 0 -> EDGE0
120 2 = pip1 == 1 -> EDGE1
121 4 = pip2 == 2 -> EDGE2
122 5 = pip2 = pip0 == 0 -> NOEE0
123 3 = pip0 = pip1 == 0 -> NODE1
124 6 = pip1 = pip2 == 0 -> NODE2 */
125 const intersection_type type_list[] = { INTERIOR, EDGE0, EDGE1, NODE1, EDGE2, NODE0, NODE2 };
126
127 bool plucker_ray_tri_intersect( const CartVect vertices[3],
128 const CartVect& ray_point,
129 const CartVect& ray_unit_direction,
130 double& dist_out,
131 const double* nonneg_ray_length = 0,
132 const double* neg_ray_length = 0,
133 const int* orientation = 0,
134 intersection_type* int_type = 0 );
135 double plucker_edge_test( const CartVect& vertexa,
136 const CartVect& vertexb,
137 const CartVect& ray,
138 const CartVect& ray_normal );
139
140 //! Find range of overlap between ray and axis-aligned box.
141 //!
142 //!\param box_min Box corner with minimum coordinate values
143 //!\param box_max Box corner with minimum coordinate values
144 //!\param ray_pt Coordinates of start point of ray
145 //!\param ray_dir Directionion vector for ray such that
146 //! the ray is defined by r(t) = ray_point + t * ray_vect
147 //! for t > 0.
148 //!\param t_enter Output: if return value is true, this value
149 //! is the parameter location along the ray at which
150 //! the ray entered the leaf. If return value is false,
151 //! then this value is undefined.
152 //!\param t_exit Output: if return value is true, this value
153 //! is the parameter location along the ray at which
154 //! the ray exited the leaf. If return value is false,
155 //! then this value is undefined.
156 //!\return true if ray intersects leaf, false otherwise.
157 bool ray_box_intersect( const CartVect& box_min,
158 const CartVect& box_max,
159 const CartVect& ray_pt,
160 const CartVect& ray_dir,
161 double& t_enter,
162 double& t_exit );
163
164 /**\brief Test if plane intersects axis-aligned box
165 *
166 * Test for intersection between an unbounded plane and
167 * an axis-aligned box.
168 *\param plane_normal Vector in plane normal direction (need *not*
169 * be a unit vector). The N in
170 * the plane equation: N . X + D = 0
171 *\param plane_coeff The scalar 'D' term in the plane equation:
172 * N . X + D = 0
173 *\param box_min_corner The smallest coordinates of the box along each
174 * axis. The corner of the box for which all three
175 * coordinate values are smaller than those of any
176 * other corner. The X, Y, Z values for the planes
177 * normal to those axes and bounding the box on the
178 * -X, -Y, and -Z sides respectively.
179 *\param box_max_corner The largest coordinates of the box along each
180 * axis. The corner of the box for which all three
181 * coordinate values are larger than those of any
182 * other corner. The X, Y, Z values for the planes
183 * normal to those axes and bounding the box on the
184 * +X, +Y, and +Z sides respectively.
185 *\return true if overlap, false otherwise.
186 */
187 bool box_plane_overlap( const CartVect& plane_normal,
188 double plane_coeff,
189 CartVect box_min_corner,
190 CartVect box_max_corner );
191
192 /**\brief Test if triangle intersects axis-aligned box
193 *
194 * Test if a triangle intersects an axis-aligned box.
195 *\param triangle_corners The corners of the triangle.
196 *\param box_min_corner The smallest coordinates of the box along each
197 * axis. The corner of the box for which all three
198 * coordinate values are smaller than those of any
199 * other corner. The X, Y, Z values for the planes
200 * normal to those axes and bounding the box on the
201 * -X, -Y, and -Z sides respectively.
202 *\param box_max_corner The largest coordinates of the box along each
203 * axis. The corner of the box for which all three
204 * coordinate values are larger than those of any
205 * other corner. The X, Y, Z values for the planes
206 * normal to those axes and bounding the box on the
207 * +X, +Y, and +Z sides respectively.
208 *\param tolerance The tolerance used in the intersection test. The box
209 * size is increased by this amount before the intersection
210 * test.
211 *\return true if overlap, false otherwise.
212 */
213 bool box_tri_overlap( const CartVect triangle_corners[3],
214 const CartVect& box_min_corner,
215 const CartVect& box_max_corner,
216 double tolerance );
217
218 /**\brief Test if triangle intersects axis-aligned box
219 *
220 * Test if a triangle intersects an axis-aligned box.
221 *\param triangle_corners The corners of the triangle.
222 *\param box_center The center of the box.
223 *\param box_hanf_dims The distance along each axis, respectively, from the
224 * box_center to the boundary of the box.
225 *\return true if overlap, false otherwise.
226 */
227 bool box_tri_overlap( const CartVect triangle_corners[3],
228 const CartVect& box_center,
229 const CartVect& box_half_dims );
230
231 bool box_point_overlap( const CartVect& box_min_corner,
232 const CartVect& box_max_corner,
233 const CartVect& point,
234 double tolerance );
235
236 /**\brief Test if the specified element intersects an axis-aligned box.
237 *
238 * Test if element intersects axis-aligned box. Use element-specific
239 * optimization if available, otherwise call box_general_elem_overlap.
240 *
241 *\param elem_corners The coordinates of the element vertices
242 *\param elem_type The toplogy of the element.
243 *\param box_center The center of the axis-aligned box
244 *\param box_half_dims Half of the width of the box in each axial
245 * direction.
246 */
247 bool box_elem_overlap( const CartVect* elem_corners,
248 EntityType elem_type,
249 const CartVect& box_center,
250 const CartVect& box_half_dims,
251 int nodecount = 0 );
252
253 /**\brief Test if the specified element intersects an axis-aligned box.
254 *
255 * Uses MBCN and separating axis theorem for general algorithm that
256 * works for all fixed-size elements (not poly*).
257 *
258 *\param elem_corners The coordinates of the element vertices
259 *\param elem_type The toplogy of the element.
260 *\param box_center The center of the axis-aligned box
261 *\param box_half_dims Half of the width of the box in each axial
262 * direction.
263 */
264 bool box_linear_elem_overlap( const CartVect* elem_corners,
265 EntityType elem_type,
266 const CartVect& box_center,
267 const CartVect& box_half_dims );
268
269 /**\brief Test if the specified element intersects an axis-aligned box.
270 *
271 * Uses MBCN and separating axis theorem for general algorithm that
272 * works for all fixed-size elements (not poly*). Box and element
273 * vertices must be translated such that box center is at origin.
274 *
275 *\param elem_corners The coordinates of the element vertices, in
276 * local coordinate system of box.
277 *\param elem_type The toplogy of the element.
278 *\param box_half_dims Half of the width of the box in each axial
279 * direction.
280 */
281 bool box_linear_elem_overlap( const CartVect* elem_corners, EntityType elem_type, const CartVect& box_half_dims );
282
283 void closest_location_on_box( const CartVect& box_min_corner,
284 const CartVect& box_max_corner,
285 const CartVect& point,
286 CartVect& closest );
287
288 /**\brief find closest location on triangle
289 *
290 * Find closest location on linear triangle.
291 *\param location Input position to evaluate from
292 *\param vertices Array of three corner vertex coordinates.
293 *\param closest_out Result position
294 */
295 void closest_location_on_tri( const CartVect& location, const CartVect* vertices, CartVect& closest_out );
296
297 /**\brief find closest location on polygon
298 *
299 * Find closest location on polygon
300 *\param location Input position to evaluate from
301 *\param vertices Array of corner vertex coordinates.
302 *\param num_vertices Length of 'vertices' array.
303 *\param closest_out Result position
304 */
305 void closest_location_on_polygon( const CartVect& location,
306 const CartVect* vertices,
307 int num_vertices,
308 CartVect& closest_out );
309
310 /**\brief find closest topological location on triangle
311 *
312 * Find closest location on linear triangle.
313 *\param location Input position to evaluate from
314 *\param vertices Array of three corner vertex coordinates.
315 *\param tolerance Tolerance to use when comparing to corners and edges
316 *\param closest_out Result position
317 *\param closest_topo Closest topological entity
318 * 0-2 : vertex index
319 * 3-5 : edge beginning at closest_topo - 3
320 * 6 : triangle interior
321 */
322 void closest_location_on_tri( const CartVect& location,
323 const CartVect* vertices,
324 double tolerance,
325 CartVect& closest_out,
326 int& closest_topo );
327
328 // Finds whether or not a box defined by the center and the half
329 // width intersects a trilinear hex defined by its eight vertices.
330 bool box_hex_overlap( const CartVect hexv[8], const CartVect& box_center, const CartVect& box_dims );
331
332 // Finds whether or not a box defined by the center and the half
333 // width intersects a linear tetrahedron defined by its four vertices.
334 bool box_tet_overlap( const CartVect tet_corners[4], const CartVect& box_center, const CartVect& box_dims );
335
336 // tests if 2 boxes overlap within a tolerance
337 // assume that data is valid, box_min1 << box_max1, and box_min2<< box_max2
338 // they are overlapping if they are overlapping in all directions
339 // for example in x direction:
340 // box_max2[0] + tolerance < box_min1[0] -- false
341 bool boxes_overlap( const CartVect& box_min1,
342 const CartVect& box_max1,
343 const CartVect& box_min2,
344 const CartVect& box_max2,
345 double tolerance );
346
347 // see if boxes formed by 2 lists of "CartVect"s overlap
348 bool bounding_boxes_overlap( const CartVect* list1, int num1, const CartVect* list2, int num2, double tolerance );
349
350 // see if boxes from vertices in 2d overlap (used in gnomonic planes right now)
351 bool bounding_boxes_overlap_2d( const double* list1, int num1, const double* list2, int num2, double tolerance );
352 // point_in_trilinear_hex
353 // Tests if a point in xyz space is within a hex element defined with
354 // its eight vertex points forming a trilinear basis function. Computes
355 // the natural coordinates with respect to the hex of the xyz point
356 // and checks if each are between +/-1. If anyone is outside the range
357 // the function returns false, otherwise it returns true.
358 //
359 bool point_in_trilinear_hex( const CartVect* hex, const CartVect& xyz, double etol );
360
361 bool nat_coords_trilinear_hex( const CartVect* corner_coords, const CartVect& x, CartVect& xi, double tol );
362 } // namespace GeomUtil
363
364 } // namespace moab
365
366 #endif