Mesh Oriented datABase  (version 5.5.1)
An array-based unstructured mesh library
GeomUtil.hpp
Go to the documentation of this file.
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 ([email protected])
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  */
106  {
107  NONE = 0,
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 */
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