[MOAB-dev] r1527 - MOAB/trunk
kraftche at mcs.anl.gov
kraftche at mcs.anl.gov
Tue Jan 15 10:14:09 CST 2008
Author: kraftche
Date: 2008-01-15 10:14:09 -0600 (Tue, 15 Jan 2008)
New Revision: 1527
Modified:
MOAB/trunk/MBAdaptiveKDTree.cpp
MOAB/trunk/MBAdaptiveKDTree.hpp
MOAB/trunk/kd_tree_test.cpp
Log:
Add a few features to MBAdaptiveKDTree:
o Make MBAdaptiveKDTreeIter bidirectional.
o Add function to get iterator for leaf containing a point.
Modified: MOAB/trunk/MBAdaptiveKDTree.cpp
===================================================================
--- MOAB/trunk/MBAdaptiveKDTree.cpp 2008-01-14 23:35:35 UTC (rev 1526)
+++ MOAB/trunk/MBAdaptiveKDTree.cpp 2008-01-15 16:14:09 UTC (rev 1527)
@@ -264,7 +264,7 @@
}
MBErrorCode MBAdaptiveKDTree::get_tree_iterator( MBEntityHandle root,
- MBAdaptiveKDTreeIter& iter )
+ MBAdaptiveKDTreeIter& iter )
{
double box[6];
MBErrorCode rval = moab()->tag_get_data( rootTag, &root, 1, box );
@@ -274,12 +274,23 @@
return get_sub_tree_iterator( root, box, box+3, iter );
}
+MBErrorCode MBAdaptiveKDTree::get_last_iterator( MBEntityHandle root,
+ MBAdaptiveKDTreeIter& iter )
+{
+ double box[6];
+ MBErrorCode rval = moab()->tag_get_data( rootTag, &root, 1, box );
+ if (MB_SUCCESS != rval)
+ return rval;
+
+ return iter.initialize( this, root, box, box+3, MBAdaptiveKDTreeIter::RIGHT );
+}
+
MBErrorCode MBAdaptiveKDTree::get_sub_tree_iterator( MBEntityHandle root,
- const double min[3],
- const double max[3],
- MBAdaptiveKDTreeIter& result )
+ const double min[3],
+ const double max[3],
+ MBAdaptiveKDTreeIter& result )
{
- return result.initialize( this, root, min, max );
+ return result.initialize( this, root, min, max, MBAdaptiveKDTreeIter::LEFT );
}
MBErrorCode MBAdaptiveKDTree::split_leaf( MBAdaptiveKDTreeIter& leaf,
@@ -302,7 +313,7 @@
if (MB_SUCCESS != set_split_plane( leaf.handle(), plane ) ||
MB_SUCCESS != moab()->add_child_meshset( leaf.handle(), left ) ||
MB_SUCCESS != moab()->add_child_meshset( leaf.handle(), right) ||
- MB_SUCCESS != leaf.step_to_first_leaf()) {
+ MB_SUCCESS != leaf.step_to_first_leaf(MBAdaptiveKDTreeIter::LEFT)) {
MBEntityHandle children[] = { left, right };
moab()->delete_entities( children, 2 );
return MB_FAILURE;
@@ -427,7 +438,8 @@
MBErrorCode MBAdaptiveKDTreeIter::initialize( MBAdaptiveKDTree* tool,
MBEntityHandle root,
const double box_min[3],
- const double box_max[3] )
+ const double box_max[3],
+ Direction direction )
{
mStack.clear();
treeTool = tool;
@@ -438,13 +450,14 @@
mBox[BMAX][1] = box_max[1];
mBox[BMAX][2] = box_max[2];
mStack.push_back( StackObj(root,0) );
- return step_to_first_leaf();
+ return step_to_first_leaf( direction );
}
-MBErrorCode MBAdaptiveKDTreeIter::step_to_first_leaf()
+MBErrorCode MBAdaptiveKDTreeIter::step_to_first_leaf( Direction direction )
{
MBErrorCode rval;
MBAdaptiveKDTree::Plane plane;
+ const Direction opposite = static_cast<Direction>(1-direction);
for (;;) {
childVect.clear();
@@ -458,17 +471,18 @@
if (MB_SUCCESS != rval)
return rval;
- mStack.push_back( StackObj(childVect[0],mBox[BMAX][plane.norm]) );
- mBox[BMAX][plane.norm] = plane.coord;
+ mStack.push_back( StackObj(childVect[direction],mBox[opposite][plane.norm]) );
+ mBox[opposite][plane.norm] = plane.coord;
}
return MB_SUCCESS;
}
-MBErrorCode MBAdaptiveKDTreeIter::step()
+MBErrorCode MBAdaptiveKDTreeIter::step( Direction direction )
{
StackObj node, parent;
MBErrorCode rval;
MBAdaptiveKDTree::Plane plane;
+ const Direction opposite = static_cast<Direction>(1-direction);
// If stack is empty, then either this iterator is uninitialized
// or we reached the end of the iteration (and return
@@ -494,23 +508,23 @@
return rval;
// If we're at the left child
- if (childVect[0] == node.entity) {
+ if (childVect[opposite] == node.entity) {
// change from box of left child to box of parent
- mBox[BMAX][plane.norm] = node.coord;
+ mBox[direction][plane.norm] = node.coord;
// push right child on stack
- node.entity = childVect[1];
- node.coord = mBox[BMIN][plane.norm];
+ node.entity = childVect[direction];
+ node.coord = mBox[opposite][plane.norm];
mStack.push_back( node );
// change from box of parent to box of right child
- mBox[BMIN][plane.norm] = plane.coord;
+ mBox[opposite][plane.norm] = plane.coord;
// descend to left-most leaf of the right child
- return step_to_first_leaf();
+ return step_to_first_leaf(opposite);
}
// The current node is the right child of the parent,
// continue up the tree.
- assert( childVect[1] == node.entity );
- mBox[BMIN][plane.norm] = node.coord;
+ assert( childVect[direction] == node.entity );
+ mBox[opposite][plane.norm] = node.coord;
node = parent;
mStack.pop_back();
}
@@ -1090,7 +1104,7 @@
return rval;
MBAdaptiveKDTreeIter iter;
- iter.initialize( this, root_set_out, bmin.array(), bmax.array() );
+ iter.initialize( this, root_set_out, bmin.array(), bmax.array(), MBAdaptiveKDTreeIter::LEFT );
for (;;) {
@@ -1199,6 +1213,55 @@
return MB_SUCCESS;
}
+MBErrorCode MBAdaptiveKDTree::leaf_containing_point( MBEntityHandle root,
+ const double point[3],
+ MBAdaptiveKDTreeIter& result )
+{
+ // get bounding box of tree
+ MBErrorCode rval = moab()->tag_get_data( rootTag, &root, 1, result.mBox );
+ if (MB_SUCCESS != rval)
+ return rval;
+
+ // test that point is inside tree
+ if (point[0] < result.box_min()[0] || point[0] > result.box_max()[0] ||
+ point[1] < result.box_min()[1] || point[1] > result.box_max()[1] ||
+ point[2] < result.box_min()[2] || point[2] > result.box_max()[2])
+ return MB_ENTITY_NOT_FOUND;
+
+ // initialize iterator at tree root
+ result.treeTool = this;
+ result.mStack.clear();
+ result.mStack.push_back( MBAdaptiveKDTreeIter::StackObj(root,0) );
+
+ // loop until we reach a leaf
+ MBAdaptiveKDTree::Plane plane;
+ for(;;) {
+ // get children
+ result.childVect.clear();
+ rval = moab()->get_child_meshsets( result.handle(), result.childVect );
+ if (MB_SUCCESS != rval)
+ return rval;
+
+ // if no children, then at leaf (done)
+ if (result.childVect.empty())
+ break;
+
+ // get split plane
+ rval = get_split_plane( result.handle(), plane );
+ if (MB_SUCCESS != rval)
+ return rval;
+
+ // step iterator to appropriate child
+ // idx: 0->left, 1->right
+ const int idx = (point[plane.norm] > plane.coord);
+ result.mStack.push_back( MBAdaptiveKDTreeIter::StackObj( result.childVect[idx],
+ result.mBox[1-idx][plane.norm] ) );
+ result.mBox[1-idx][plane.norm] = plane.coord;
+ }
+
+ return MB_SUCCESS;
+}
+
struct NodeDistance {
MBEntityHandle handle;
MBCartVect dist; // from_point - closest_point_on_box
Modified: MOAB/trunk/MBAdaptiveKDTree.hpp
===================================================================
--- MOAB/trunk/MBAdaptiveKDTree.hpp 2008-01-14 23:35:35 UTC (rev 1526)
+++ MOAB/trunk/MBAdaptiveKDTree.hpp 2008-01-15 16:14:09 UTC (rev 1527)
@@ -85,6 +85,10 @@
//! Get iterator for tree
MBErrorCode get_tree_iterator( MBEntityHandle tree_root,
MBAdaptiveKDTreeIter& result );
+
+ //! Get iterator at right-most ('last') leaf.
+ MBErrorCode get_last_iterator( MBEntityHandle tree_root,
+ MBAdaptiveKDTreeIter& result );
//! Get iterator for tree or subtree
MBErrorCode get_sub_tree_iterator( MBEntityHandle tree_root,
@@ -183,6 +187,14 @@
const double point[3],
MBEntityHandle& leaf_out );
+ //! Get iterator at leaf containing input position.
+ //!
+ //! Returns MB_ENTITY_NOT_FOUND if point is not within
+ //! bounding box of tree.
+ MBErrorCode leaf_containing_point( MBEntityHandle tree_root,
+ const double xyz[3],
+ MBAdaptiveKDTreeIter& result );
+
//! Find all leaves within a given distance from point in space.
MBErrorCode leaves_within_distance( MBEntityHandle tree_root,
const double from_point[3],
@@ -195,6 +207,10 @@
//! Iterate over leaves of an adapative kD-tree
class MBAdaptiveKDTreeIter
{
+public:
+
+ enum Direction { LEFT = 0, RIGHT = 1 };
+
private:
struct StackObj {
@@ -213,7 +229,7 @@
//! Descend tree to left most leaf from current position
//! No-op if at leaf.
- MBErrorCode step_to_first_leaf();
+ MBErrorCode step_to_first_leaf( Direction direction );
friend class MBAdaptiveKDTree;
public:
@@ -223,7 +239,8 @@
MBErrorCode initialize( MBAdaptiveKDTree* tool,
MBEntityHandle root,
const double box_min[3],
- const double box_max[3] );
+ const double box_max[3],
+ Direction direction );
MBAdaptiveKDTree* tool() const
{ return treeTool; }
@@ -243,11 +260,26 @@
//! Get depth in tree. root is at depth of 1.
unsigned depth() const
{ return mStack.size(); }
+
+ //! Advance the iterator either left or right in the tree
+ //! Note: stepping past the end of the tree will invalidate
+ //! the iterator. It will *not* be work step the
+ //! other direction.
+ MBErrorCode step( Direction direction );
//! Advance to next leaf
//! Returns MB_ENTITY_NOT_FOUND if at end.
- MBErrorCode step();
+ //! Note: steping past the end of the tree will invalidate
+ //! the iterator. Calling back() will not work.
+ MBErrorCode step() { return step(RIGHT); }
+
+ //! Move back to previous leaf
+ //! Returns MB_ENTITY_NOT_FOUND if at beginning.
+ //! Note: steping past the start of the tree will invalidate
+ //! the iterator. Calling step() will not work.
+ MBErrorCode back() { return step(LEFT); }
+
//! Return the side of the box bounding this tree node
//! that is shared with the immediately adjacent sibling
//! (the tree node that shares a common parent node with
Modified: MOAB/trunk/kd_tree_test.cpp
===================================================================
--- MOAB/trunk/kd_tree_test.cpp 2008-01-14 23:35:35 UTC (rev 1526)
+++ MOAB/trunk/kd_tree_test.cpp 2008-01-15 16:14:09 UTC (rev 1527)
@@ -5,10 +5,15 @@
#include <assert.h>
#include <float.h>
+#include "MBCartVect.hpp"
+
const unsigned INTERVALS = 4;
const unsigned DEPTH = 7; // 3*log2(INTERVALS)+1
const char* TAG_NAME = "TEST_DATA";
+void test_iterator_back( MBAdaptiveKDTree& tool, MBEntityHandle root );
+void test_point_search( MBAdaptiveKDTree& tool, MBEntityHandle root );
+
int main()
{
// Initialize MOAB & create tree tool
@@ -72,6 +77,9 @@
const int num_leaves = INTERVALS*INTERVALS*INTERVALS;
assert( num_leaves == counter );
+ test_iterator_back( tool, root );
+ test_point_search( tool, root );
+
// reduce tree depth to DEPTH-1 by merging adjacent leaf pairs,
// make new "leaf" have smaller of two data values on original pair
for (err = tool.get_tree_iterator( root, iter ); !err; err = iter.step()) {
@@ -193,3 +201,123 @@
return 0;
}
+
+void test_iterator_back( MBAdaptiveKDTree& tool, MBEntityHandle root )
+{
+ MBAdaptiveKDTreeIter iter;
+ MBErrorCode rval = tool.get_tree_iterator( root, iter );
+ assert( MB_SUCCESS == rval );
+
+ MBCartVect min( iter.box_min() );
+ MBCartVect max( iter.box_max() );
+ MBEntityHandle leaf = iter.handle();
+
+ // going back from first location should fail.
+ rval = iter.back();
+ assert( MB_ENTITY_NOT_FOUND == rval );
+ rval = tool.get_tree_iterator( root, iter );
+ assert( MB_SUCCESS == rval );
+
+ // make sure iterator is valid
+ assert( iter.box_min()[0] == min[0] && iter.box_min()[1] == min[1] && iter.box_min()[2] == min[2] );
+ assert( iter.box_max()[0] == max[0] && iter.box_max()[1] == max[1] && iter.box_max()[2] == max[2] );
+ assert( iter.handle() == leaf );
+
+ while (MB_SUCCESS == iter.step()) {
+ // Get values at current iterator location
+ MBCartVect next_min( iter.box_min() );
+ MBCartVect next_max( iter.box_max() );
+ MBEntityHandle next_leaf = iter.handle();
+
+ // step back to previous location
+ rval = iter.back();
+ assert( MB_SUCCESS == rval );
+
+ // check expected values for previous location
+ assert( iter.box_min()[0] == min[0] && iter.box_min()[1] == min[1] && iter.box_min()[2] == min[2] );
+ assert( iter.box_max()[0] == max[0] && iter.box_max()[1] == max[1] && iter.box_max()[2] == max[2] );
+ assert( iter.handle() == leaf );
+
+ // advance iterator to 'current' location
+ rval = iter.step();
+ assert( MB_SUCCESS == rval );
+
+ // check that iterator values are correct
+ assert( iter.box_min()[0] == next_min[0] && iter.box_min()[1] == next_min[1] && iter.box_min()[2] == next_min[2] );
+ assert( iter.box_max()[0] == next_max[0] && iter.box_max()[1] == next_max[1] && iter.box_max()[2] == next_max[2] );
+ assert( iter.handle() == next_leaf );
+
+ // store values for next iteration
+ min = next_min;
+ max = next_max;
+ leaf = next_leaf;
+ }
+}
+
+void test_point_search( MBAdaptiveKDTree& tool, MBEntityHandle root )
+{
+ MBErrorCode rval;
+ MBEntityHandle leaf;
+ MBAdaptiveKDTreeIter iter, iter2;
+
+ // points to test
+ MBCartVect left( 0.5 );
+ MBCartVect right( MBCartVect(INTERVALS) - left );
+
+ // compare leaf search to iterator search
+ rval = tool.leaf_containing_point( root, left.array(), leaf );
+ assert( MB_SUCCESS == rval );
+ rval = tool.leaf_containing_point( root, left.array(), iter );
+ assert( MB_SUCCESS == rval );
+ assert( iter.handle() == leaf );
+
+ // iterator should be at 'first' leaf
+ rval = tool.get_tree_iterator( root, iter2 );
+ assert( MB_SUCCESS == rval );
+ for (;;) {
+ assert( iter.handle() == iter2.handle() );
+ assert( iter.depth() == iter2.depth() );
+ assert( iter.box_min()[0] == iter2.box_min()[0] );
+ assert( iter.box_min()[1] == iter2.box_min()[1] );
+ assert( iter.box_min()[2] == iter2.box_min()[2] );
+ assert( iter.box_max()[0] == iter2.box_max()[0] );
+ assert( iter.box_max()[1] == iter2.box_max()[1] );
+ assert( iter.box_max()[2] == iter2.box_max()[2] );
+
+ rval = iter2.step();
+ if (MB_SUCCESS != rval)
+ break;
+ rval = iter.step();
+ assert( MB_SUCCESS == rval );
+ }
+
+ // compare leaf search to iterator search
+ rval = tool.leaf_containing_point( root, right.array(), leaf );
+ assert( MB_SUCCESS == rval );
+ rval = tool.leaf_containing_point( root, right.array(), iter );
+ assert( MB_SUCCESS == rval );
+ assert( iter.handle() == leaf );
+
+ // iterator should be at 'last' leaf
+ rval = tool.get_last_iterator( root, iter2 );
+ assert( MB_SUCCESS == rval );
+ for (;;) {
+ assert( iter.handle() == iter2.handle() );
+ assert( iter.depth() == iter2.depth() );
+ assert( iter.box_min()[0] == iter2.box_min()[0] );
+ assert( iter.box_min()[1] == iter2.box_min()[1] );
+ assert( iter.box_min()[2] == iter2.box_min()[2] );
+ assert( iter.box_max()[0] == iter2.box_max()[0] );
+ assert( iter.box_max()[1] == iter2.box_max()[1] );
+ assert( iter.box_max()[2] == iter2.box_max()[2] );
+
+ rval = iter2.back();
+ if (MB_SUCCESS != rval)
+ break;
+ rval = iter.back();
+ assert( MB_SUCCESS == rval );
+ }
+}
+
+
+
More information about the moab-dev
mailing list