Mesh Oriented datABase  (version 5.5.1)
An array-based unstructured mesh library
moab::AdaptiveKDTree Class Reference

Adaptive KD tree, for sorting and searching entities spatially. More...

#include <AdaptiveKDTree.hpp>

+ Inheritance diagram for moab::AdaptiveKDTree:
+ Collaboration diagram for moab::AdaptiveKDTree:

Classes

struct  Plane
 Split plane. More...
 

Public Types

enum  Axis { X = 0 , Y = 1 , Z = 2 }
 Enumeriate split plane directions. More...
 
enum  CandidatePlaneSet { SUBDIVISION = 0 , SUBDIVISION_SNAP , VERTEX_MEDIAN , VERTEX_SAMPLE }
 methods for selecting candidate split planes More...
 

Public Member Functions

 AdaptiveKDTree (Interface *iface)
 
 AdaptiveKDTree (Interface *iface, const Range &entities, EntityHandle *tree_root_set=NULL, FileOptions *opts=NULL)
 Constructor (build the tree on construction) Construct a tree object, and build the tree with entities input. See comments for build_tree() for detailed description of arguments. More...
 
 ~AdaptiveKDTree ()
 
ErrorCode parse_options (FileOptions &options)
 Parse options for tree creation. More...
 
virtual ErrorCode build_tree (const Range &entities, EntityHandle *tree_root_set=NULL, FileOptions *options=NULL)
 
virtual ErrorCode reset_tree ()
 Reset the tree, optionally checking we have the right root. More...
 
virtual ErrorCode point_search (const double *point, EntityHandle &leaf_out, const double iter_tol=1.0e-10, const double inside_tol=1.0e-6, bool *multiple_leaves=NULL, EntityHandle *start_node=NULL, CartVect *params=NULL)
 Get leaf containing input position. More...
 
ErrorCode point_search (const double *point, AdaptiveKDTreeIter &leaf_it, const double iter_tol=1.0e-10, const double inside_tol=1.0e-6, bool *multiple_leaves=NULL, EntityHandle *start_node=NULL)
 Get leaf containing input position. More...
 
virtual ErrorCode distance_search (const double *point, const double distance, std::vector< EntityHandle > &leaves_out, const double iter_tol=1.0e-10, const double inside_tol=1.0e-6, std::vector< double > *dists_out=NULL, std::vector< CartVect > *params_out=NULL, EntityHandle *start_node=NULL)
 Find all leaves within a given distance from point If dists_out input non-NULL, also returns distances from each leaf; if point i is inside leaf, 0 is given as dists_out[i]. If params_out is non-NULL and myEval is non-NULL, will evaluate individual entities in tree nodes and return containing entities in leaves_out. In those cases, if params_out is also non-NULL, will return parameters in those elements in that vector. More...
 
ErrorCode get_info (EntityHandle root, double min[3], double max[3], unsigned int &dep)
 Return some basic information about the tree Stats are returned for tree starting from input node or tree root (root = 0) More...
 
ErrorCode get_split_plane (EntityHandle node, Plane &plane)
 Get split plane for tree node. More...
 
ErrorCode set_split_plane (EntityHandle node, const Plane &plane)
 Set split plane for tree node. More...
 
ErrorCode get_tree_iterator (EntityHandle tree_root, AdaptiveKDTreeIter &result)
 Get iterator for tree. More...
 
ErrorCode get_last_iterator (EntityHandle tree_root, AdaptiveKDTreeIter &result)
 Get iterator at right-most ('last') leaf. More...
 
ErrorCode get_sub_tree_iterator (EntityHandle tree_root, const double box_min[3], const double box_max[3], AdaptiveKDTreeIter &result)
 Get iterator for tree or subtree. More...
 
ErrorCode split_leaf (AdaptiveKDTreeIter &leaf, Plane plane)
 Split leaf of tree Updates iterator location to point to first new leaf node. More...
 
ErrorCode split_leaf (AdaptiveKDTreeIter &leaf, Plane plane, EntityHandle &left_child, EntityHandle &right_child)
 Split leaf of tree Updates iterator location to point to first new leaf node. More...
 
ErrorCode split_leaf (AdaptiveKDTreeIter &leaf, Plane plane, const Range &left_entities, const Range &right_entities)
 Split leaf of tree Updates iterator location to point to first new leaf node. More...
 
ErrorCode split_leaf (AdaptiveKDTreeIter &leaf, Plane plane, const std::vector< EntityHandle > &left_entities, const std::vector< EntityHandle > &right_entities)
 Split leaf of tree Updates iterator location to point to first new leaf node. More...
 
ErrorCode merge_leaf (AdaptiveKDTreeIter &iter)
 Merge the leaf pointed to by the current iterator with it's sibling. If the sibling is not a leaf, multiple merges may be done. More...
 
ErrorCode closest_triangle (EntityHandle tree_root, const double from_coords[3], double closest_point_out[3], EntityHandle &triangle_out)
 Find triangle closest to input position. More...
 
ErrorCode sphere_intersect_triangles (EntityHandle tree_root, const double center[3], double radius, std::vector< EntityHandle > &triangles)
 
ErrorCode ray_intersect_triangles (EntityHandle tree_root, const double tolerance, const double ray_unit_dir[3], const double ray_base_pt[3], std::vector< EntityHandle > &triangles_out, std::vector< double > &distance_out, int result_count_limit=0, double distance_limit=-1.0)
 
ErrorCode compute_depth (EntityHandle root, unsigned int &min_depth, unsigned int &max_depth)
 
virtual ErrorCode print ()
 print various things about this tree More...
 
- Public Member Functions inherited from moab::Tree
 Tree (Interface *iface)
 Constructor (bare) More...
 
virtual ~Tree ()
 Destructor. More...
 
ErrorCode delete_tree_sets ()
 Delete the entity sets associated with the tree, starting with the root and traversing children. More...
 
virtual ErrorCode get_bounding_box (BoundBox &box, EntityHandle *tree_node=NULL) const
 Get bounding box for tree below tree_node, or entire tree If no tree has been built yet, returns +/- DBL_MAX for all dimensions. Note for some tree types, boxes are not available for non-root nodes, and this function will return failure if non-root is passed in. More...
 
ErrorCode find_all_trees (Range &results)
 Find all trees, by bounding box tag. More...
 
Interfacemoab ()
 Return the MOAB interface associated with this tree. More...
 
const Interfacemoab () const
 Return the MOAB interface associated with this tree. More...
 
double get_max_depth ()
 Get max depth set on tree. More...
 
double get_max_per_leaf ()
 Get max entities per leaf set on tree. More...
 
TreeStatstree_stats ()
 Get tree traversal stats object. More...
 
const TreeStatstree_stats () const
 Get tree traversal stats object. More...
 
ErrorCode create_root (const double box_min[3], const double box_max[3], EntityHandle &root_handle)
 Create tree root and tag with bounding box. More...
 
ElemEvaluatorget_eval ()
 get/set the ElemEvaluator More...
 
void set_eval (ElemEvaluator *eval)
 get/set the ElemEvaluator More...
 

Private Member Functions

ErrorCode init ()
 
ErrorCode find_close_triangle (EntityHandle root, const double from_point[3], double pt[3], EntityHandle &triangle)
 find a triangle near the input point More...
 
ErrorCode make_tag (Interface *iface, std::string name, TagType storage, DataType type, int count, void *default_val, Tag &tag_handle, std::vector< Tag > &created_tags)
 
ErrorCode intersect_children_with_elems (const Range &elems, AdaptiveKDTree::Plane plane, double eps, CartVect box_min, CartVect box_max, Range &left_tris, Range &right_tris, Range &both_tris, double &metric_value)
 
ErrorCode best_subdivision_snap_plane (int num_planes, const AdaptiveKDTreeIter &iter, Range &best_left, Range &best_right, Range &best_both, AdaptiveKDTree::Plane &best_plane, std::vector< double > &tmp_data, double eps)
 
ErrorCode best_subdivision_plane (int num_planes, const AdaptiveKDTreeIter &iter, Range &best_left, Range &best_right, Range &best_both, AdaptiveKDTree::Plane &best_plane, double eps)
 
ErrorCode best_vertex_median_plane (int num_planes, const AdaptiveKDTreeIter &iter, Range &best_left, Range &best_right, Range &best_both, AdaptiveKDTree::Plane &best_plane, std::vector< double > &coords, double eps)
 
ErrorCode best_vertex_sample_plane (int num_planes, const AdaptiveKDTreeIter &iter, Range &best_left, Range &best_right, Range &best_both, AdaptiveKDTree::Plane &best_plane, std::vector< double > &coords, std::vector< EntityHandle > &indices, double eps)
 

Private Attributes

Tag planeTag
 
Tag axisTag
 
unsigned splitsPerDir
 
CandidatePlaneSet planeSet
 
bool spherical
 
double radius
 

Static Private Attributes

static const char * treeName = "AKDTree"
 

Friends

class AdaptiveKDTreeIter
 

Additional Inherited Members

- Protected Member Functions inherited from moab::Tree
ErrorCode parse_common_options (FileOptions &options)
 Parse options common to all trees. More...
 
Tag get_box_tag (bool create_if_missing=true)
 Get the box tag, possibly constructing it first. More...
 
- Protected Attributes inherited from moab::Tree
InterfacembImpl
 
BoundBox boundBox
 
int maxPerLeaf
 
int maxDepth
 
int treeDepth
 
double minWidth
 
unsigned int meshsetFlags
 
bool cleanUp
 
EntityHandle myRoot
 
Tag boxTag
 
std::string boxTagName
 
TreeStats treeStats
 
ElemEvaluatormyEval
 

Detailed Description

Adaptive KD tree, for sorting and searching entities spatially.

Definition at line 23 of file AdaptiveKDTree.hpp.

Member Enumeration Documentation

◆ Axis

Enumeriate split plane directions.

Enumerator

Definition at line 148 of file AdaptiveKDTree.hpp.

149  {
150  X = 0,
151  Y = 1,
152  Z = 2
153  };

◆ CandidatePlaneSet

methods for selecting candidate split planes

Enumerator
SUBDIVISION 

Candidiate planes at evenly spaced intervals.

SUBDIVISION_SNAP 

Like SUBDIVISION, except snap to closest vertex coordinate.

VERTEX_MEDIAN 

Median vertex coodinate values.

VERTEX_SAMPLE 

Random sampling of vertex coordinate values.

Definition at line 248 of file AdaptiveKDTree.hpp.

249  {
250  //! Candidiate planes at evenly spaced intervals
251  SUBDIVISION = 0,
252  //! Like SUBDIVISION, except snap to closest vertex coordinate
253  SUBDIVISION_SNAP, // = 1
254  //! Median vertex coodinate values
255  VERTEX_MEDIAN, // = 2
256  //! Random sampling of vertex coordinate values
257  VERTEX_SAMPLE // = 3
258  };

Constructor & Destructor Documentation

◆ AdaptiveKDTree() [1/2]

moab::AdaptiveKDTree::AdaptiveKDTree ( Interface iface)

Definition at line 36 of file AdaptiveKDTree.cpp.

37  : Tree( iface ), planeTag( 0 ), axisTag( 0 ), splitsPerDir( 3 ), planeSet( SUBDIVISION_SNAP ), spherical( false ),
38  radius( 1.0 )
39 {
41 
42  ErrorCode rval = init();
43  if( MB_SUCCESS != rval ) throw rval;
44 }

References moab::Tree::boxTagName, ErrorCode, init(), MB_SUCCESS, and treeName.

◆ AdaptiveKDTree() [2/2]

moab::AdaptiveKDTree::AdaptiveKDTree ( Interface iface,
const Range entities,
EntityHandle tree_root_set = NULL,
FileOptions opts = NULL 
)

Constructor (build the tree on construction) Construct a tree object, and build the tree with entities input. See comments for build_tree() for detailed description of arguments.

Parameters
ifaceMOAB instance
entitiesEntities to build tree around
tree_rootRoot set for tree (see function description)
optsOptions for tree (see function description)

Definition at line 46 of file AdaptiveKDTree.cpp.

50  : Tree( iface ), planeTag( 0 ), axisTag( 0 ), splitsPerDir( 3 ), planeSet( SUBDIVISION_SNAP ), spherical( false ),
51  radius( 1.0 )
52 {
54 
55  ErrorCode rval;
56  if( opts )
57  {
58  rval = parse_options( *opts );
59  if( MB_SUCCESS != rval ) throw rval;
60  }
61 
62  rval = init();
63  if( MB_SUCCESS != rval ) throw rval;
64 
65  rval = build_tree( entities, tree_root_set, opts );
66  if( MB_SUCCESS != rval ) throw rval;
67 }

References moab::Tree::boxTagName, build_tree(), entities, ErrorCode, init(), MB_SUCCESS, parse_options(), and treeName.

◆ ~AdaptiveKDTree()

moab::AdaptiveKDTree::~AdaptiveKDTree ( )

Definition at line 69 of file AdaptiveKDTree.cpp.

70 {
71  if( !cleanUp ) return;
72 
73  if( myRoot )
74  {
75  reset_tree();
76  myRoot = 0;
77  }
78 }

References moab::Tree::cleanUp, moab::Tree::myRoot, and reset_tree().

Member Function Documentation

◆ best_subdivision_plane()

ErrorCode moab::AdaptiveKDTree::best_subdivision_plane ( int  num_planes,
const AdaptiveKDTreeIter iter,
Range best_left,
Range best_right,
Range best_both,
AdaptiveKDTree::Plane best_plane,
double  eps 
)
private

Definition at line 942 of file AdaptiveKDTree.cpp.

949 {
950  double metric_val = std::numeric_limits< unsigned >::max();
951 
952  ErrorCode r;
953  const CartVect box_min( iter.box_min() );
954  const CartVect box_max( iter.box_max() );
955  const CartVect diff( box_max - box_min );
956 
957  Range entities;
958  r = iter.tool()->moab()->get_entities_by_handle( iter.handle(), entities );
959  if( MB_SUCCESS != r ) return r;
960  const size_t p_count = entities.size();
961 
962  for( int axis = 0; axis < 3; ++axis )
963  {
964  int plane_count = num_planes;
965  if( ( num_planes + 1 ) * eps >= diff[axis] ) plane_count = (int)( diff[axis] / eps ) - 1;
966 
967  for( int p = 1; p <= plane_count; ++p )
968  {
969  AdaptiveKDTree::Plane plane = { box_min[axis] + ( p / ( 1.0 + plane_count ) ) * diff[axis], axis };
970  Range left, right, both;
971  double val;
972  r = intersect_children_with_elems( entities, plane, eps, box_min, box_max, left, right, both, val );
973  if( MB_SUCCESS != r ) return r;
974  const size_t sdiff = p_count - both.size();
975  if( left.size() == sdiff || right.size() == sdiff ) continue;
976 
977  if( val >= metric_val ) continue;
978 
979  metric_val = val;
980  best_plane = plane;
981  best_left.swap( left );
982  best_right.swap( right );
983  best_both.swap( both );
984  }
985  }
986 
987  return MB_SUCCESS;
988 }

References moab::AdaptiveKDTreeIter::box_max(), box_max(), moab::AdaptiveKDTreeIter::box_min(), box_min(), entities, ErrorCode, moab::Interface::get_entities_by_handle(), moab::AdaptiveKDTreeIter::handle(), intersect_children_with_elems(), MB_SUCCESS, moab::Tree::moab(), moab::Range::size(), moab::Range::swap(), and moab::AdaptiveKDTreeIter::tool().

Referenced by build_tree().

◆ best_subdivision_snap_plane()

ErrorCode moab::AdaptiveKDTree::best_subdivision_snap_plane ( int  num_planes,
const AdaptiveKDTreeIter iter,
Range best_left,
Range best_right,
Range best_both,
AdaptiveKDTree::Plane best_plane,
std::vector< double > &  tmp_data,
double  eps 
)
private

Definition at line 990 of file AdaptiveKDTree.cpp.

998 {
999  double metric_val = std::numeric_limits< unsigned >::max();
1000 
1001  ErrorCode r;
1002  // const CartVect tol(eps*diff);
1003 
1004  Range entities, vertices;
1005  r = iter.tool()->moab()->get_entities_by_handle( iter.handle(), entities );
1006  if( MB_SUCCESS != r ) return r;
1007  const size_t p_count = entities.size();
1008  r = iter.tool()->moab()->get_adjacencies( entities, 0, false, vertices, Interface::UNION );
1009  if( MB_SUCCESS != r ) return r;
1010 
1011  unsigned int nverts = vertices.size();
1012  tmp_data.resize( 3 * nverts );
1013  r = iter.tool()->moab()->get_coords( vertices, &tmp_data[0], &tmp_data[nverts], &tmp_data[2 * nverts] );
1014  if( MB_SUCCESS != r ) return r;
1015 
1016  // calculate bounding box of vertices
1017  // decide based on the actual box the splitting plane
1018  // do not decide based on iterator box.
1019  // it could be too big
1020  // BoundBox box;
1021  // r = box.update(*moab(), vertices);
1022  CartVect box_min;
1023  CartVect box_max;
1024  for( int dir = 0; dir < 3; dir++ )
1025  {
1026  double amin = tmp_data[dir * nverts];
1027  double amax = amin;
1028  double* p = &tmp_data[dir * nverts + 1];
1029  for( unsigned int i = 1; i < nverts; i++ )
1030  {
1031  if( *p < amin ) amin = *p;
1032  if( *p > amax ) amax = *p;
1033  p++;
1034  }
1035  box_min[dir] = amin;
1036  box_max[dir] = amax;
1037  }
1038  CartVect diff( box_max - box_min );
1039 
1040  for( int axis = 0; axis < 3; ++axis )
1041  {
1042  int plane_count = num_planes;
1043 
1044  // if num_planes results in width < eps, reset the plane count
1045  if( ( num_planes + 1 ) * eps >= diff[axis] ) plane_count = (int)( diff[axis] / eps ) - 1;
1046 
1047  for( int p = 1; p <= plane_count; ++p )
1048  {
1049 
1050  // coord of this plane on axis
1051  double coord = box_min[axis] + ( p / ( 1.0 + plane_count ) ) * diff[axis];
1052 
1053  // find closest vertex coordinate to this plane position
1054  unsigned int istrt = axis * nverts;
1055  double closest_coord = tmp_data[istrt];
1056  for( unsigned i = 1; i < nverts; ++i )
1057  if( fabs( coord - tmp_data[istrt + i] ) < fabs( coord - closest_coord ) )
1058  closest_coord = tmp_data[istrt + i];
1059  if( closest_coord - box_min[axis] <= eps || box_max[axis] - closest_coord <= eps ) continue;
1060 
1061  // seprate elems into left/right/both, and compute separating metric
1062  AdaptiveKDTree::Plane plane = { closest_coord, axis };
1063  Range left, right, both;
1064  double val;
1065  r = intersect_children_with_elems( entities, plane, eps, box_min, box_max, left, right, both, val );
1066  if( MB_SUCCESS != r ) return r;
1067  const size_t d = p_count - both.size();
1068  if( left.size() == d || right.size() == d ) continue;
1069 
1070  if( val >= metric_val ) continue;
1071 
1072  metric_val = val;
1073  best_plane = plane;
1074  best_left.swap( left );
1075  best_right.swap( right );
1076  best_both.swap( both );
1077  }
1078  }
1079 
1080  return MB_SUCCESS;
1081 }

References box_max(), box_min(), entities, ErrorCode, moab::Interface::get_adjacencies(), moab::Interface::get_coords(), moab::Interface::get_entities_by_handle(), moab::AdaptiveKDTreeIter::handle(), intersect_children_with_elems(), MB_SUCCESS, moab::Tree::moab(), moab::Range::size(), moab::Range::swap(), moab::AdaptiveKDTreeIter::tool(), and moab::Interface::UNION.

Referenced by build_tree().

◆ best_vertex_median_plane()

ErrorCode moab::AdaptiveKDTree::best_vertex_median_plane ( int  num_planes,
const AdaptiveKDTreeIter iter,
Range best_left,
Range best_right,
Range best_both,
AdaptiveKDTree::Plane best_plane,
std::vector< double > &  coords,
double  eps 
)
private

Definition at line 1083 of file AdaptiveKDTree.cpp.

1091 {
1092  double metric_val = std::numeric_limits< unsigned >::max();
1093 
1094  ErrorCode r;
1095  const CartVect box_min( iter.box_min() );
1096  const CartVect box_max( iter.box_max() );
1097 
1098  Range entities, vertices;
1099  r = iter.tool()->moab()->get_entities_by_handle( iter.handle(), entities );
1100  if( MB_SUCCESS != r ) return r;
1101  const size_t p_count = entities.size();
1102  r = iter.tool()->moab()->get_adjacencies( entities, 0, false, vertices, Interface::UNION );
1103  if( MB_SUCCESS != r ) return r;
1104 
1105  coords.resize( vertices.size() );
1106  for( int axis = 0; axis < 3; ++axis )
1107  {
1108  if( box_max[axis] - box_min[axis] <= 2 * eps ) continue;
1109 
1110  double* ptrs[] = { 0, 0, 0 };
1111  ptrs[axis] = &coords[0];
1112  r = iter.tool()->moab()->get_coords( vertices, ptrs[0], ptrs[1], ptrs[2] );
1113  if( MB_SUCCESS != r ) return r;
1114 
1115  std::sort( coords.begin(), coords.end() );
1116  std::vector< double >::iterator citer;
1117  citer = std::upper_bound( coords.begin(), coords.end(), box_min[axis] + eps );
1118  const size_t count = std::upper_bound( citer, coords.end(), box_max[axis] - eps ) - citer;
1119  size_t step;
1120  int np = num_planes;
1121  if( count < 2 * (size_t)num_planes )
1122  {
1123  step = 1;
1124  np = count - 1;
1125  }
1126  else
1127  {
1128  step = count / ( num_planes + 1 );
1129  }
1130 
1131  for( int p = 1; p <= np; ++p )
1132  {
1133 
1134  citer += step;
1135  AdaptiveKDTree::Plane plane = { *citer, axis };
1136  Range left, right, both;
1137  double val;
1138  r = intersect_children_with_elems( entities, plane, eps, box_min, box_max, left, right, both, val );
1139  if( MB_SUCCESS != r ) return r;
1140  const size_t diff = p_count - both.size();
1141  if( left.size() == diff || right.size() == diff ) continue;
1142 
1143  if( val >= metric_val ) continue;
1144 
1145  metric_val = val;
1146  best_plane = plane;
1147  best_left.swap( left );
1148  best_right.swap( right );
1149  best_both.swap( both );
1150  }
1151  }
1152 
1153  return MB_SUCCESS;
1154 }

References moab::AdaptiveKDTreeIter::box_max(), box_max(), moab::AdaptiveKDTreeIter::box_min(), box_min(), entities, ErrorCode, moab::Interface::get_adjacencies(), moab::Interface::get_coords(), moab::Interface::get_entities_by_handle(), moab::AdaptiveKDTreeIter::handle(), intersect_children_with_elems(), MB_SUCCESS, moab::Tree::moab(), moab::Range::size(), moab::Range::swap(), moab::AdaptiveKDTreeIter::tool(), and moab::Interface::UNION.

Referenced by build_tree().

◆ best_vertex_sample_plane()

ErrorCode moab::AdaptiveKDTree::best_vertex_sample_plane ( int  num_planes,
const AdaptiveKDTreeIter iter,
Range best_left,
Range best_right,
Range best_both,
AdaptiveKDTree::Plane best_plane,
std::vector< double > &  coords,
std::vector< EntityHandle > &  indices,
double  eps 
)
private

Definition at line 1156 of file AdaptiveKDTree.cpp.

1165 {
1166  const size_t random_elem_threshold = 20 * num_planes;
1167  double metric_val = std::numeric_limits< unsigned >::max();
1168 
1169  ErrorCode r;
1170  const CartVect box_min( iter.box_min() );
1171  const CartVect box_max( iter.box_max() );
1172 
1173  Range entities, vertices;
1174  r = iter.tool()->moab()->get_entities_by_handle( iter.handle(), entities );
1175  if( MB_SUCCESS != r ) return r;
1176 
1177  // We are selecting random vertex coordinates to use for candidate split
1178  // planes. So if element list is large, begin by selecting random elements.
1179  const size_t p_count = entities.size();
1180  coords.resize( 3 * num_planes );
1181  if( p_count < random_elem_threshold )
1182  {
1183  r = iter.tool()->moab()->get_adjacencies( entities, 0, false, vertices, Interface::UNION );
1184  if( MB_SUCCESS != r ) return r;
1185  }
1186  else
1187  {
1188  indices.resize( random_elem_threshold );
1189  const int num_rand = p_count / RAND_MAX + 1;
1190  for( size_t j = 0; j < random_elem_threshold; ++j )
1191  {
1192  size_t rnd = rand();
1193  for( int i = num_rand; i > 1; --i )
1194  rnd *= rand();
1195  rnd %= p_count;
1196  indices[j] = entities[rnd];
1197  }
1198  r = iter.tool()->moab()->get_adjacencies( &indices[0], random_elem_threshold, 0, false, vertices,
1199  Interface::UNION );
1200  if( MB_SUCCESS != r ) return r;
1201  }
1202 
1203  coords.resize( vertices.size() );
1204  for( int axis = 0; axis < 3; ++axis )
1205  {
1206  if( box_max[axis] - box_min[axis] <= 2 * eps ) continue;
1207 
1208  double* ptrs[] = { 0, 0, 0 };
1209  ptrs[axis] = &coords[0];
1210  r = iter.tool()->moab()->get_coords( vertices, ptrs[0], ptrs[1], ptrs[2] );
1211  if( MB_SUCCESS != r ) return r;
1212 
1213  size_t num_valid_coords = 0;
1214  for( size_t i = 0; i < coords.size(); ++i )
1215  if( coords[i] > box_min[axis] + eps && coords[i] < box_max[axis] - eps ) ++num_valid_coords;
1216 
1217  if( 2 * (size_t)num_planes > num_valid_coords )
1218  {
1219  indices.clear();
1220  for( size_t i = 0; i < coords.size(); ++i )
1221  if( coords[i] > box_min[axis] + eps && coords[i] < box_max[axis] - eps ) indices.push_back( i );
1222  }
1223  else
1224  {
1225  indices.resize( num_planes );
1226  // make sure random indices are sufficient to cover entire range
1227  const int num_rand = coords.size() / RAND_MAX + 1;
1228  for( int j = 0; j < num_planes; ++j )
1229  {
1230  size_t rnd;
1231  do
1232  {
1233  rnd = rand();
1234  for( int i = num_rand; i > 1; --i )
1235  rnd *= rand();
1236  rnd %= coords.size();
1237  } while( coords[rnd] <= box_min[axis] + eps || coords[rnd] >= box_max[axis] - eps );
1238  indices[j] = rnd;
1239  }
1240  }
1241 
1242  for( unsigned p = 0; p < indices.size(); ++p )
1243  {
1244 
1245  AdaptiveKDTree::Plane plane = { coords[indices[p]], axis };
1246  Range left, right, both;
1247  double val;
1248  r = intersect_children_with_elems( entities, plane, eps, box_min, box_max, left, right, both, val );
1249  if( MB_SUCCESS != r ) return r;
1250  const size_t diff = p_count - both.size();
1251  if( left.size() == diff || right.size() == diff ) continue;
1252 
1253  if( val >= metric_val ) continue;
1254 
1255  metric_val = val;
1256  best_plane = plane;
1257  best_left.swap( left );
1258  best_right.swap( right );
1259  best_both.swap( both );
1260  }
1261  }
1262 
1263  return MB_SUCCESS;
1264 }

References moab::AdaptiveKDTreeIter::box_max(), box_max(), moab::AdaptiveKDTreeIter::box_min(), box_min(), entities, ErrorCode, moab::Interface::get_adjacencies(), moab::Interface::get_coords(), moab::Interface::get_entities_by_handle(), moab::AdaptiveKDTreeIter::handle(), intersect_children_with_elems(), MB_SUCCESS, moab::Tree::moab(), moab::Range::size(), moab::Range::swap(), moab::AdaptiveKDTreeIter::tool(), and moab::Interface::UNION.

Referenced by build_tree().

◆ build_tree()

ErrorCode moab::AdaptiveKDTree::build_tree ( const Range entities,
EntityHandle tree_root_set = NULL,
FileOptions options = NULL 
)
virtual

Build the tree Build a tree with the entities input. If a non-NULL tree_root_set pointer is input, use the pointed-to set as the root of this tree (*tree_root_set!=0) otherwise construct a new root set and pass its handle back in *tree_root_set. Options vary by tree type; see Tree.hpp for common options; options specific to AdaptiveKDTree: SPLITS_PER_DIR: number of candidate splits considered per direction; default = 3 PLANE_SET: method used to decide split planes; see CandidatePlaneSet enum (below) for possible values; default = 1 (SUBDIVISION_SNAP)

Parameters
entitiesEntities with which to build the tree
tree_rootRoot set for tree (see function description)
optsOptions for tree (see function description)
Returns
Error is returned only on build failure

Implements moab::Tree.

Definition at line 80 of file AdaptiveKDTree.cpp.

81 {
82  ErrorCode rval;
83  CpuTimer cp;
84 
85  if( options )
86  {
87  rval = parse_options( *options );
88  if( MB_SUCCESS != rval ) return rval;
89 
90  if( !options->all_seen() ) return MB_FAILURE;
91  }
92 
93  // calculate bounding box of elements
94  BoundBox box;
95  rval = box.update( *moab(), entities, spherical, radius );
96  if( MB_SUCCESS != rval ) return rval;
97 
98  // create tree root
99  EntityHandle tmp_root;
100  if( !tree_root_set ) tree_root_set = &tmp_root;
101  rval = create_root( box.bMin.array(), box.bMax.array(), *tree_root_set );
102  if( MB_SUCCESS != rval ) return rval;
103  rval = moab()->add_entities( *tree_root_set, entities );
104  if( MB_SUCCESS != rval ) return rval;
105 
106  AdaptiveKDTreeIter iter;
107  iter.initialize( this, *tree_root_set, box.bMin.array(), box.bMax.array(), AdaptiveKDTreeIter::LEFT );
108 
109  std::vector< double > tmp_data;
110  std::vector< EntityHandle > tmp_data2;
111  for( ;; )
112  {
113 
114  int pcount;
115  rval = moab()->get_number_entities_by_handle( iter.handle(), pcount );
116  if( MB_SUCCESS != rval ) break;
117 
118  const size_t p_count = pcount;
119  Range best_left, best_right, best_both;
120  Plane best_plane = { HUGE_VAL, -1 };
121  if( (int)p_count > maxPerLeaf && (int)iter.depth() < maxDepth )
122  {
123  switch( planeSet )
124  {
126  rval = best_subdivision_plane( splitsPerDir, iter, best_left, best_right, best_both, best_plane,
127  minWidth );
128  break;
130  rval = best_subdivision_snap_plane( splitsPerDir, iter, best_left, best_right, best_both,
131  best_plane, tmp_data, minWidth );
132  break;
134  rval = best_vertex_median_plane( splitsPerDir, iter, best_left, best_right, best_both, best_plane,
135  tmp_data, minWidth );
136  break;
138  rval = best_vertex_sample_plane( splitsPerDir, iter, best_left, best_right, best_both, best_plane,
139  tmp_data, tmp_data2, minWidth );
140  break;
141  default:
142  rval = MB_FAILURE;
143  }
144 
145  if( MB_SUCCESS != rval ) return rval;
146  }
147 
148  if( best_plane.norm >= 0 )
149  {
150  best_left.merge( best_both );
151  best_right.merge( best_both );
152  rval = split_leaf( iter, best_plane, best_left, best_right );
153  if( MB_SUCCESS != rval ) return rval;
154  }
155  else
156  {
157  rval = iter.step();
158  if( MB_ENTITY_NOT_FOUND == rval )
159  {
161  treeStats.initTime = cp.time_elapsed();
162  return rval; // at end
163  }
164  else if( MB_SUCCESS != rval )
165  break;
166  }
167  }
168 
169  reset_tree();
170 
171  treeStats.reset();
172 
173  return rval;
174 }

References moab::Interface::add_entities(), moab::FileOptions::all_seen(), moab::CartVect::array(), best_subdivision_plane(), best_subdivision_snap_plane(), best_vertex_median_plane(), best_vertex_sample_plane(), moab::BoundBox::bMax, moab::BoundBox::bMin, moab::TreeStats::compute_stats(), moab::Tree::create_root(), moab::AdaptiveKDTreeIter::depth(), entities, ErrorCode, moab::Interface::get_number_entities_by_handle(), moab::AdaptiveKDTreeIter::handle(), moab::AdaptiveKDTreeIter::initialize(), moab::TreeStats::initTime, moab::AdaptiveKDTreeIter::LEFT, moab::Tree::maxDepth, moab::Tree::maxPerLeaf, MB_ENTITY_NOT_FOUND, MB_SUCCESS, moab::Tree::mbImpl, moab::Range::merge(), moab::Tree::minWidth, moab::Tree::moab(), moab::Tree::myRoot, moab::AdaptiveKDTree::Plane::norm, parse_options(), planeSet, radius, moab::TreeStats::reset(), reset_tree(), spherical, split_leaf(), splitsPerDir, moab::AdaptiveKDTreeIter::step(), SUBDIVISION, SUBDIVISION_SNAP, moab::CpuTimer::time_elapsed(), moab::Tree::treeStats, moab::BoundBox::update(), VERTEX_MEDIAN, and VERTEX_SAMPLE.

Referenced by AdaptiveKDTree(), moab::TempestRemapper::ConstructCoveringSet(), moab::Coupler::initialize_tree(), moab::Intx2Mesh::intersect_meshes(), moab::Intx2Mesh::intersect_meshes_kdtree(), main(), moab::MergeMesh::merge_all(), merge_duplicate_vertices(), moab::MergeMesh::merge_entities(), and moab::ReadNCDF::update().

◆ closest_triangle()

ErrorCode moab::AdaptiveKDTree::closest_triangle ( EntityHandle  tree_root,
const double  from_coords[3],
double  closest_point_out[3],
EntityHandle triangle_out 
)

Find triangle closest to input position.

Parameters
from_coordsThe input position to test against
closest_point_outThe closest point on the set of triangles in the tree
triangle_outThe triangle closest to the input position

Find the triangles in a set that are closer to the input position than any triangles in the 'closest_tris' list.

closest_tris is assumed to contain a list of triangles for which the first is the closest known triangle to the input position and the first entry in 'closest_pts' is the closest location on that triangle. Any other values in the lists must be other triangles for which the closest point is within the input tolerance of the closest closest point. This function will update the lists as appropriate if any closer triangles or triangles within the tolerance of the current closest location are found. The first entry is maintained as the closest of the list of triangles.

Definition at line 1763 of file AdaptiveKDTree.cpp.

1767 {
1768  ErrorCode rval;
1769  double shortest_dist_sqr = HUGE_VAL;
1770  std::vector< EntityHandle > leaves;
1771  const CartVect from( from_coords );
1772  CartVect closest_pt;
1773 
1774  // Find the leaf containing the input point
1775  // This search does not take into account any bounding box for the
1776  // tree, so it always returns one leaf.
1777  assert( tree_root );
1778  rval = find_close_triangle( tree_root, from_coords, closest_pt.array(), triangle_out );
1779  if( MB_SUCCESS != rval ) return rval;
1780 
1781  // Find any other leaves for which the bounding box is within
1782  // the same distance from the input point as the current closest
1783  // point is.
1784  CartVect diff = closest_pt - from;
1785  rval = distance_search( from_coords, sqrt( diff % diff ), leaves, 1.0e-10, 1.0e-6, NULL, NULL, &tree_root );
1786  if( MB_SUCCESS != rval ) return rval;
1787 
1788  // Check any close leaves to see if they contain triangles that
1789  // are as close to or closer than the current closest triangle(s).
1790  for( unsigned i = 0; i < leaves.size(); ++i )
1791  {
1792  rval = closest_to_triangles( moab(), leaves[i], from, shortest_dist_sqr, closest_pt, triangle_out );
1793  if( MB_SUCCESS != rval ) return rval;
1794  }
1795 
1796  // pass back resulting position
1797  closest_pt.get( closest_point_out );
1798  return MB_SUCCESS;
1799 }

References moab::CartVect::array(), moab::closest_to_triangles(), distance_search(), ErrorCode, find_close_triangle(), moab::CartVect::get(), MB_SUCCESS, and moab::Tree::moab().

◆ compute_depth()

ErrorCode moab::AdaptiveKDTree::compute_depth ( EntityHandle  root,
unsigned int &  min_depth,
unsigned int &  max_depth 
)

Definition at line 2031 of file AdaptiveKDTree.cpp.

2032 {
2033  AdaptiveKDTreeIter iter;
2034  get_tree_iterator( root, iter );
2035  iter.step_to_first_leaf( AdaptiveKDTreeIter::LEFT );
2036  min_depth = max_depth = iter.depth();
2037 
2038  int num_of_elements = 0, max, min;
2039  moab()->get_number_entities_by_handle( iter.handle(), num_of_elements );
2040  max = min = num_of_elements;
2041  while( MB_SUCCESS == iter.step() )
2042  {
2043  int temp = 0;
2044  moab()->get_number_entities_by_handle( iter.handle(), temp );
2045  max = std::max( max, temp );
2046  min = std::min( min, temp );
2047  if( iter.depth() > max_depth )
2048  max_depth = iter.depth();
2049  else if( iter.depth() < min_depth )
2050  min_depth = iter.depth();
2051  }
2052  return MB_SUCCESS;
2053 }

References moab::AdaptiveKDTreeIter::depth(), moab::Interface::get_number_entities_by_handle(), get_tree_iterator(), moab::AdaptiveKDTreeIter::handle(), moab::AdaptiveKDTreeIter::LEFT, MB_SUCCESS, moab::Tree::moab(), moab::AdaptiveKDTreeIter::step(), and moab::AdaptiveKDTreeIter::step_to_first_leaf().

Referenced by get_info().

◆ distance_search()

ErrorCode moab::AdaptiveKDTree::distance_search ( const double *  point,
const double  distance,
std::vector< EntityHandle > &  leaves_out,
const double  iter_tol = 1.0e-10,
const double  inside_tol = 1.0e-6,
std::vector< double > *  dists_out = NULL,
std::vector< CartVect > *  params_out = NULL,
EntityHandle start_node = NULL 
)
virtual

Find all leaves within a given distance from point If dists_out input non-NULL, also returns distances from each leaf; if point i is inside leaf, 0 is given as dists_out[i]. If params_out is non-NULL and myEval is non-NULL, will evaluate individual entities in tree nodes and return containing entities in leaves_out. In those cases, if params_out is also non-NULL, will return parameters in those elements in that vector.

Parameters
pointPoint to be located in tree
distanceDistance within which to query
leaves_outLeaves within distance or containing point
iter_tolTolerance for convergence of point search
inside_tolTolerance for inside element calculation
dists_outIf non-NULL, will contain distsances to leaves
params_outIf non-NULL, will contain parameters of the point in the ents in leaves_out
start_nodeStart from this tree node (non-NULL) instead of tree root (NULL)

when we pop it from the list.)

when we pop it from the list.)

Implements moab::Tree.

Definition at line 1388 of file AdaptiveKDTree.cpp.

1396 {
1398  const double dist_sqr = distance * distance;
1399  const CartVect from( from_point );
1400  std::vector< NodeDistance > list,
1401  result_list_nodes; // list of subtrees to traverse, and results
1402  // pre-allocate space for default max tree depth
1403  list.reserve( maxDepth );
1404 
1405  // misc temporary values
1406  Plane plane;
1407  NodeDistance node;
1408  ErrorCode rval;
1409  std::vector< EntityHandle > children;
1410 
1411  // Get distance from input position to bounding box of tree
1412  // (zero if inside box)
1413  BoundBox box;
1414  rval = get_bounding_box( box );
1415  if( MB_SUCCESS == rval && !box.contains_point( from_point, iter_tol ) )
1416  {
1418  return MB_SUCCESS;
1419  }
1420 
1421  // if bounding box is not available (e.g. not starting from true root)
1422  // just start with zero. Less efficient, but will work.
1423  node.dist = CartVect( 0.0 );
1424  if( MB_SUCCESS == rval )
1425  {
1426  for( int i = 0; i < 3; ++i )
1427  {
1428  if( from_point[i] < box.bMin[i] )
1429  node.dist[i] = box.bMin[i] - from_point[i];
1430  else if( from_point[i] > box.bMax[i] )
1431  node.dist[i] = from_point[i] - box.bMax[i];
1432  }
1433  if( node.dist % node.dist > dist_sqr )
1434  {
1436  return MB_SUCCESS;
1437  }
1438  }
1439 
1440  // begin with root in list
1441  node.handle = ( tree_root ? *tree_root : myRoot );
1442  list.push_back( node );
1443 
1444  while( !list.empty() )
1445  {
1446 
1447  node = list.back();
1448  list.pop_back();
1450 
1451  // If leaf node, test contained triangles
1452  children.clear();
1453  rval = moab()->get_child_meshsets( node.handle, children );
1454  if( children.empty() )
1455  {
1457  if( myEval && result_params )
1458  {
1459  EntityHandle ent;
1460  CartVect params;
1461  rval = myEval->find_containing_entity( node.handle, from_point, iter_tol, inside_tol, ent,
1462  params.array(), &treeStats.traversalLeafObjectTests );
1463  if( MB_SUCCESS != rval )
1464  return rval;
1465  else if( ent )
1466  {
1467  result_list.push_back( ent );
1468  result_params->push_back( params );
1469  if( result_dists ) result_dists->push_back( 0.0 );
1470  }
1471  }
1472  else
1473  {
1474  result_list_nodes.push_back( node );
1475  continue;
1476  }
1477  }
1478 
1479  // If not leaf node, add children to working list
1480  rval = get_split_plane( node.handle, plane );
1481  if( MB_SUCCESS != rval ) return rval;
1482 
1483  const double d = from[plane.norm] - plane.coord;
1484 
1485  // right of plane?
1486  if( d > 0 )
1487  {
1488  node.handle = children[1];
1489  list.push_back( node );
1490  // if the split plane is close to the input point, add
1491  // the left child also (we'll check the exact distance
1492  /// when we pop it from the list.)
1493  if( d <= distance )
1494  {
1495  node.dist[plane.norm] = d;
1496  if( node.dist % node.dist <= dist_sqr )
1497  {
1498  node.handle = children[0];
1499  list.push_back( node );
1500  }
1501  }
1502  }
1503  // left of plane
1504  else
1505  {
1506  node.handle = children[0];
1507  list.push_back( node );
1508  // if the split plane is close to the input point, add
1509  // the right child also (we'll check the exact distance
1510  /// when we pop it from the list.)
1511  if( -d <= distance )
1512  {
1513  node.dist[plane.norm] = -d;
1514  if( node.dist % node.dist <= dist_sqr )
1515  {
1516  node.handle = children[1];
1517  list.push_back( node );
1518  }
1519  }
1520  }
1521  }
1522 
1523  if( myEval && result_params ) return MB_SUCCESS;
1524 
1525  // separate loops to avoid if test inside loop
1526 
1527  result_list.reserve( result_list_nodes.size() );
1528  for( std::vector< NodeDistance >::iterator vit = result_list_nodes.begin(); vit != result_list_nodes.end(); ++vit )
1529  result_list.push_back( ( *vit ).handle );
1530 
1531  if( result_dists && distance > 0.0 )
1532  {
1533  result_dists->reserve( result_list_nodes.size() );
1534  for( std::vector< NodeDistance >::iterator vit = result_list_nodes.begin(); vit != result_list_nodes.end();
1535  ++vit )
1536  result_dists->push_back( ( *vit ).dist.length() );
1537  }
1538 
1539  return MB_SUCCESS;
1540 }

References moab::CartVect::array(), moab::BoundBox::bMax, moab::BoundBox::bMin, children, moab::BoundBox::contains_point(), moab::AdaptiveKDTree::Plane::coord, moab::NodeDistance::dist, ErrorCode, moab::ElemEvaluator::find_containing_entity(), moab::Tree::get_bounding_box(), moab::Interface::get_child_meshsets(), get_split_plane(), moab::NodeDistance::handle, moab::TreeStats::leavesVisited, moab::Tree::maxDepth, MB_SUCCESS, moab::Tree::moab(), moab::Tree::myEval, moab::Tree::myRoot, moab::TreeStats::nodesVisited, moab::AdaptiveKDTree::Plane::norm, moab::TreeStats::numTraversals, moab::TreeStats::traversalLeafObjectTests, and moab::Tree::treeStats.

Referenced by closest_triangle(), moab::MergeMesh::find_merged_to(), moab::Intx2Mesh::intersect_meshes(), moab::Intx2Mesh::intersect_meshes_kdtree(), merge_duplicate_vertices(), moab::Coupler::nat_param(), sphere_intersect_triangles(), and moab::ReadNCDF::update().

◆ find_close_triangle()

ErrorCode moab::AdaptiveKDTree::find_close_triangle ( EntityHandle  root,
const double  from_point[3],
double  pt[3],
EntityHandle triangle 
)
private

find a triangle near the input point

Definition at line 1593 of file AdaptiveKDTree.cpp.

1597 {
1598  ErrorCode rval;
1599  Range tris;
1600  Plane split;
1601  std::vector< EntityHandle > stack;
1602  std::vector< EntityHandle > children( 2 );
1603  stack.reserve( 30 );
1604  assert( root );
1605  stack.push_back( root );
1606 
1607  while( !stack.empty() )
1608  {
1609  EntityHandle node = stack.back();
1610  stack.pop_back();
1611 
1612  for( ;; )
1613  { // loop until we find a leaf
1614 
1615  children.clear();
1616  rval = moab()->get_child_meshsets( node, children );
1617  if( MB_SUCCESS != rval ) return rval;
1618 
1619  // loop termination criterion
1620  if( children.empty() ) break;
1621 
1622  // if not a leaf, get split plane
1623  rval = get_split_plane( node, split );
1624  if( MB_SUCCESS != rval ) return rval;
1625 
1626  // continue down the side that contains the point,
1627  // and push the other side onto the stack in case
1628  // we need to check it later.
1629  int rs = split.right_side( from );
1630  node = children[rs];
1631  stack.push_back( children[1 - rs] );
1632  }
1633 
1634  // We should now be at a leaf.
1635  // If it has some triangles, we're done.
1636  // If not, continue searching for another leaf.
1637  tris.clear();
1638  rval = moab()->get_entities_by_type( node, MBTRI, tris );
1639  if( !tris.empty() )
1640  {
1641  double dist_sqr = HUGE_VAL;
1642  CartVect point( pt );
1643  rval = closest_to_triangles( moab(), tris, CartVect( from ), dist_sqr, point, triangle );
1644  point.get( pt );
1645  return rval;
1646  }
1647  }
1648 
1649  // If we got here, then we traversed the entire tree
1650  // and all the leaves were empty.
1651  return MB_ENTITY_NOT_FOUND;
1652 }

References children, moab::Range::clear(), moab::closest_to_triangles(), moab::Range::empty(), ErrorCode, moab::CartVect::get(), moab::Interface::get_child_meshsets(), moab::Interface::get_entities_by_type(), get_split_plane(), MB_ENTITY_NOT_FOUND, MB_SUCCESS, MBTRI, moab::Tree::moab(), and split.

Referenced by closest_triangle().

◆ get_info()

ErrorCode moab::AdaptiveKDTree::get_info ( EntityHandle  root,
double  min[3],
double  max[3],
unsigned int &  max_dep 
)
virtual

Return some basic information about the tree Stats are returned for tree starting from input node or tree root (root = 0)

Parameters
rootIf non-0, give stats below and including root
minMinimum corner of bounding box
maxMaximum corner of bounding box
max_depMaximum depth of tree below root

Reimplemented from moab::Tree.

Definition at line 2055 of file AdaptiveKDTree.cpp.

2056 {
2057  BoundBox box;
2058  ErrorCode result = get_bounding_box( box, &root );
2059  if( MB_SUCCESS != result ) return result;
2060  box.bMin.get( bmin );
2061  box.bMax.get( bmax );
2062 
2063  unsigned min_depth;
2064  return compute_depth( root, min_depth, dep );
2065 }

References moab::BoundBox::bMax, moab::BoundBox::bMin, compute_depth(), ErrorCode, moab::CartVect::get(), moab::Tree::get_bounding_box(), and MB_SUCCESS.

Referenced by moab::Coupler::initialize_tree().

◆ get_last_iterator()

ErrorCode moab::AdaptiveKDTree::get_last_iterator ( EntityHandle  tree_root,
AdaptiveKDTreeIter result 
)

Get iterator at right-most ('last') leaf.

Definition at line 327 of file AdaptiveKDTree.cpp.

328 {
329  double box[6];
330  ErrorCode rval = moab()->tag_get_data( boxTag, &root, 1, box );
331  if( MB_SUCCESS != rval ) return rval;
332 
333  return iter.initialize( this, root, box, box + 3, AdaptiveKDTreeIter::RIGHT );
334 }

References moab::Tree::boxTag, ErrorCode, moab::AdaptiveKDTreeIter::initialize(), MB_SUCCESS, moab::Tree::moab(), moab::AdaptiveKDTreeIter::RIGHT, and moab::Interface::tag_get_data().

◆ get_split_plane()

ErrorCode moab::AdaptiveKDTree::get_split_plane ( EntityHandle  node,
Plane plane 
)

Get split plane for tree node.

Definition at line 285 of file AdaptiveKDTree.cpp.

286 {
287 #ifndef MB_AD_KD_TREE_USE_SINGLE_TAG
288  ErrorCode r1, r2;
289  r1 = moab()->tag_get_data( planeTag, &entity, 1, &plane.coord );
290  r2 = moab()->tag_get_data( axisTag, &entity, 1, &plane.norm );
291  return MB_SUCCESS == r1 ? r2 : r1;
292 #elif defined( MB_AD_KD_TREE_USE_TWO_DOUBLE_TAG )
293  double values[2];
294  ErrorCode rval = moab()->tag_get_data( planeTag, &entity, 1, values );
295  plane.coord = values[0];
296  plane.norm = (int)values[1];
297  return rval;
298 #else
299  return moab()->tag_get_data( planeTag, &entity, 1, &plane );
300 #endif
301 }

References axisTag, moab::AdaptiveKDTree::Plane::coord, ErrorCode, MB_SUCCESS, moab::Tree::moab(), moab::AdaptiveKDTree::Plane::norm, planeTag, and moab::Interface::tag_get_data().

Referenced by distance_search(), find_close_triangle(), moab::AdaptiveKDTreeIter::get_neighbors(), moab::AdaptiveKDTreeIter::get_parent_split_plane(), merge_leaf(), point_search(), ray_intersect_triangles(), moab::AdaptiveKDTreeIter::sibling_side(), moab::AdaptiveKDTreeIter::step(), and moab::AdaptiveKDTreeIter::step_to_first_leaf().

◆ get_sub_tree_iterator()

ErrorCode moab::AdaptiveKDTree::get_sub_tree_iterator ( EntityHandle  tree_root,
const double  box_min[3],
const double  box_max[3],
AdaptiveKDTreeIter result 
)

Get iterator for tree or subtree.

Definition at line 336 of file AdaptiveKDTree.cpp.

340 {
341  return result.initialize( this, root, min, max, AdaptiveKDTreeIter::LEFT );
342 }

References moab::AdaptiveKDTreeIter::initialize(), and moab::AdaptiveKDTreeIter::LEFT.

Referenced by get_tree_iterator().

◆ get_tree_iterator()

ErrorCode moab::AdaptiveKDTree::get_tree_iterator ( EntityHandle  tree_root,
AdaptiveKDTreeIter result 
)

Get iterator for tree.

Definition at line 318 of file AdaptiveKDTree.cpp.

319 {
320  double box[6];
321  ErrorCode rval = moab()->tag_get_data( boxTag, &root, 1, box );
322  if( MB_SUCCESS != rval ) return rval;
323 
324  return get_sub_tree_iterator( root, box, box + 3, iter );
325 }

References moab::Tree::boxTag, ErrorCode, get_sub_tree_iterator(), MB_SUCCESS, moab::Tree::moab(), and moab::Interface::tag_get_data().

Referenced by compute_depth(), moab::MergeMesh::find_merged_to(), moab::Coupler::nat_param(), print(), and moab::ReadNCDF::update().

◆ init()

ErrorCode moab::AdaptiveKDTree::init ( )
private

Definition at line 240 of file AdaptiveKDTree.cpp.

241 {
242  std::vector< Tag > ctl;
243 
244 #ifndef MB_AD_KD_TREE_USE_SINGLE_TAG
245  // create two tags, one for axis direction and one for axis coordinate
246  std::string n1( treeName ), n2( treeName );
247  n1 += "_coord";
248  n2 += "_norm";
249  ErrorCode rval = make_tag( moab(), n1, MB_TAG_DENSE, MB_TYPE_DOUBLE, 1, 0, planeTag, ctl );
250  if( MB_SUCCESS != rval ) return rval;
251  rval = make_tag( moab(), n2, MB_TAG_DENSE, MB_TYPE_INT, 1, 0, axisTag, ctl );
252  if( MB_SUCCESS != rval ) return rval;
253 
254 #elif defined( MB_AD_KD_TREE_USE_TWO_DOUBLE_TAG )
255  // create tag to hold two doubles, one for location and one for axis
256  std::string double_tag_name = std::string( treeName ) + std::string( "_coord_norm" );
257  ErrorCode rval = make_tag( moab(), double_tag_name, MB_TAG_DENSE, MB_TYPE_DOUBLE, 2, 0, planeTag, ctl );
258  if( MB_SUCCESS != rval ) return rval;
259 #else
260  // create opaque tag to hold struct Plane
261  ErrorCode rval = make_tag( moab(), tagname, MB_TAG_DENSE, MB_TYPE_OPAQUE, sizeof( Plane ), 0, planeTag, ctl );
262  if( MB_SUCCESS != rval ) return rval;
263 
264 #ifdef MOAB_HAVE_HDF5
265  // create a mesh tag holding the HDF5 type for a struct Plane
266  Tag type_tag;
267  std::string type_tag_name = "__hdf5_tag_type_";
268  type_tag_name += boxTagName;
269  rval = make_tag( moab(), type_tag_name, MB_TAG_MESH, MB_TYPE_OPAQUE, sizeof( hid_t ), 0, type_tag, ctl );
270  if( MB_SUCCESS != rval ) return rval;
271  // create HDF5 type object describing struct Plane
272  Plane p;
273  hid_t handle = H5Tcreate( H5T_COMPOUND, sizeof( Plane ) );
274  H5Tinsert( handle, "coord", &( p.coord ) - &p, H5T_NATIVE_DOUBLE );
275  H5Tinsert( handle, "norm", &( p.axis ) - &p, H5T_NATIVE_INT );
276  EntityHandle root = 0;
277  rval = mbImpl->tag_set_data( type_tag, &root, 1, &handle );
278  if( MB_SUCCESS != rval ) return rval;
279 #endif
280 #endif
281 
282  return rval;
283 }

References axisTag, moab::Tree::boxTagName, moab::AdaptiveKDTree::Plane::coord, ErrorCode, make_tag(), MB_SUCCESS, MB_TAG_DENSE, MB_TAG_MESH, MB_TYPE_DOUBLE, MB_TYPE_OPAQUE, moab::Tree::mbImpl, moab::Tree::moab(), planeTag, moab::Interface::tag_set_data(), and treeName.

Referenced by AdaptiveKDTree().

◆ intersect_children_with_elems()

ErrorCode moab::AdaptiveKDTree::intersect_children_with_elems ( const Range elems,
AdaptiveKDTree::Plane  plane,
double  eps,
CartVect  box_min,
CartVect  box_max,
Range left_tris,
Range right_tris,
Range both_tris,
double &  metric_value 
)
private

Definition at line 777 of file AdaptiveKDTree.cpp.

786 {
787  left_tris.clear();
788  right_tris.clear();
789  both_tris.clear();
790  CartVect coords[16];
791 
792  // get extents of boxes for left and right sides
793  BoundBox left_box( box_min, box_max ), right_box( box_min, box_max );
794  right_box.bMin = box_min;
795  left_box.bMax = box_max;
796  right_box.bMin[plane.norm] = left_box.bMax[plane.norm] = plane.coord;
797  const CartVect left_cen = 0.5 * ( left_box.bMax + box_min );
798  const CartVect left_dim = 0.5 * ( left_box.bMax - box_min );
799  const CartVect right_cen = 0.5 * ( box_max + right_box.bMin );
800  const CartVect right_dim = 0.5 * ( box_max - right_box.bMin );
801  const CartVect dim = box_max - box_min;
802  const double max_tol = std::max( dim[0], std::max( dim[1], dim[2] ) ) / 10;
803 
804  // test each entity
805  ErrorCode rval;
806  int count, count2;
807  const EntityHandle *conn, *conn2;
808 
809  const Range::const_iterator elem_begin = elems.lower_bound( MBEDGE );
810  const Range::const_iterator poly_begin = elems.lower_bound( MBPOLYHEDRON, elem_begin );
811  const Range::const_iterator set_begin = elems.lower_bound( MBENTITYSET, poly_begin );
812  Range::iterator left_ins = left_tris.begin();
813  Range::iterator right_ins = right_tris.begin();
814  Range::iterator both_ins = both_tris.begin();
815  Range::const_iterator i;
816 
817  // vertices
818  for( i = elems.begin(); i != elem_begin; ++i )
819  {
821  rval = moab()->get_coords( &*i, 1, coords[0].array() );
822  if( MB_SUCCESS != rval ) return rval;
823 
824  bool lo = false, ro = false;
825  if( coords[0][plane.norm] <= plane.coord ) lo = true;
826  if( coords[0][plane.norm] >= plane.coord ) ro = true;
827 
828  if( lo && ro )
829  both_ins = both_tris.insert( both_ins, *i, *i );
830  else if( lo )
831  left_ins = left_tris.insert( left_ins, *i, *i );
832  else // if (ro)
833  right_ins = right_tris.insert( right_ins, *i, *i );
834  }
835 
836  // non-polyhedron elements
837  std::vector< EntityHandle > dum_vector;
838  for( i = elem_begin; i != poly_begin; ++i )
839  {
841  rval = moab()->get_connectivity( *i, conn, count, true, &dum_vector );
842  if( MB_SUCCESS != rval ) return rval;
843  if( count > (int)( sizeof( coords ) / sizeof( coords[0] ) ) ) return MB_FAILURE;
844  rval = moab()->get_coords( &conn[0], count, coords[0].array() );
845  if( MB_SUCCESS != rval ) return rval;
846 
847  bool lo = false, ro = false;
848  for( int j = 0; j < count; ++j )
849  {
850  if( coords[j][plane.norm] <= plane.coord ) lo = true;
851  if( coords[j][plane.norm] >= plane.coord ) ro = true;
852  }
853 
854  // Triangle must be in at least one leaf. If test against plane
855  // identified that leaf, then we're done. If triangle is on both
856  // sides of plane, do more precise test to ensure that it is really
857  // in both.
858  // BoundBox box;
859  // box.update(*moab(), *i);
860  if( lo && ro )
861  {
862  double tol = eps;
863  lo = ro = false;
864  while( !lo && !ro && tol <= max_tol )
865  {
866  tree_stats().boxElemTests += 2;
867  lo = GeomUtil::box_elem_overlap( coords, TYPE_FROM_HANDLE( *i ), left_cen, left_dim + CartVect( tol ),
868  count );
869  ro = GeomUtil::box_elem_overlap( coords, TYPE_FROM_HANDLE( *i ), right_cen, right_dim + CartVect( tol ),
870  count );
871 
872  tol *= 10.0;
873  }
874  }
875  if( lo && ro )
876  both_ins = both_tris.insert( both_ins, *i, *i );
877  else if( lo )
878  left_ins = left_tris.insert( left_ins, *i, *i );
879  else if( ro )
880  right_ins = right_tris.insert( right_ins, *i, *i );
881  }
882 
883  // polyhedra
884  for( i = poly_begin; i != set_begin; ++i )
885  {
887  rval = moab()->get_connectivity( *i, conn, count, true );
888  if( MB_SUCCESS != rval ) return rval;
889 
890  // just check the bounding box of the polyhedron
891  bool lo = false, ro = false;
892  for( int j = 0; j < count; ++j )
893  {
894  rval = moab()->get_connectivity( conn[j], conn2, count2, true );
895  if( MB_SUCCESS != rval ) return rval;
896 
897  for( int k = 0; k < count2; ++k )
898  {
899  rval = moab()->get_coords( conn2 + k, 1, coords[0].array() );
900  if( MB_SUCCESS != rval ) return rval;
901  if( coords[0][plane.norm] <= plane.coord ) lo = true;
902  if( coords[0][plane.norm] >= plane.coord ) ro = true;
903  }
904  }
905 
906  if( lo && ro )
907  both_ins = both_tris.insert( both_ins, *i, *i );
908  else if( lo )
909  left_ins = left_tris.insert( left_ins, *i, *i );
910  else if( ro )
911  right_ins = right_tris.insert( right_ins, *i, *i );
912  }
913 
914  // sets
915  BoundBox tbox;
916  for( i = set_begin; i != elems.end(); ++i )
917  {
919  rval = tbox.update( *moab(), *i, spherical, radius );
920  if( MB_SUCCESS != rval ) return rval;
921 
922  bool lo = false, ro = false;
923  if( tbox.bMin[plane.norm] <= plane.coord ) lo = true;
924  if( tbox.bMax[plane.norm] >= plane.coord ) ro = true;
925 
926  if( lo && ro )
927  both_ins = both_tris.insert( both_ins, *i, *i );
928  else if( lo )
929  left_ins = left_tris.insert( left_ins, *i, *i );
930  else // if (ro)
931  right_ins = right_tris.insert( right_ins, *i, *i );
932  }
933 
934  CartVect box_dim = box_max - box_min;
935  double area_left = left_dim[0] * left_dim[1] + left_dim[1] * left_dim[2] + left_dim[2] * left_dim[0];
936  double area_right = right_dim[0] * right_dim[1] + right_dim[1] * right_dim[2] + right_dim[2] * right_dim[0];
937  double area_both = box_dim[0] * box_dim[1] + box_dim[1] * box_dim[2] + box_dim[2] * box_dim[0];
938  metric_value = ( area_left * left_tris.size() + area_right * right_tris.size() ) / area_both + both_tris.size();
939  return MB_SUCCESS;
940 }

References moab::Range::begin(), moab::BoundBox::bMax, moab::BoundBox::bMin, moab::GeomUtil::box_elem_overlap(), box_max(), box_min(), moab::TreeStats::boxElemTests, moab::Range::clear(), moab::TreeStats::constructLeafObjectTests, moab::AdaptiveKDTree::Plane::coord, dim, moab::Range::end(), ErrorCode, moab::Interface::get_connectivity(), moab::Interface::get_coords(), moab::Range::insert(), left_box, moab::Range::lower_bound(), MB_SUCCESS, MBEDGE, MBENTITYSET, MBPOLYHEDRON, moab::Tree::moab(), moab::AdaptiveKDTree::Plane::norm, radius, right_box, moab::Range::size(), spherical, moab::Tree::tree_stats(), moab::TYPE_FROM_HANDLE(), and moab::BoundBox::update().

Referenced by best_subdivision_plane(), best_subdivision_snap_plane(), best_vertex_median_plane(), and best_vertex_sample_plane().

◆ make_tag()

ErrorCode moab::AdaptiveKDTree::make_tag ( Interface iface,
std::string  name,
TagType  storage,
DataType  type,
int  count,
void *  default_val,
Tag tag_handle,
std::vector< Tag > &  created_tags 
)
private

Definition at line 209 of file AdaptiveKDTree.cpp.

217 {
218  ErrorCode rval =
219  iface->tag_get_handle( name.c_str(), count, type, tag_handle, MB_TAG_CREAT | storage, default_val );
220 
221  if( MB_SUCCESS == rval )
222  {
223  if( std::find( created_tags.begin(), created_tags.end(), tag_handle ) == created_tags.end() )
224  created_tags.push_back( tag_handle );
225  }
226  else
227  {
228  while( !created_tags.empty() )
229  {
230  iface->tag_delete( created_tags.back() );
231  created_tags.pop_back();
232  }
233 
234  planeTag = axisTag = (Tag)-1;
235  }
236 
237  return rval;
238 }

References axisTag, ErrorCode, iface, MB_SUCCESS, MB_TAG_CREAT, and planeTag.

Referenced by init().

◆ merge_leaf()

ErrorCode moab::AdaptiveKDTree::merge_leaf ( AdaptiveKDTreeIter iter)

Merge the leaf pointed to by the current iterator with it's sibling. If the sibling is not a leaf, multiple merges may be done.

Definition at line 419 of file AdaptiveKDTree.cpp.

420 {
421  ErrorCode rval;
422  if( iter.depth() == 1 ) // at root
423  return MB_FAILURE;
424 
425  // Move iter to parent
426 
427  AdaptiveKDTreeIter::StackObj node = iter.mStack.back();
428  iter.mStack.pop_back();
429 
430  iter.childVect.clear();
431  rval = moab()->get_child_meshsets( iter.mStack.back().entity, iter.childVect );
432  if( MB_SUCCESS != rval ) return rval;
433  Plane plane;
434  rval = get_split_plane( iter.mStack.back().entity, plane );
435  if( MB_SUCCESS != rval ) return rval;
436 
437  int child_idx = iter.childVect[0] == node.entity ? 0 : 1;
438  assert( iter.childVect[child_idx] == node.entity );
439  iter.mBox[1 - child_idx][plane.norm] = node.coord;
440 
441  // Get all entities from children and put them in parent
442  EntityHandle parent = iter.handle();
443  moab()->remove_child_meshset( parent, iter.childVect[0] );
444  moab()->remove_child_meshset( parent, iter.childVect[1] );
445  std::vector< EntityHandle > stack( iter.childVect );
446 
447  Range range;
448  while( !stack.empty() )
449  {
450  EntityHandle h = stack.back();
451  stack.pop_back();
452  range.clear();
453  rval = moab()->get_entities_by_handle( h, range );
454  if( MB_SUCCESS != rval ) return rval;
455  rval = moab()->add_entities( parent, range );
456  if( MB_SUCCESS != rval ) return rval;
457 
458  iter.childVect.clear();
459  rval = moab()->get_child_meshsets( h, iter.childVect );MB_CHK_ERR( rval );
460  if( !iter.childVect.empty() )
461  {
462  moab()->remove_child_meshset( h, iter.childVect[0] );
463  moab()->remove_child_meshset( h, iter.childVect[1] );
464  stack.push_back( iter.childVect[0] );
465  stack.push_back( iter.childVect[1] );
466  }
467 
468  rval = moab()->delete_entities( &h, 1 );
469  if( MB_SUCCESS != rval ) return rval;
470  }
471 
472  return MB_SUCCESS;
473 }

References moab::Interface::add_entities(), moab::Range::back(), moab::AdaptiveKDTreeIter::childVect, moab::Range::clear(), moab::AdaptiveKDTreeIter::StackObj::coord, moab::Interface::delete_entities(), moab::AdaptiveKDTreeIter::depth(), moab::AdaptiveKDTreeIter::StackObj::entity, ErrorCode, moab::Interface::get_child_meshsets(), moab::Interface::get_entities_by_handle(), get_split_plane(), moab::AdaptiveKDTreeIter::handle(), MB_CHK_ERR, MB_SUCCESS, moab::AdaptiveKDTreeIter::mBox, moab::Tree::moab(), moab::AdaptiveKDTreeIter::mStack, moab::AdaptiveKDTree::Plane::norm, and moab::Interface::remove_child_meshset().

◆ parse_options()

ErrorCode moab::AdaptiveKDTree::parse_options ( FileOptions options)
virtual

Parse options for tree creation.

Parameters
optionsOptions passed in by application
Returns
Failure is returned if any options were passed in and not interpreted; could mean inappropriate options for a particular tree type

Implements moab::Tree.

Definition at line 176 of file AdaptiveKDTree.cpp.

177 {
178  ErrorCode rval = parse_common_options( opts );
179  if( MB_SUCCESS != rval ) return rval;
180 
181  // SPLITS_PER_DIR: number of candidate splits considered per direction; default = 3
182  int tmp_int;
183  rval = opts.get_int_option( "SPLITS_PER_DIR", tmp_int );
184  if( MB_SUCCESS == rval ) splitsPerDir = tmp_int;
185 
186  // PLANE_SET: method used to decide split planes; see CandidatePlaneSet enum (below)
187  // for possible values; default = 1 (SUBDIVISION_SNAP)
188  rval = opts.get_int_option( "PLANE_SET", tmp_int );
189  if( MB_SUCCESS == rval && ( tmp_int < SUBDIVISION || tmp_int > VERTEX_SAMPLE ) )
190  return MB_FAILURE;
191  else if( MB_ENTITY_NOT_FOUND == rval )
193  else
194  planeSet = (CandidatePlaneSet)( tmp_int );
195 
196  rval = opts.get_toggle_option( "SPHERICAL", false, spherical );
197  if( MB_SUCCESS != rval ) spherical = false;
198 
199  double tmp = 1.0;
200  rval = opts.get_real_option( "RADIUS", tmp );
201  if( MB_SUCCESS != rval )
202  radius = 1.0;
203  else
204  radius = tmp;
205 
206  return MB_SUCCESS;
207 }

References ErrorCode, moab::FileOptions::get_int_option(), moab::FileOptions::get_real_option(), moab::FileOptions::get_toggle_option(), MB_ENTITY_NOT_FOUND, MB_SUCCESS, moab::Tree::parse_common_options(), planeSet, radius, spherical, splitsPerDir, SUBDIVISION, and VERTEX_SAMPLE.

Referenced by AdaptiveKDTree(), build_tree(), moab::Intx2Mesh::intersect_meshes(), and moab::Intx2Mesh::intersect_meshes_kdtree().

◆ point_search() [1/2]

ErrorCode moab::AdaptiveKDTree::point_search ( const double *  point,
AdaptiveKDTreeIter leaf_it,
const double  iter_tol = 1.0e-10,
const double  inside_tol = 1.0e-6,
bool *  multiple_leaves = NULL,
EntityHandle start_node = NULL 
)

Get leaf containing input position.

Does not take into account global bounding box of tree.

  • Therefore there is always one leaf containing the point.
  • If caller wants to account for global bounding box, then caller can test against that box and not call this method at all if the point is outside the box, as there is no leaf containing the point in that case.
    Parameters
    pointPoint to be located in tree
    leaf_itIterator to leaf containing point
    iter_tolTolerance for convergence of point search
    inside_tolTolerance for inside element calculation
    multiple_leavesSome tree types can have multiple leaves containing a point; if non-NULL, this parameter is returned true if multiple leaves contain the input point
    start_nodeStart from this tree node (non-NULL) instead of tree root (NULL)
    Returns
    Non-success returned only in case of failure; not-found indicated by leaf_out=0

Definition at line 1321 of file AdaptiveKDTree.cpp.

1327 {
1328  ErrorCode rval;
1330 
1331  // kdtrees never have multiple leaves containing a pt
1332  if( multiple_leaves ) *multiple_leaves = false;
1333 
1334  leaf_it.mBox[0] = boundBox.bMin;
1335  leaf_it.mBox[1] = boundBox.bMax;
1336 
1337  // test that point is inside tree
1338  if( !boundBox.contains_point( point, iter_tol ) )
1339  {
1341  return MB_ENTITY_NOT_FOUND;
1342  }
1343 
1344  // initialize iterator at tree root
1345  leaf_it.treeTool = this;
1346  leaf_it.mStack.clear();
1347  leaf_it.mStack.push_back( AdaptiveKDTreeIter::StackObj( ( start_node ? *start_node : myRoot ), 0 ) );
1348 
1349  // loop until we reach a leaf
1350  AdaptiveKDTree::Plane plane;
1351  for( ;; )
1352  {
1354 
1355  // get children
1356  leaf_it.childVect.clear();
1357  rval = moab()->get_child_meshsets( leaf_it.handle(), leaf_it.childVect );
1358  if( MB_SUCCESS != rval ) return rval;
1359 
1360  // if no children, then at leaf (done)
1361  if( leaf_it.childVect.empty() )
1362  {
1364  break;
1365  }
1366 
1367  // get split plane
1368  rval = get_split_plane( leaf_it.handle(), plane );
1369  if( MB_SUCCESS != rval ) return rval;
1370 
1371  // step iterator to appropriate child
1372  // idx: 0->left, 1->right
1373  const int idx = ( point[plane.norm] > plane.coord );
1374  leaf_it.mStack.push_back(
1375  AdaptiveKDTreeIter::StackObj( leaf_it.childVect[idx], leaf_it.mBox[1 - idx][plane.norm] ) );
1376  leaf_it.mBox[1 - idx][plane.norm] = plane.coord;
1377  }
1378 
1379  return MB_SUCCESS;
1380 }

References moab::BoundBox::bMax, moab::BoundBox::bMin, moab::Tree::boundBox, moab::AdaptiveKDTreeIter::childVect, moab::BoundBox::contains_point(), moab::AdaptiveKDTree::Plane::coord, ErrorCode, moab::Interface::get_child_meshsets(), get_split_plane(), moab::AdaptiveKDTreeIter::handle(), moab::TreeStats::leavesVisited, MB_ENTITY_NOT_FOUND, MB_SUCCESS, moab::AdaptiveKDTreeIter::mBox, moab::Tree::moab(), moab::AdaptiveKDTreeIter::mStack, moab::Tree::myRoot, moab::TreeStats::nodesVisited, moab::AdaptiveKDTree::Plane::norm, moab::TreeStats::numTraversals, moab::Tree::treeStats, and moab::AdaptiveKDTreeIter::treeTool.

◆ point_search() [2/2]

ErrorCode moab::AdaptiveKDTree::point_search ( const double *  point,
EntityHandle leaf_out,
const double  iter_tol = 1.0e-10,
const double  inside_tol = 1.0e-6,
bool *  multiple_leaves = NULL,
EntityHandle start_node = NULL,
CartVect params = NULL 
)
virtual

Get leaf containing input position.

Does not take into account global bounding box of tree.

  • Therefore there is always one leaf containing the point.
  • If caller wants to account for global bounding box, then caller can test against that box and not call this method at all if the point is outside the box, as there is no leaf containing the point in that case.
    Parameters
    pointPoint to be located in tree
    leaf_outLeaf containing point
    iter_tolTolerance for convergence of point search
    inside_tolTolerance for inside element calculation
    multiple_leavesSome tree types can have multiple leaves containing a point; if non-NULL, this parameter is returned true if multiple leaves contain the input point
    start_nodeStart from this tree node (non-NULL) instead of tree root (NULL)
    Returns
    Non-success returned only in case of failure; not-found indicated by leaf_out=0

Implements moab::Tree.

Definition at line 1266 of file AdaptiveKDTree.cpp.

1273 {
1274  std::vector< EntityHandle > children;
1275  Plane plane;
1276 
1278  leaf_out = 0;
1279  BoundBox box;
1280  // kdtrees never have multiple leaves containing a pt
1281  if( multiple_leaves ) *multiple_leaves = false;
1282 
1283  EntityHandle node = ( start_node ? *start_node : myRoot );
1284 
1286  ErrorCode rval = get_bounding_box( box, &node );
1287  if( MB_SUCCESS != rval ) return rval;
1288  if( !box.contains_point( point, iter_tol ) ) return MB_SUCCESS;
1289 
1290  rval = moab()->get_child_meshsets( node, children );
1291  if( MB_SUCCESS != rval ) return rval;
1292 
1293  while( !children.empty() )
1294  {
1296 
1297  rval = get_split_plane( node, plane );
1298  if( MB_SUCCESS != rval ) return rval;
1299 
1300  const double d = point[plane.norm] - plane.coord;
1301  node = children[( d > 0.0 )];
1302 
1303  children.clear();
1304  rval = moab()->get_child_meshsets( node, children );
1305  if( MB_SUCCESS != rval ) return rval;
1306  }
1307 
1309  if( myEval && params )
1310  {
1311  rval = myEval->find_containing_entity( node, point, iter_tol, inside_tol, leaf_out, params->array(),
1313  if( MB_SUCCESS != rval ) return rval;
1314  }
1315  else
1316  leaf_out = node;
1317 
1318  return MB_SUCCESS;
1319 }

References moab::CartVect::array(), children, moab::BoundBox::contains_point(), moab::AdaptiveKDTree::Plane::coord, ErrorCode, moab::ElemEvaluator::find_containing_entity(), moab::Tree::get_bounding_box(), moab::Interface::get_child_meshsets(), get_split_plane(), moab::TreeStats::leavesVisited, MB_SUCCESS, moab::Tree::moab(), moab::Tree::myEval, moab::Tree::myRoot, moab::TreeStats::nodesVisited, moab::AdaptiveKDTree::Plane::norm, moab::TreeStats::numTraversals, moab::TreeStats::traversalLeafObjectTests, and moab::Tree::treeStats.

Referenced by moab::TempestRemapper::ConstructCoveringSet(), main(), and moab::Coupler::nat_param().

◆ print()

ErrorCode moab::AdaptiveKDTree::print ( )
virtual

print various things about this tree

Implements moab::Tree.

Definition at line 2129 of file AdaptiveKDTree.cpp.

2130 {
2131  Range range;
2132 
2133  Range tree_sets, elem2d, elem3d, verts, all;
2134  moab()->get_child_meshsets( myRoot, tree_sets, 0 );
2135 
2136  {
2137  std::vector<EntityHandle> elem2d_vec, elem3d_vec, verts_vec;
2138  for( Range::iterator rit = tree_sets.begin(); rit != tree_sets.end(); ++rit )
2139  {
2140  moab()->get_entities_by_dimension( *rit, 2, elem2d_vec );
2141  moab()->get_entities_by_dimension( *rit, 3, elem3d_vec );
2142  moab()->get_entities_by_type( *rit, MBVERTEX, verts_vec );
2143  }
2144  std::sort(elem2d_vec.begin(), elem2d_vec.end());
2145  std::copy( elem2d_vec.rbegin(), elem2d_vec.rend(), range_inserter( elem2d ) );
2146  std::sort(elem3d_vec.begin(), elem3d_vec.end());
2147  std::copy( elem3d_vec.rbegin(), elem3d_vec.rend(), range_inserter( elem3d ) );
2148  std::sort(verts_vec.begin(), verts_vec.end());
2149  std::copy( verts_vec.rbegin(), verts_vec.rend(), range_inserter( verts ) );
2150  }
2151 
2152  all.merge( verts );
2153  all.merge( elem2d );
2154  all.merge( elem3d );
2155  tree_sets.insert( myRoot );
2156  unsigned long long set_used, set_amortized, set_store_used, set_store_amortized, set_tag_used, set_tag_amortized,
2157  elem_used, elem_amortized;
2158  moab()->estimated_memory_use( tree_sets, &set_used, &set_amortized, &set_store_used, &set_store_amortized, 0, 0, 0,
2159  0, &set_tag_used, &set_tag_amortized );
2160  moab()->estimated_memory_use( all, &elem_used, &elem_amortized );
2161 
2162  int num_2d = 0, num_3d = 0;
2163  ;
2164  moab()->get_number_entities_by_dimension( 0, 2, num_2d );
2165  moab()->get_number_entities_by_dimension( 0, 3, num_3d );
2166 
2167  BoundBox box;
2168  ErrorCode rval = get_bounding_box( box, &myRoot );
2169  if( MB_SUCCESS != rval || box == BoundBox() ) throw rval;
2170  double diff[3] = { box.bMax[0] - box.bMin[0], box.bMax[1] - box.bMin[1], box.bMax[2] - box.bMin[2] };
2171  double tree_vol = diff[0] * diff[1] * diff[2];
2172  double tree_surf_area = 2 * ( diff[0] * diff[1] + diff[1] * diff[2] + diff[2] * diff[0] );
2173 
2174  SimpleStat< unsigned > depth, size;
2175  SimpleStat< double > vol, surf;
2176 
2177  AdaptiveKDTreeIter iter;
2178  get_tree_iterator( myRoot, iter );
2179  do
2180  {
2181  depth.add( iter.depth() );
2182 
2183  int num_leaf_elem;
2184  moab()->get_number_entities_by_handle( iter.handle(), num_leaf_elem );
2185  size.add( num_leaf_elem );
2186 
2187  const double* n = iter.box_min();
2188  const double* x = iter.box_max();
2189  double dims[3] = { x[0] - n[0], x[1] - n[1], x[2] - n[2] };
2190 
2191  double leaf_vol = dims[0] * dims[1] * dims[2];
2192  vol.add( leaf_vol );
2193 
2194  double area = 2.0 * ( dims[0] * dims[1] + dims[1] * dims[2] + dims[2] * dims[0] );
2195  surf.add( area );
2196 
2197  } while( MB_SUCCESS == iter.step() );
2198 
2199  printf( "------------------------------------------------------------------\n" );
2200  printf( "tree volume: %f\n", tree_vol );
2201  printf( "total elements: %d\n", num_2d + num_3d );
2202  printf( "number of leaves: %lu\n", (unsigned long)depth.count );
2203  printf( "number of nodes: %lu\n", (unsigned long)tree_sets.size() );
2204  printf( "volume ratio: %0.2f%%\n", 100 * ( vol.sum / tree_vol ) );
2205  printf( "surface ratio: %0.2f%%\n", 100 * ( surf.sum / tree_surf_area ) );
2206  printf( "\nmemory: used amortized\n" );
2207  printf( " ---------- ----------\n" );
2208  printf( "elements %10s %10s\n", mem_to_string( elem_used ).c_str(), mem_to_string( elem_amortized ).c_str() );
2209  printf( "sets (total)%10s %10s\n", mem_to_string( set_used ).c_str(), mem_to_string( set_amortized ).c_str() );
2210  printf( "sets %10s %10s\n", mem_to_string( set_store_used ).c_str(),
2211  mem_to_string( set_store_amortized ).c_str() );
2212  printf( "set tags %10s %10s\n", mem_to_string( set_tag_used ).c_str(),
2213  mem_to_string( set_tag_amortized ).c_str() );
2214  printf( "\nleaf stats: min avg rms max std.dev\n" );
2215  printf( " ---------- ---------- ---------- ---------- ----------\n" );
2216  printf( "depth %10u %10.1f %10.1f %10u %10.2f\n", depth.min, depth.avg(), depth.rms(), depth.max,
2217  depth.dev() );
2218  printf( "triangles %10u %10.1f %10.1f %10u %10.2f\n", size.min, size.avg(), size.rms(), size.max, size.dev() );
2219  printf( "volume %10.2g %10.2g %10.2g %10.2g %10.2g\n", vol.min, vol.avg(), vol.rms(), vol.max, vol.dev() );
2220  printf( "surf. area %10.2g %10.2g %10.2g %10.2g %10.2g\n", surf.min, surf.avg(), surf.rms(), surf.max,
2221  surf.dev() );
2222  printf( "------------------------------------------------------------------\n" );
2223 
2224  return MB_SUCCESS;
2225 }

References moab::SimpleStat< T >::add(), moab::SimpleStat< T >::avg(), moab::Range::begin(), moab::BoundBox::bMax, moab::BoundBox::bMin, moab::AdaptiveKDTreeIter::box_max(), moab::AdaptiveKDTreeIter::box_min(), moab::SimpleStat< T >::count, moab::AdaptiveKDTreeIter::depth(), moab::SimpleStat< T >::dev(), moab::Range::end(), ErrorCode, moab::Interface::estimated_memory_use(), moab::Tree::get_bounding_box(), moab::Interface::get_child_meshsets(), moab::Interface::get_entities_by_dimension(), moab::Interface::get_entities_by_type(), moab::Interface::get_number_entities_by_dimension(), moab::Interface::get_number_entities_by_handle(), get_tree_iterator(), moab::AdaptiveKDTreeIter::handle(), moab::Range::insert(), moab::SimpleStat< T >::max, MB_SUCCESS, MBVERTEX, moab::mem_to_string(), moab::Range::merge(), moab::SimpleStat< T >::min, moab::Tree::moab(), moab::Tree::myRoot, moab::SimpleStat< T >::rms(), size, moab::Range::size(), moab::AdaptiveKDTreeIter::step(), and moab::SimpleStat< T >::sum.

◆ ray_intersect_triangles()

ErrorCode moab::AdaptiveKDTree::ray_intersect_triangles ( EntityHandle  tree_root,
const double  tolerance,
const double  ray_unit_dir[3],
const double  ray_base_pt[3],
std::vector< EntityHandle > &  triangles_out,
std::vector< double > &  distance_out,
int  result_count_limit = 0,
double  distance_limit = -1.0 
)

Definition at line 1851 of file AdaptiveKDTree.cpp.

1859 {
1860  ErrorCode rval;
1861  double ray_beg = 0.0;
1862  if( ray_end < 0.0 ) ray_end = HUGE_VAL;
1863 
1864  // if root has bounding box, trim ray to that box
1865  CartVect tvec( tol );
1866  BoundBox box;
1867  const CartVect ray_pt( ray_pt_in ), ray_dir( ray_dir_in );
1868  rval = get_bounding_box( box );
1869  if( MB_SUCCESS == rval )
1870  {
1871  if( !GeomUtil::segment_box_intersect( box.bMin - tvec, box.bMax + tvec, ray_pt, ray_dir, ray_beg, ray_end ) )
1872  return MB_SUCCESS; // ray misses entire tree.
1873  }
1874 
1875  Range tris;
1876  Range::iterator iter;
1877  CartVect tri_coords[3];
1878  const EntityHandle* tri_conn;
1879  int conn_len;
1880  double tri_t;
1881 
1882  Plane plane;
1883  std::vector< EntityHandle > children;
1884  std::vector< NodeSeg > list;
1885  NodeSeg seg( root, ray_beg, ray_end );
1886  list.push_back( seg );
1887 
1888  while( !list.empty() )
1889  {
1890  seg = list.back();
1891  list.pop_back();
1892 
1893  // If we are limited to a certain number of intersections
1894  // (max_ints != 0), then ray_end will contain the distance
1895  // to the furthest intersection we have so far. If the
1896  // tree node is further than that, skip it.
1897  if( seg.beg > ray_end ) continue;
1898 
1899  // Check if at a leaf
1900  children.clear();
1901  rval = moab()->get_child_meshsets( seg.handle, children );
1902  if( MB_SUCCESS != rval ) return rval;
1903  if( children.empty() )
1904  { // leaf
1905 
1906  tris.clear();
1907  rval = moab()->get_entities_by_type( seg.handle, MBTRI, tris );
1908  if( MB_SUCCESS != rval ) return rval;
1909 
1910  for( iter = tris.begin(); iter != tris.end(); ++iter )
1911  {
1912  rval = moab()->get_connectivity( *iter, tri_conn, conn_len );
1913  if( MB_SUCCESS != rval ) return rval;
1914  rval = moab()->get_coords( tri_conn, 3, tri_coords[0].array() );
1915  if( MB_SUCCESS != rval ) return rval;
1916 
1917  if( GeomUtil::ray_tri_intersect( tri_coords, ray_pt, ray_dir, tri_t, &ray_end ) )
1918  {
1919  if( !max_ints )
1920  {
1921  if( std::find( tris_out.begin(), tris_out.end(), *iter ) == tris_out.end() )
1922  {
1923  tris_out.push_back( *iter );
1924  dists_out.push_back( tri_t );
1925  }
1926  }
1927  else if( tri_t < ray_end )
1928  {
1929  if( std::find( tris_out.begin(), tris_out.end(), *iter ) == tris_out.end() )
1930  {
1931  if( tris_out.size() < (unsigned)max_ints )
1932  {
1933  tris_out.resize( tris_out.size() + 1 );
1934  dists_out.resize( dists_out.size() + 1 );
1935  }
1936  int w = tris_out.size() - 1;
1937  for( ; w > 0 && tri_t < dists_out[w - 1]; --w )
1938  {
1939  tris_out[w] = tris_out[w - 1];
1940  dists_out[w] = dists_out[w - 1];
1941  }
1942  tris_out[w] = *iter;
1943  dists_out[w] = tri_t;
1944  if( tris_out.size() >= (unsigned)max_ints )
1945  // when we have already reached the max intx points, we cans safely
1946  // reset ray_end, because we will accept new points only "closer"
1947  // than the last one
1948  ray_end = dists_out.back();
1949  }
1950  }
1951  }
1952  }
1953 
1954  continue;
1955  }
1956 
1957  rval = get_split_plane( seg.handle, plane );
1958  if( MB_SUCCESS != rval ) return rval;
1959 
1960  // Consider two planes that are the split plane +/- the tolerance.
1961  // Calculate the segment parameter at which the line segment intersects
1962  // the true plane, and also the difference between that value and the
1963  // intersection with either of the +/- tol planes.
1964  const double inv_dir = 1.0 / ray_dir[plane.norm]; // only do division once
1965  const double t = ( plane.coord - ray_pt[plane.norm] ) * inv_dir; // intersection with plane
1966  const double diff = tol * inv_dir; // t adjustment for +tol plane
1967  // const double t0 = t - diff; // intersection with -tol plane
1968  // const double t1 = t + diff; // intersection with +tol plane
1969 
1970  // The index of the child tree node (0 or 1) that is on the
1971  // side of the plane to which the ray direction points. That is,
1972  // if the ray direction is opposite the plane normal, the index
1973  // of the child corresponding to the side beneath the plane. If
1974  // the ray direction is the same as the plane normal, the index
1975  // of the child corresponding to the side above the plane.
1976  const int fwd_child = ( ray_dir[plane.norm] > 0.0 );
1977 
1978  // Note: we maintain seg.beg <= seg.end at all times, so assume that here.
1979 
1980  // If segment is parallel to plane
1981  if( !Util::is_finite( t ) )
1982  {
1983  if( ray_pt[plane.norm] - tol <= plane.coord ) list.push_back( NodeSeg( children[0], seg.beg, seg.end ) );
1984  if( ray_pt[plane.norm] + tol >= plane.coord ) list.push_back( NodeSeg( children[1], seg.beg, seg.end ) );
1985  }
1986  // If segment is entirely to one side of plane such that the
1987  // intersection with the split plane is past the end of the segment
1988  else if( seg.end + diff < t )
1989  {
1990  // If segment direction is opposite that of plane normal, then
1991  // being past the end of the segment means that we are to the
1992  // right (or above) the plane and what the right child (index == 1).
1993  // Otherwise we want the left child (index == 0);
1994  list.push_back( NodeSeg( children[1 - fwd_child], seg.beg, seg.end ) );
1995  }
1996  // If the segment is entirely to one side of the plane such that
1997  // the intersection with the split plane is before the start of the
1998  // segment
1999  else if( seg.beg - diff > t )
2000  {
2001  // If segment direction is opposite that of plane normal, then
2002  // being before the start of the segment means that we are to the
2003  // left (or below) the plane and what the left child (index == 0).
2004  // Otherwise we want the right child (index == 1);
2005  list.push_back( NodeSeg( children[fwd_child], seg.beg, seg.end ) );
2006  }
2007  // Otherwise we must intersect the plane.
2008  // Note: be careful not to grow the segment if t is slightly
2009  // outside the current segment, as doing so would effectively
2010  // increase the tolerance as we descend the tree.
2011  else if( t <= seg.beg )
2012  {
2013  list.push_back( NodeSeg( children[1 - fwd_child], seg.beg, seg.beg ) );
2014  list.push_back( NodeSeg( children[fwd_child], seg.beg, seg.end ) );
2015  }
2016  else if( t >= seg.end )
2017  {
2018  list.push_back( NodeSeg( children[1 - fwd_child], seg.beg, seg.end ) );
2019  list.push_back( NodeSeg( children[fwd_child], seg.end, seg.end ) );
2020  }
2021  else
2022  {
2023  list.push_back( NodeSeg( children[1 - fwd_child], seg.beg, t ) );
2024  list.push_back( NodeSeg( children[fwd_child], t, seg.end ) );
2025  }
2026  }
2027 
2028  return MB_SUCCESS;
2029 }

References moab::NodeSeg::beg, moab::Range::begin(), moab::BoundBox::bMax, moab::BoundBox::bMin, children, moab::Range::clear(), moab::AdaptiveKDTree::Plane::coord, moab::NodeSeg::end, moab::Range::end(), ErrorCode, moab::Tree::get_bounding_box(), moab::Interface::get_child_meshsets(), moab::Interface::get_connectivity(), moab::Interface::get_coords(), moab::Interface::get_entities_by_type(), get_split_plane(), moab::NodeSeg::handle, moab::Util::is_finite(), MB_SUCCESS, MBTRI, moab::Tree::moab(), moab::AdaptiveKDTree::Plane::norm, moab::GeomUtil::ray_tri_intersect(), and moab::GeomUtil::segment_box_intersect().

◆ reset_tree()

ErrorCode moab::AdaptiveKDTree::reset_tree ( )
inlinevirtual

Reset the tree, optionally checking we have the right root.

Implements moab::Tree.

Definition at line 528 of file AdaptiveKDTree.hpp.

529 {
530  return delete_tree_sets();
531 }

References moab::Tree::delete_tree_sets().

Referenced by build_tree(), main(), and ~AdaptiveKDTree().

◆ set_split_plane()

ErrorCode moab::AdaptiveKDTree::set_split_plane ( EntityHandle  node,
const Plane plane 
)

Set split plane for tree node.

Definition at line 303 of file AdaptiveKDTree.cpp.

304 {
305 #ifndef MB_AD_KD_TREE_USE_SINGLE_TAG
306  ErrorCode r1, r2;
307  r1 = moab()->tag_set_data( planeTag, &entity, 1, &plane.coord );
308  r2 = moab()->tag_set_data( axisTag, &entity, 1, &plane.norm );
309  return MB_SUCCESS == r1 ? r2 : r1;
310 #elif defined( MB_AD_KD_TREE_USE_TWO_DOUBLE_TAG )
311  double values[2] = { plane.coord, static_cast< double >( plane.norm ) };
312  return moab()->tag_set_data( planeTag, &entity, 1, values );
313 #else
314  return moab()->tag_set_data( planeTag, &entity, 1, &plane );
315 #endif
316 }

References axisTag, moab::AdaptiveKDTree::Plane::coord, ErrorCode, MB_SUCCESS, moab::Tree::moab(), moab::AdaptiveKDTree::Plane::norm, planeTag, and moab::Interface::tag_set_data().

Referenced by split_leaf().

◆ sphere_intersect_triangles()

ErrorCode moab::AdaptiveKDTree::sphere_intersect_triangles ( EntityHandle  tree_root,
const double  center[3],
double  radius,
std::vector< EntityHandle > &  triangles 
)

Definition at line 1801 of file AdaptiveKDTree.cpp.

1805 {
1806  ErrorCode rval;
1807  std::vector< EntityHandle > leaves;
1808  const CartVect from( center );
1809  CartVect closest_pt;
1810  const EntityHandle* conn;
1811  CartVect coords[3];
1812  int conn_len;
1813 
1814  // get leaves of tree that intersect sphere
1815  assert( tree_root );
1816  rval = distance_search( center, rad, leaves, 1.0e-10, 1.0e-6, NULL, NULL, &tree_root );
1817  if( MB_SUCCESS != rval ) return rval;
1818 
1819  // search each leaf for triangles intersecting sphere
1820  for( unsigned i = 0; i < leaves.size(); ++i )
1821  {
1822  Range tris;
1823  rval = moab()->get_entities_by_type( leaves[i], MBTRI, tris );
1824  if( MB_SUCCESS != rval ) return rval;
1825 
1826  for( Range::iterator j = tris.begin(); j != tris.end(); ++j )
1827  {
1828  rval = moab()->get_connectivity( *j, conn, conn_len );
1829  if( MB_SUCCESS != rval ) return rval;
1830  rval = moab()->get_coords( conn, 3, coords[0].array() );
1831  if( MB_SUCCESS != rval ) return rval;
1832  GeomUtil::closest_location_on_tri( from, coords, closest_pt );
1833  closest_pt -= from;
1834  if( ( closest_pt % closest_pt ) <= ( rad * rad ) ) triangles.push_back( *j );
1835  }
1836  }
1837 
1838  // remove duplicates from triangle list
1839  std::sort( triangles.begin(), triangles.end() );
1840  triangles.erase( std::unique( triangles.begin(), triangles.end() ), triangles.end() );
1841  return MB_SUCCESS;
1842 }

References moab::Range::begin(), center(), moab::GeomUtil::closest_location_on_tri(), distance_search(), moab::Range::end(), ErrorCode, moab::Interface::get_connectivity(), moab::Interface::get_coords(), moab::Interface::get_entities_by_type(), MB_SUCCESS, MBTRI, and moab::Tree::moab().

◆ split_leaf() [1/4]

ErrorCode moab::AdaptiveKDTree::split_leaf ( AdaptiveKDTreeIter leaf,
Plane  plane 
)

Split leaf of tree Updates iterator location to point to first new leaf node.

Definition at line 371 of file AdaptiveKDTree.cpp.

372 {
373  EntityHandle left, right;
374  return split_leaf( leaf, plane, left, right );
375 }

Referenced by build_tree(), and split_leaf().

◆ split_leaf() [2/4]

ErrorCode moab::AdaptiveKDTree::split_leaf ( AdaptiveKDTreeIter leaf,
Plane  plane,
const Range left_entities,
const Range right_entities 
)

Split leaf of tree Updates iterator location to point to first new leaf node.

Definition at line 377 of file AdaptiveKDTree.cpp.

381 {
382  EntityHandle left, right, parent = leaf.handle();
383  ErrorCode rval = split_leaf( leaf, plane, left, right );
384  if( MB_SUCCESS != rval ) return rval;
385 
386  if( MB_SUCCESS == moab()->add_entities( left, left_entities ) &&
387  MB_SUCCESS == moab()->add_entities( right, right_entities ) &&
388  MB_SUCCESS == moab()->clear_meshset( &parent, 1 ) )
389  return MB_SUCCESS;
390 
391  moab()->remove_child_meshset( parent, left );
392  moab()->remove_child_meshset( parent, right );
393  EntityHandle children[] = { left, right };
394  moab()->delete_entities( children, 2 );
395  return MB_FAILURE;
396 }

References children, moab::Interface::delete_entities(), ErrorCode, moab::AdaptiveKDTreeIter::handle(), MB_SUCCESS, moab::Tree::moab(), moab::Interface::remove_child_meshset(), and split_leaf().

◆ split_leaf() [3/4]

ErrorCode moab::AdaptiveKDTree::split_leaf ( AdaptiveKDTreeIter leaf,
Plane  plane,
const std::vector< EntityHandle > &  left_entities,
const std::vector< EntityHandle > &  right_entities 
)

Split leaf of tree Updates iterator location to point to first new leaf node.

Definition at line 398 of file AdaptiveKDTree.cpp.

402 {
403  EntityHandle left, right, parent = leaf.handle();
404  ErrorCode rval = split_leaf( leaf, plane, left, right );
405  if( MB_SUCCESS != rval ) return rval;
406 
407  if( MB_SUCCESS == moab()->add_entities( left, &left_entities[0], left_entities.size() ) &&
408  MB_SUCCESS == moab()->add_entities( right, &right_entities[0], right_entities.size() ) &&
409  MB_SUCCESS == moab()->clear_meshset( &parent, 1 ) )
410  return MB_SUCCESS;
411 
412  moab()->remove_child_meshset( parent, left );
413  moab()->remove_child_meshset( parent, right );
414  EntityHandle children[] = { left, right };
415  moab()->delete_entities( children, 2 );
416  return MB_FAILURE;
417 }

References children, moab::Interface::delete_entities(), ErrorCode, moab::AdaptiveKDTreeIter::handle(), MB_SUCCESS, moab::Tree::moab(), moab::Interface::remove_child_meshset(), and split_leaf().

◆ split_leaf() [4/4]

ErrorCode moab::AdaptiveKDTree::split_leaf ( AdaptiveKDTreeIter leaf,
Plane  plane,
EntityHandle left_child,
EntityHandle right_child 
)

Split leaf of tree Updates iterator location to point to first new leaf node.

Definition at line 344 of file AdaptiveKDTree.cpp.

345 {
346  ErrorCode rval;
347 
348  rval = moab()->create_meshset( meshsetFlags, left );
349  if( MB_SUCCESS != rval ) return rval;
350 
351  rval = moab()->create_meshset( meshsetFlags, right );
352  if( MB_SUCCESS != rval )
353  {
354  moab()->delete_entities( &left, 1 );
355  return rval;
356  }
357 
358  if( MB_SUCCESS != set_split_plane( leaf.handle(), plane ) ||
359  MB_SUCCESS != moab()->add_child_meshset( leaf.handle(), left ) ||
360  MB_SUCCESS != moab()->add_child_meshset( leaf.handle(), right ) ||
361  MB_SUCCESS != leaf.step_to_first_leaf( AdaptiveKDTreeIter::LEFT ) )
362  {
363  EntityHandle children[] = { left, right };
364  moab()->delete_entities( children, 2 );
365  return MB_FAILURE;
366  }
367 
368  return MB_SUCCESS;
369 }

References children, moab::Interface::create_meshset(), moab::Interface::delete_entities(), ErrorCode, moab::AdaptiveKDTreeIter::handle(), moab::AdaptiveKDTreeIter::LEFT, MB_SUCCESS, moab::Tree::meshsetFlags, moab::Tree::moab(), set_split_plane(), and moab::AdaptiveKDTreeIter::step_to_first_leaf().

Friends And Related Function Documentation

◆ AdaptiveKDTreeIter

friend class AdaptiveKDTreeIter
friend

Definition at line 264 of file AdaptiveKDTree.hpp.

Member Data Documentation

◆ axisTag

Tag moab::AdaptiveKDTree::axisTag
private

Definition at line 331 of file AdaptiveKDTree.hpp.

Referenced by get_split_plane(), init(), make_tag(), and set_split_plane().

◆ planeSet

CandidatePlaneSet moab::AdaptiveKDTree::planeSet
private

Definition at line 335 of file AdaptiveKDTree.hpp.

Referenced by build_tree(), and parse_options().

◆ planeTag

Tag moab::AdaptiveKDTree::planeTag
private

Definition at line 331 of file AdaptiveKDTree.hpp.

Referenced by get_split_plane(), init(), make_tag(), and set_split_plane().

◆ radius

double moab::AdaptiveKDTree::radius
private

Definition at line 338 of file AdaptiveKDTree.hpp.

Referenced by build_tree(), intersect_children_with_elems(), and parse_options().

◆ spherical

bool moab::AdaptiveKDTree::spherical
private

Definition at line 337 of file AdaptiveKDTree.hpp.

Referenced by build_tree(), intersect_children_with_elems(), and parse_options().

◆ splitsPerDir

unsigned moab::AdaptiveKDTree::splitsPerDir
private

Definition at line 333 of file AdaptiveKDTree.hpp.

Referenced by build_tree(), and parse_options().

◆ treeName

const char * moab::AdaptiveKDTree::treeName = "AKDTree"
staticprivate

Definition at line 329 of file AdaptiveKDTree.hpp.

Referenced by AdaptiveKDTree(), and init().


The documentation for this class was generated from the following files: