Mesh Oriented datABase  (version 5.4.1)
Array-based unstructured mesh datastructure
GeomUtil.hpp
Go to the documentation of this file.
00001 /*
00002  * MOAB, a Mesh-Oriented datABase, is a software component for creating,
00003  * storing and accessing finite element mesh data.
00004  *
00005  * Copyright 2004 Sandia Corporation.  Under the terms of Contract
00006  * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government
00007  * retains certain rights in this software.
00008  *
00009  * This library is free software; you can redistribute it and/or
00010  * modify it under the terms of the GNU Lesser General Public
00011  * License as published by the Free Software Foundation; either
00012  * version 2.1 of the License, or (at your option) any later version.
00013  *
00014  */
00015 
00016 /**\file Geometry.hpp
00017  *\author Jason Kraftcheck ([email protected])
00018  *\date 2006-07-27
00019  */
00020 
00021 #ifndef MB_GEOM_UTIL_HPP
00022 #define MB_GEOM_UTIL_HPP
00023 
00024 #include "moab/CartVect.hpp"
00025 #include <cmath>
00026 
00027 namespace moab
00028 {
00029 
00030 /** \class GeomUtil
00031  * \brief Functions for computational geometry on triangles, rays, and boxes
00032  */
00033 namespace GeomUtil
00034 {
00035 
00036     /** Given a line segment and an axis-aligned box,
00037      *  return the sub-segment of the line segment that
00038      *  itersects the box.
00039      *
00040      *  Can be used to intersect ray with box by passing seg_end
00041      *  as HUGE_VAL or std::numeric_limits<double>::maximum().
00042      *
00043      *\param box_min   Minimum corner of axis-aligned box
00044      *\param box_max   Maximum corner of axis-aligned box
00045      *\param seg_pt    A point in the line containing the segement
00046      *\param seg_unit_dir A unit vector in the direction of the line
00047      *                 containing the semgent.
00048      *\param seg_start The distance from seg_pt in the direction of
00049      *                 seg_unit_dir at which the segment begins.
00050      *                 As input, the start of the original segment, as output, the
00051      *                 start of the sub-segment intersecting the box.
00052      *                 Note:  seg_start must be less than seg_end
00053      *\param seg_end   The distance from seg_pt in the direction of
00054      *                 seg_unit_dir at which the segment ends.
00055      *                 As input, the end of the original segment, as output, the
00056      *                 end of the sub-segment intersecting the box.
00057      *                 Note:  seg_start must be less than seg_end
00058      *\return true if line semgent intersects box, false otherwise.
00059      */
00060     bool segment_box_intersect( CartVect box_min,
00061                                 CartVect box_max,
00062                                 const CartVect& seg_pt,
00063                                 const CartVect& seg_unit_dir,
00064                                 double& seg_start,
00065                                 double& seg_end );
00066 
00067     /**\brief Test for intersection between a ray and a triangle.
00068      *\param ray_point  The start point of the ray.
00069      *\param ray_unit_direciton  The direction of the ray. Must be a unit vector.
00070      *\param t_out Output: The distance along the ray from ray_point in the
00071      *                  direction of ray_unit_direction at which the ray
00072      *                  itersected the triangle.
00073      *\param ray_length Optional:  If non-null, a pointer a maximum length
00074      *                  for the ray, converting this function to a segment-tri-
00075      *                  intersect test.
00076      *\return true if intersection, false otherwise.
00077      */
00078     bool ray_tri_intersect( const CartVect vertices[3],
00079                             const CartVect& ray_point,
00080                             const CartVect& ray_unit_direction,
00081                             double& t_out,
00082                             const double* ray_length = 0 );
00083 
00084     /**\brief Plücker test for intersection between a ray and a triangle.
00085      *\param vertices            Nodes of the triangle.
00086      *\param ray_point           The start point of the ray.
00087      *\param ray_unit_direction  The direction of the ray. Must be a unit vector.
00088      *\param t_out               Output: The distance along the ray from ray_point in the
00089      *                           direction of ray_unit_direction at which the ray
00090      *                           intersected the triangle.
00091      *\param nonneg_ray_length   Optional: If non-null, a maximum length for the ray,
00092      *                           converting this function to a segment-tri-intersect
00093      *                           test.
00094      *\param neg_ray_length      Optional: If non-null, a maximum length for the ray
00095      *                           behind the origin, converting this function to a
00096      *                           segment-tri-intersect test.
00097      *\param orientation         Optional: Reject intersections without the specified
00098      *                           orientation of ray with respect to triangle normal
00099      *                           vector. Indicate desired orientation by passing
00100      *                           1 (forward), -1 (reverse), or 0 (no preference).
00101      *\param int_type            Optional Output: The type of intersection; used to
00102      *                           identify edge/node intersections.
00103      *\return true if intersection, false otherwise.
00104      */
00105     enum intersection_type
00106     {
00107         NONE = 0,
00108         INTERIOR,
00109         NODE0,
00110         NODE1,
00111         NODE2,
00112         EDGE0,
00113         EDGE1,
00114         EDGE2
00115     };
00116     /* intersection type is determined by which of the intermediate values are == 0.  There
00117        are three such values that can therefore be encoded in a 3 bit integer.
00118        0 = none are == 0 -> interior type
00119        1 = pip0 == 0 -> EDGE0
00120        2 = pip1 == 1 -> EDGE1
00121        4 = pip2 == 2 -> EDGE2
00122        5 = pip2 = pip0 == 0 -> NOEE0
00123        3 = pip0 = pip1 == 0 -> NODE1
00124        6 = pip1 = pip2 == 0 -> NODE2 */
00125     const intersection_type type_list[] = { INTERIOR, EDGE0, EDGE1, NODE1, EDGE2, NODE0, NODE2 };
00126 
00127     bool plucker_ray_tri_intersect( const CartVect vertices[3],
00128                                     const CartVect& ray_point,
00129                                     const CartVect& ray_unit_direction,
00130                                     double& dist_out,
00131                                     const double* nonneg_ray_length = 0,
00132                                     const double* neg_ray_length    = 0,
00133                                     const int* orientation          = 0,
00134                                     intersection_type* int_type     = 0 );
00135     double plucker_edge_test( const CartVect& vertexa,
00136                               const CartVect& vertexb,
00137                               const CartVect& ray,
00138                               const CartVect& ray_normal );
00139 
00140     //! Find range of overlap between ray and axis-aligned box.
00141     //!
00142     //!\param box_min   Box corner with minimum coordinate values
00143     //!\param box_max   Box corner with minimum coordinate values
00144     //!\param ray_pt    Coordinates of start point of ray
00145     //!\param ray_dir   Directionion vector for ray such that
00146     //!                 the ray is defined by r(t) = ray_point + t * ray_vect
00147     //!                 for t > 0.
00148     //!\param t_enter   Output: if return value is true, this value
00149     //!                 is the parameter location along the ray at which
00150     //!                 the ray entered the leaf.  If return value is false,
00151     //!                 then this value is undefined.
00152     //!\param t_exit    Output: if return value is true, this value
00153     //!                 is the parameter location along the ray at which
00154     //!                 the ray exited the leaf.  If return value is false,
00155     //!                 then this value is undefined.
00156     //!\return true if ray intersects leaf, false otherwise.
00157     bool ray_box_intersect( const CartVect& box_min,
00158                             const CartVect& box_max,
00159                             const CartVect& ray_pt,
00160                             const CartVect& ray_dir,
00161                             double& t_enter,
00162                             double& t_exit );
00163 
00164     /**\brief Test if plane intersects axis-aligned box
00165      *
00166      * Test for intersection between an unbounded plane and
00167      * an axis-aligned box.
00168      *\param plane_normal Vector in plane normal direction (need *not*
00169      *                    be a unit vector).  The N in
00170      *                    the plane equation: N . X + D = 0
00171      *\param plane_coeff  The scalar 'D' term in the plane equation:
00172      *                    N . X + D = 0
00173      *\param box_min_corner The smallest coordinates of the box along each
00174      *                    axis.  The corner of the box for which all three
00175      *                    coordinate values are smaller than those of any
00176      *                    other corner.  The X, Y, Z values for the planes
00177      *                    normal to those axes and bounding the box on the
00178      *                    -X, -Y, and -Z sides respectively.
00179      *\param box_max_corner The largest coordinates of the box along each
00180      *                    axis.  The corner of the box for which all three
00181      *                    coordinate values are larger than those of any
00182      *                    other corner.  The X, Y, Z values for the planes
00183      *                    normal to those axes and bounding the box on the
00184      *                    +X, +Y, and +Z sides respectively.
00185      *\return true if overlap, false otherwise.
00186      */
00187     bool box_plane_overlap( const CartVect& plane_normal,
00188                             double plane_coeff,
00189                             CartVect box_min_corner,
00190                             CartVect box_max_corner );
00191 
00192     /**\brief Test if triangle intersects axis-aligned box
00193      *
00194      * Test if a triangle intersects an axis-aligned box.
00195      *\param triangle_corners  The corners of the triangle.
00196      *\param box_min_corner The smallest coordinates of the box along each
00197      *                    axis.  The corner of the box for which all three
00198      *                    coordinate values are smaller than those of any
00199      *                    other corner.  The X, Y, Z values for the planes
00200      *                    normal to those axes and bounding the box on the
00201      *                    -X, -Y, and -Z sides respectively.
00202      *\param box_max_corner The largest coordinates of the box along each
00203      *                    axis.  The corner of the box for which all three
00204      *                    coordinate values are larger than those of any
00205      *                    other corner.  The X, Y, Z values for the planes
00206      *                    normal to those axes and bounding the box on the
00207      *                    +X, +Y, and +Z sides respectively.
00208      *\param tolerance    The tolerance used in the intersection test.  The box
00209      *                    size is increased by this amount before the intersection
00210      *                    test.
00211      *\return true if overlap, false otherwise.
00212      */
00213     bool box_tri_overlap( const CartVect triangle_corners[3],
00214                           const CartVect& box_min_corner,
00215                           const CartVect& box_max_corner,
00216                           double tolerance );
00217 
00218     /**\brief Test if triangle intersects axis-aligned box
00219      *
00220      * Test if a triangle intersects an axis-aligned box.
00221      *\param triangle_corners  The corners of the triangle.
00222      *\param box_center   The center of the box.
00223      *\param box_hanf_dims The distance along each axis, respectively, from the
00224      *                    box_center to the boundary of the box.
00225      *\return true if overlap, false otherwise.
00226      */
00227     bool box_tri_overlap( const CartVect triangle_corners[3],
00228                           const CartVect& box_center,
00229                           const CartVect& box_half_dims );
00230 
00231     bool box_point_overlap( const CartVect& box_min_corner,
00232                             const CartVect& box_max_corner,
00233                             const CartVect& point,
00234                             double tolerance );
00235 
00236     /**\brief Test if the specified element intersects an axis-aligned box.
00237      *
00238      * Test if element intersects axis-aligned box.  Use element-specific
00239      * optimization if available, otherwise call box_general_elem_overlap.
00240      *
00241      *\param elem_corners The coordinates of the element vertices
00242      *\param elem_type    The toplogy of the element.
00243      *\param box_center   The center of the axis-aligned box
00244      *\param box_half_dims Half of the width of the box in each axial
00245      *                     direction.
00246      */
00247     bool box_elem_overlap( const CartVect* elem_corners,
00248                            EntityType elem_type,
00249                            const CartVect& box_center,
00250                            const CartVect& box_half_dims,
00251                            int nodecount = 0 );
00252 
00253     /**\brief Test if the specified element intersects an axis-aligned box.
00254      *
00255      * Uses MBCN and separating axis theorem for general algorithm that
00256      * works for all fixed-size elements (not poly*).
00257      *
00258      *\param elem_corners The coordinates of the element vertices
00259      *\param elem_type    The toplogy of the element.
00260      *\param box_center   The center of the axis-aligned box
00261      *\param box_half_dims Half of the width of the box in each axial
00262      *                     direction.
00263      */
00264     bool box_linear_elem_overlap( const CartVect* elem_corners,
00265                                   EntityType elem_type,
00266                                   const CartVect& box_center,
00267                                   const CartVect& box_half_dims );
00268 
00269     /**\brief Test if the specified element intersects an axis-aligned box.
00270      *
00271      * Uses MBCN and separating axis theorem for general algorithm that
00272      * works for all fixed-size elements (not poly*).  Box and element
00273      * vertices must be translated such that box center is at origin.
00274      *
00275      *\param elem_corners The coordinates of the element vertices, in
00276      *                    local coordinate system of box.
00277      *\param elem_type    The toplogy of the element.
00278      *\param box_half_dims Half of the width of the box in each axial
00279      *                     direction.
00280      */
00281     bool box_linear_elem_overlap( const CartVect* elem_corners, EntityType elem_type, const CartVect& box_half_dims );
00282 
00283     void closest_location_on_box( const CartVect& box_min_corner,
00284                                   const CartVect& box_max_corner,
00285                                   const CartVect& point,
00286                                   CartVect& closest );
00287 
00288     /**\brief find closest location on triangle
00289      *
00290      * Find closest location on linear triangle.
00291      *\param location  Input position to evaluate from
00292      *\param vertices  Array of three corner vertex coordinates.
00293      *\param closest_out Result position
00294      */
00295     void closest_location_on_tri( const CartVect& location, const CartVect* vertices, CartVect& closest_out );
00296 
00297     /**\brief find closest location on polygon
00298      *
00299      * Find closest location on polygon
00300      *\param location  Input position to evaluate from
00301      *\param vertices  Array of corner vertex coordinates.
00302      *\param num_vertices Length of 'vertices' array.
00303      *\param closest_out Result position
00304      */
00305     void closest_location_on_polygon( const CartVect& location,
00306                                       const CartVect* vertices,
00307                                       int num_vertices,
00308                                       CartVect& closest_out );
00309 
00310     /**\brief find closest topological location on triangle
00311      *
00312      * Find closest location on linear triangle.
00313      *\param location  Input position to evaluate from
00314      *\param vertices  Array of three corner vertex coordinates.
00315      *\param tolerance Tolerance to use when comparing to corners and edges
00316      *\param closest_out Result position
00317      *\param closest_topo Closest topological entity
00318      *                     0-2 : vertex index
00319      *                     3-5 : edge beginning at closest_topo - 3
00320      *                       6 : triangle interior
00321      */
00322     void closest_location_on_tri( const CartVect& location,
00323                                   const CartVect* vertices,
00324                                   double tolerance,
00325                                   CartVect& closest_out,
00326                                   int& closest_topo );
00327 
00328     // Finds whether or not a box defined by the center and the half
00329     // width intersects a trilinear hex defined by its eight vertices.
00330     bool box_hex_overlap( const CartVect hexv[8], const CartVect& box_center, const CartVect& box_dims );
00331 
00332     // Finds whether or not a box defined by the center and the half
00333     // width intersects a linear tetrahedron defined by its four vertices.
00334     bool box_tet_overlap( const CartVect tet_corners[4], const CartVect& box_center, const CartVect& box_dims );
00335 
00336     // tests if 2 boxes overlap within a tolerance
00337     // assume that data is valid, box_min1 << box_max1, and box_min2<< box_max2
00338     // they are overlapping if they are overlapping in all directions
00339     // for example in x direction:
00340     //   box_max2[0] + tolerance < box_min1[0] -- false
00341     bool boxes_overlap( const CartVect& box_min1,
00342                         const CartVect& box_max1,
00343                         const CartVect& box_min2,
00344                         const CartVect& box_max2,
00345                         double tolerance );
00346 
00347     // see if boxes formed by 2 lists of "CartVect"s overlap
00348     bool bounding_boxes_overlap( const CartVect* list1, int num1, const CartVect* list2, int num2, double tolerance );
00349 
00350     // see if boxes from vertices in 2d overlap (used in gnomonic planes right now)
00351     bool bounding_boxes_overlap_2d( const double* list1, int num1, const double* list2, int num2, double tolerance );
00352     // point_in_trilinear_hex
00353     // Tests if a point in xyz space is within a hex element defined with
00354     // its eight vertex points forming a trilinear basis function.  Computes
00355     // the natural coordinates with respect to the hex of the xyz point
00356     // and checks if each are between +/-1.  If anyone is outside the range
00357     // the function returns false, otherwise it returns true.
00358     //
00359     bool point_in_trilinear_hex( const CartVect* hex, const CartVect& xyz, double etol );
00360 
00361     bool nat_coords_trilinear_hex( const CartVect* corner_coords, const CartVect& x, CartVect& xi, double tol );
00362 }  // namespace GeomUtil
00363 
00364 }  // namespace moab
00365 
00366 #endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines