Class for constructing and querying Hierarchical Oriented Bounding Box trees. More...

#include <OrientedBoxTreeTool.hpp>

Collaboration diagram for moab::OrientedBoxTreeTool:

Classes
class	IntRegCtxt
	Default/Base class to provide a context for registering intersections. More...
class	Op
	Implement this and pass instance to preorder_traverse. More...
struct	SetData
struct	Settings
	Misc. knobs controlling tree subdivision. More...
class	TrvStats
	Traversal statistics structure. More...
Public Types
typedef std::pair< const double , const double >	IntersectSearchWindow
	This provides the search range for ray intersections, measured relative to the origin of the ray.
Public Member Functions
	OrientedBoxTreeTool (Interface i, const char tag_name=0, bool destroy_created_trees=false)
	~OrientedBoxTreeTool ()
ErrorCode	build (const Range &entities, EntityHandle &set_handle_out, const Settings *settings=0)
	Build oriented bounding box tree.
ErrorCode	join_trees (const Range &tree_roots, EntityHandle &root_set_out, const Settings *settings=0)
	Build a tree of sets, where each set contains triangles.
ErrorCode	ray_intersect_triangles (std::vector< double > &distances_out, std::vector< EntityHandle > &facets_out, EntityHandle root_set, double tolerance, const double ray_point[3], const double unit_ray_dir[3], const double ray_length=0, TrvStats accum=0)
	Intersect a ray with the triangles contained within the tree.
ErrorCode	ray_intersect_boxes (Range &boxes_out, EntityHandle root_set, double tolerance, const double ray_point[3], const double unit_ray_dir[3], const double ray_length=0, TrvStats accum=0)
	Intersect ray with tree.
ErrorCode	ray_intersect_triangles (std::vector< double > &intersection_distances_out, std::vector< EntityHandle > &intersection_facets_out, const Range &leaf_boxes_containing_tris, double tolerance, const double ray_point[3], const double unit_ray_dir[3], const double ray_length=0, unsigned int raytri_test_count=0)
	Intersect ray with triangles contained in passed MBENTITYSETs.
ErrorCode	ray_intersect_sets (std::vector< double > &distances_out, std::vector< EntityHandle > &sets_out, std::vector< EntityHandle > &facets_out, EntityHandle root_set, double tolerance, const double ray_point[3], const double unit_ray_dir[3], IntersectSearchWindow &search_win, IntRegCtxt &register_intersection, TrvStats *accum=0)
	Intersect a ray with the triangles contained within the tree.
ErrorCode	ray_intersect_sets (std::vector< double > &distances_out, std::vector< EntityHandle > &sets_out, std::vector< EntityHandle > &facets_out, EntityHandle root_set, double tolerance, const double ray_point[3], const double unit_ray_dir[3], const double ray_length=0, TrvStats accum=0)
ErrorCode	ray_intersect_sets (EntityHandle root_set, double tolerance, const double ray_point[3], const double unit_ray_dir[3], IntersectSearchWindow &search_win, IntRegCtxt &register_intersection, TrvStats *accum=0)
ErrorCode	closest_to_location (const double point, EntityHandle tree_root, double point_out, EntityHandle &facet_out, EntityHandle set_out=0, TrvStats accum=0)
	Find closest surface, facet in surface, and location on facet.
ErrorCode	closest_to_location (const double point, EntityHandle tree_root, double tolerance, std::vector< EntityHandle > &facets_out, std::vector< EntityHandle > sets_out=0, TrvStats *accum=0)
	Find closest facet(s) to input position.
ErrorCode	sphere_intersect_triangles (const double center, double radius, EntityHandle tree_root, std::vector< EntityHandle > &facets_out, std::vector< EntityHandle > sets_out=0, TrvStats *accum=0)
	Find facets intersected by a sphere.
ErrorCode	get_close_tris (CartVect int_pt, double tol, const EntityHandle rootSet, const EntityHandle geomVol, const Tag *senseTag, std::vector< EntityHandle > &close_tris, std::vector< int > &close_senses)
ErrorCode	box (EntityHandle node_set, double center[3], double axis1[3], double axis2[3], double axis3[3])
	Get oriented box at node in tree.
ErrorCode	delete_tree (EntityHandle root_set)
ErrorCode	remove_root (EntityHandle root_set)
	Remove obb tree root from the Oriented Box Tree Tool data structure.
void	print (EntityHandle tree_root_set, std::ostream &stream, bool list_contents=false, const char *id_tag_name=0)
	Print out tree.
ErrorCode	stats (EntityHandle tree_root_set, std::ostream &stream)
	Print tree statistics.
ErrorCode	stats (EntityHandle set, unsigned &entities_in_tree, double &root_volume, double &tot_node_volume, double &tot_to_root_volume, unsigned &tree_height, unsigned &node_count, unsigned &num_leaves)
	Get tree statistics.
ErrorCode	preorder_traverse (EntityHandle root_set, Op &operation, TrvStats *accum=0)
	Visitor pattern - do operation for each tree node.
Interface *	get_moab_instance () const
ErrorCode	box (EntityHandle node_set, OrientedBox &box)
	Get oriented box at node in tree.
Private Member Functions
ErrorCode	build_tree (const Range &entities, EntityHandle &set, int depth, const Settings &settings)
ErrorCode	build_sets (std::list< SetData > &sets, EntityHandle &node_set, int depth, const Settings &settings)
ErrorCode	recursive_stats (OrientedBoxTreeTool tool, Interface instance, EntityHandle set, int depth, StatData &data, unsigned &count_out, CartVect &dimensions_out)
Private Attributes
Interface *	instance
Tag	tagHandle
bool	cleanUpTrees
std::vector< EntityHandle >	createdTrees

Detailed Description

Class for constructing and querying Hierarchical Oriented Bounding Box trees.

Definition at line 42 of file OrientedBoxTreeTool.hpp.

Member Typedef Documentation

typedef std::pair< const double*, const double* > moab::OrientedBoxTreeTool::IntersectSearchWindow

This provides the search range for ray intersections, measured relative to the origin of the ray.

first: nonnegative limit for search second: negative limit for search

These are const double* so that the window is always defined by pointing to other quantities, but it is not posible to change those quantities via the window.

Definition at line 55 of file OrientedBoxTreeTool.hpp.

Constructor & Destructor Documentation

moab::OrientedBoxTreeTool::OrientedBoxTreeTool	(	Interface *	i,
		const char *	tag_name = `0`,
		bool	destroy_created_trees = `false`
	)

Definition at line 49 of file OrientedBoxTreeTool.cpp.

References moab::DEFAULT_TAG_NAME, ErrorCode, instance, MB_SUCCESS, moab::OrientedBox::tag_handle(), and tagHandle.

    : instance( i ), cleanUpTrees( destroy_created_trees )
{
    if( !tag_name ) tag_name = DEFAULT_TAG_NAME;
    ErrorCode rval = OrientedBox::tag_handle( tagHandle, instance, tag_name );
    if( MB_SUCCESS != rval ) tagHandle = 0;
}

moab::OrientedBoxTreeTool::~OrientedBoxTreeTool ( )

Definition at line 57 of file OrientedBoxTreeTool.cpp.

References cleanUpTrees, createdTrees, delete_tree(), ErrorCode, instance, MB_SUCCESS, moab::Interface::tag_get_by_ptr(), and tagHandle.

{
    if( !cleanUpTrees ) return;

    while( !createdTrees.empty() )
    {
        EntityHandle tree = createdTrees.back();
        // make sure this is a tree (rather than some other, stale handle)
        const void* data_ptr = 0;
        ErrorCode rval       = instance->tag_get_by_ptr( tagHandle, &tree, 1, &data_ptr );
        if( MB_SUCCESS == rval ) rval = delete_tree( tree );
        if( MB_SUCCESS != rval ) createdTrees.pop_back();
    }
}

Member Function Documentation

ErrorCode moab::OrientedBoxTreeTool::box	(	EntityHandle	node_set,
		double	center[3],
		double	axis1[3],
		double	axis2[3],
		double	axis3[3]
	)

Get oriented box at node in tree.

Get the oriented box for a node in an oriented bounding box tree.

Definition at line 91 of file OrientedBoxTreeTool.cpp.

References moab::OrientedBox::center, ErrorCode, moab::CartVect::get(), and moab::OrientedBox::scaled_axis().

Referenced by closest_to_location(), FBiGeom_getEntBoundBox(), moab::GeomTopoTool::get_obb(), moab::FBEngine::getEntBoundBox(), main(), moab::TreeNodePrinter::print_geometry(), recursive_stats(), sphere_intersect_triangles(), moab::RayIntersector::visit(), and moab::RayIntersectSets::visit().

{
    OrientedBox obb;
    ErrorCode rval = this->box( set, obb );
    obb.center.get( center );
    obb.scaled_axis( 0 ).get( axis1 );
    obb.scaled_axis( 1 ).get( axis2 );
    obb.scaled_axis( 2 ).get( axis3 );
    return rval;
}

ErrorCode moab::OrientedBoxTreeTool::box	(	EntityHandle	node_set,
		OrientedBox &	box
	)

Get oriented box at node in tree.

Get the oriented box for a node in an oriented bounding box tree.

NOTE: This function is provided for internal MOAB use only. The definition of OrientedBox is not available as a part of the MOAB API

Definition at line 86 of file OrientedBoxTreeTool.cpp.

References instance, moab::Interface::tag_get_data(), and tagHandle.

{
    return instance->tag_get_data( tagHandle, &set, 1, &obb );
}

ErrorCode moab::OrientedBoxTreeTool::build	(	const Range &	entities,
		EntityHandle &	set_handle_out,
		const Settings *	settings = `0`
	)

Build oriented bounding box tree.

Build an oriented bounding box tree.

Parameters:

entities	A list of either vertices or 2-D elements (not both) for which to build a tree.
set_handle_out	A handle for the entity set representing the root of the tree.

Definition at line 115 of file OrientedBoxTreeTool.cpp.

References moab::Range::all_of_dimension(), build_tree(), MB_TYPE_OUT_OF_RANGE, and moab::OrientedBoxTreeTool::Settings::valid().

Referenced by moab::GeomTopoTool::construct_obb_tree(), and main().

{
    if( !entities.all_of_dimension( 2 ) ) return MB_TYPE_OUT_OF_RANGE;
    if( settings && !settings->valid() ) return MB_FAILURE;

    return build_tree( entities, set_handle_out, 0, settings ? *settings : Settings() );
}

ErrorCode moab::OrientedBoxTreeTool::build_sets	(	std::list< SetData > &	sets,
		EntityHandle &	node_set,
		int	depth,
		const Settings &	settings
	)		`[private]`

Definition at line 340 of file OrientedBoxTreeTool.cpp.

References moab::Interface::add_child_meshset(), child, moab::Range::clear(), moab::OrientedBox::compute_from_covariance_data(), moab::Interface::create_meshset(), delete_tree(), ErrorCode, moab::Interface::get_adjacencies(), moab::Interface::get_entities_by_dimension(), instance, MB_SUCCESS, moab::OrientedBoxTreeTool::Settings::set_options, moab::split_sets(), moab::Interface::tag_set_data(), tagHandle, and moab::Interface::UNION.

Referenced by join_trees().

{
    ErrorCode rval;
    int count = sets.size();
    if( 0 == count ) return MB_FAILURE;

    // calculate box
    OrientedBox obox;

    // make vector go out of scope when done, so memory is released
    {
        Range elems;
        std::vector< OrientedBox::CovarienceData > data( sets.size() );
        data.clear();
        for( std::list< SetData >::iterator i = sets.begin(); i != sets.end(); ++i )
        {
            data.push_back( i->box_data );
            rval = instance->get_entities_by_dimension( i->handle, 2, elems, true );
            if( MB_SUCCESS != rval ) return rval;
        }

        Range points;
        rval = instance->get_adjacencies( elems, 0, false, points, Interface::UNION );
        if( MB_SUCCESS != rval ) return rval;

        rval = OrientedBox::compute_from_covariance_data( obox, instance, &data[0], data.size(), points );
        if( MB_SUCCESS != rval ) return rval;
    }

    // If only one set in list...
    if( count == 1 )
    {
        node_set = sets.front().handle;
        return instance->tag_set_data( tagHandle, &node_set, 1, &obox );
    }

    // create an entity set for the tree node
    rval = instance->create_meshset( settings.set_options, node_set );
    if( MB_SUCCESS != rval ) return rval;

    rval = instance->tag_set_data( tagHandle, &node_set, 1, &obox );
    if( MB_SUCCESS != rval )
    {
        delete_tree( node_set );
        return rval;
    }

    double best_ratio = 2.0;
    std::list< SetData > best_left_list, best_right_list;
    for( int axis = 0; axis < 2; ++axis )
    {
        std::list< SetData > left_list, right_list;
        rval = split_sets( instance, obox, axis, sets, left_list, right_list );
        if( MB_SUCCESS != rval )
        {
            delete_tree( node_set );
            return rval;
        }

        double ratio = fabs( (double)right_list.size() - left_list.size() ) / sets.size();

        if( ratio < best_ratio )
        {
            best_ratio = ratio;
            best_left_list.swap( left_list );
            best_right_list.swap( right_list );
        }
    }

    // We must subdivide the list of sets, because we want to guarantee that
    // there is a node in the tree corresponding to each of the sets.  So if
    // we couldn't find a usable split plane, just split them in an arbitrary
    // manner.
    if( best_left_list.empty() || best_right_list.empty() )
    {
        best_left_list.clear();
        best_right_list.clear();
        std::list< SetData >* lists[2] = { &best_left_list, &best_right_list };
        int i                          = 0;
        while( !sets.empty() )
        {
            lists[i]->push_back( sets.front() );
            sets.pop_front();
            i = 1 - i;
        }
    }
    else
    {
        sets.clear();  // release memory before recursion
    }

    // Create child sets

    EntityHandle child = 0;

    rval = build_sets( best_left_list, child, depth + 1, settings );
    if( MB_SUCCESS != rval )
    {
        delete_tree( node_set );
        return rval;
    }
    rval = instance->add_child_meshset( node_set, child );
    if( MB_SUCCESS != rval )
    {
        delete_tree( node_set );
        delete_tree( child );
        return rval;
    }

    rval = build_sets( best_right_list, child, depth + 1, settings );
    if( MB_SUCCESS != rval )
    {
        delete_tree( node_set );
        return rval;
    }
    rval = instance->add_child_meshset( node_set, child );
    if( MB_SUCCESS != rval )
    {
        delete_tree( node_set );
        delete_tree( child );
        return rval;
    }

    return MB_SUCCESS;
}

ErrorCode moab::OrientedBoxTreeTool::build_tree	(	const Range &	entities,
		EntityHandle &	set,
		int	depth,
		const Settings &	settings
	)		`[private]`

Definition at line 204 of file OrientedBoxTreeTool.cpp.

References moab::Interface::add_child_meshset(), moab::Interface::add_entities(), moab::OrientedBoxTreeTool::Settings::best_split_ratio, child, moab::OrientedBox::compute_from_2d_cells(), moab::Interface::create_meshset(), createdTrees, delete_tree(), moab::Range::empty(), ErrorCode, instance, moab::OrientedBoxTreeTool::Settings::max_depth, moab::OrientedBoxTreeTool::Settings::max_leaf_entities, MB_SUCCESS, moab::OrientedBoxTreeTool::Settings::set_options, moab::Range::size(), moab::split_box(), moab::Range::swap(), moab::Interface::tag_set_data(), tagHandle, and moab::OrientedBoxTreeTool::Settings::worst_split_ratio.

Referenced by build().

{
    OrientedBox tmp_box;
    ErrorCode rval;

    if( entities.empty() )
    {
        Matrix3 axis;
        tmp_box = OrientedBox( axis, CartVect( 0. ) );
    }
    else
    {
        rval = OrientedBox::compute_from_2d_cells( tmp_box, instance, entities );
        if( MB_SUCCESS != rval ) return rval;
    }

    // create an entity set for the tree node
    rval = instance->create_meshset( settings.set_options, set );
    if( MB_SUCCESS != rval ) return rval;

    rval = instance->tag_set_data( tagHandle, &set, 1, &tmp_box );
    if( MB_SUCCESS != rval )
    {
        delete_tree( set );
        return rval;
    }

    // check if should create children
    bool leaf = true;
    ++depth;
    if( ( !settings.max_depth || depth < settings.max_depth ) &&
        entities.size() > (unsigned)settings.max_leaf_entities )
    {
        // try splitting with planes normal to each axis of the box
        // until we find an acceptable split
        double best_ratio = settings.worst_split_ratio;  // worst case ratio
        Range best_left_list, best_right_list;
        // Axes are sorted from shortest to longest, so search backwards
        for( int axis = 2; best_ratio > settings.best_split_ratio && axis >= 0; --axis )
        {
            Range left_list, right_list;

            rval = split_box( instance, tmp_box, axis, entities, left_list, right_list );
            if( MB_SUCCESS != rval )
            {
                delete_tree( set );
                return rval;
            }

            double ratio = fabs( (double)right_list.size() - left_list.size() ) / entities.size();

            if( ratio < best_ratio )
            {
                best_ratio = ratio;
                best_left_list.swap( left_list );
                best_right_list.swap( right_list );
            }
        }

        // create children
        if( !best_left_list.empty() )
        {
            EntityHandle child = 0;

            rval = build_tree( best_left_list, child, depth, settings );
            if( MB_SUCCESS != rval )
            {
                delete_tree( set );
                return rval;
            }
            rval = instance->add_child_meshset( set, child );
            if( MB_SUCCESS != rval )
            {
                delete_tree( set );
                delete_tree( child );
                return rval;
            }

            rval = build_tree( best_right_list, child, depth, settings );
            if( MB_SUCCESS != rval )
            {
                delete_tree( set );
                return rval;
            }
            rval = instance->add_child_meshset( set, child );
            if( MB_SUCCESS != rval )
            {
                delete_tree( set );
                delete_tree( child );
                return rval;
            }

            leaf = false;
        }
    }

    if( leaf )
    {
        rval = instance->add_entities( set, entities );
        if( MB_SUCCESS != rval )
        {
            delete_tree( set );
            return rval;
        }
    }

    createdTrees.push_back( set );
    return MB_SUCCESS;
}

ErrorCode moab::OrientedBoxTreeTool::closest_to_location	(	const double *	point,
		EntityHandle	tree_root,
		double *	point_out,
		EntityHandle &	facet_out,
		EntityHandle *	set_out = `0`,
		TrvStats *	accum = `0`
	)

Find closest surface, facet in surface, and location on facet.

Find the closest location in the tree to the specified location.

Parameters:

point	Location to search from
point_out	Closest location on closest facet
facet_out	Closest 2D element to input position
set_out	Set containing closest facet. 0 if tree was not constructed using 'set_build'

Definition at line 1118 of file OrientedBoxTreeTool.cpp.

References moab::Range::begin(), box(), children, moab::Range::clear(), moab::OrientedBox::closest_location_in_box(), moab::GeomUtil::closest_location_on_polygon(), moab::GeomUtil::closest_location_on_tri(), moab::Range::empty(), moab::Range::end(), moab::OrientedBoxTreeTool::TrvStats::end_traversal(), ErrorCode, moab::Range::front(), moab::CartVect::get(), moab::Interface::get_child_meshsets(), moab::Interface::get_connectivity(), moab::Interface::get_coords(), moab::Interface::get_entities_by_dimension(), moab::Interface::get_entities_by_type(), moab::OrientedBoxTreeTool::TrvStats::increment(), moab::OrientedBoxTreeTool::TrvStats::increment_leaf(), instance, MB_MULTIPLE_ENTITIES_FOUND, MB_SUCCESS, MBENTITYSET, and moab::Range::size().

Referenced by moab::GeomQueryTool::closest_to_location(), moab::GeomQueryTool::get_normal(), moab::FBEngine::getEntClosestPt(), moab::FBEngine::getEntNrmlXYZ(), moab::SmoothFace::project_to_facets_main(), and moab::GeomQueryTool::test_volume_boundary().

{
    ErrorCode rval;
    const CartVect loc( point );
    double smallest_dist_sqr = std::numeric_limits< double >::max();

    EntityHandle current_set = 0;
    Range sets;
    std::vector< EntityHandle > children( 2 );
    std::vector< double > coords;
    std::vector< OBBTreeCPFrame > stack;
    int max_depth = -1;

    stack.push_back( OBBTreeCPFrame( 0.0, root, current_set, 0 ) );

    while( !stack.empty() )
    {

        // pop from top of stack
        EntityHandle node = stack.back().node;
        double dist_sqr   = stack.back().dist_sqr;
        current_set       = stack.back().mset;
        int current_depth = stack.back().depth;
        stack.pop_back();

        // If current best result is closer than the box, skip this tree node.
        if( dist_sqr > smallest_dist_sqr ) continue;

        // increment traversal statistics.
        if( accum )
        {
            accum->increment( current_depth );
            max_depth = std::max( max_depth, current_depth );
        }

        // Check if this node has a set associated with it
        if( set_out && !current_set )
        {
            sets.clear();
            rval = instance->get_entities_by_type( node, MBENTITYSET, sets );
            if( MB_SUCCESS != rval ) return rval;
            if( !sets.empty() )
            {
                if( sets.size() != 1 ) return MB_MULTIPLE_ENTITIES_FOUND;
                current_set = sets.front();
            }
        }

        // Get child boxes
        children.clear();
        rval = instance->get_child_meshsets( node, children );
        if( MB_SUCCESS != rval ) return rval;

        // if not a leaf node
        if( !children.empty() )
        {
            if( children.size() != 2 ) return MB_MULTIPLE_ENTITIES_FOUND;

            // get boxes
            OrientedBox box1, box2;
            rval = box( children[0], box1 );
            if( MB_SUCCESS != rval ) return rval;
            rval = box( children[1], box2 );
            if( MB_SUCCESS != rval ) return rval;

            // get distance from each box
            CartVect pt1, pt2;
            box1.closest_location_in_box( loc, pt1 );
            box2.closest_location_in_box( loc, pt2 );
            pt1 -= loc;
            pt2 -= loc;
            const double dsqr1 = pt1 % pt1;
            const double dsqr2 = pt2 % pt2;

            // push children on tree such that closer one is on top
            if( dsqr1 < dsqr2 )
            {
                stack.push_back( OBBTreeCPFrame( dsqr2, children[1], current_set, current_depth + 1 ) );
                stack.push_back( OBBTreeCPFrame( dsqr1, children[0], current_set, current_depth + 1 ) );
            }
            else
            {
                stack.push_back( OBBTreeCPFrame( dsqr1, children[0], current_set, current_depth + 1 ) );
                stack.push_back( OBBTreeCPFrame( dsqr2, children[1], current_set, current_depth + 1 ) );
            }
        }
        else
        {  // LEAF NODE
            if( accum )
            {
                accum->increment_leaf( current_depth );
            }

            Range facets;
            rval = instance->get_entities_by_dimension( node, 2, facets );
            if( MB_SUCCESS != rval ) return rval;

            const EntityHandle* conn = NULL;
            int len                  = 0;
            CartVect tmp, diff;
            for( Range::iterator i = facets.begin(); i != facets.end(); ++i )
            {
                rval = instance->get_connectivity( *i, conn, len, true );
                if( MB_SUCCESS != rval ) return rval;

                coords.resize( 3 * len );
                rval = instance->get_coords( conn, len, &coords[0] );
                if( MB_SUCCESS != rval ) return rval;

                if( len == 3 )
                    GeomUtil::closest_location_on_tri( loc, (CartVect*)( &coords[0] ), tmp );
                else
                    GeomUtil::closest_location_on_polygon( loc, (CartVect*)( &coords[0] ), len, tmp );

                diff     = tmp - loc;
                dist_sqr = diff % diff;
                if( dist_sqr < smallest_dist_sqr )
                {
                    smallest_dist_sqr = dist_sqr;
                    facet_out         = *i;
                    tmp.get( point_out );
                    if( set_out ) *set_out = current_set;
                }
            }
        }  // LEAF NODE
    }

    if( accum )
    {
        accum->end_traversal( max_depth );
    }

    return MB_SUCCESS;
}

ErrorCode moab::OrientedBoxTreeTool::closest_to_location	(	const double *	point,
		EntityHandle	tree_root,
		double	tolerance,
		std::vector< EntityHandle > &	facets_out,
		std::vector< EntityHandle > *	sets_out = `0`,
		TrvStats *	accum = `0`
	)

Find closest facet(s) to input position.

Find the closest location(s) in the tree to the specified location.

Parameters:

point	Location to search from
facets_out	Closest 2D elements to input position are appended to this list
sets_out	If non-null, sets owning facets are appended to this list.

Definition at line 1258 of file OrientedBoxTreeTool.cpp.

References moab::Range::begin(), box(), children, moab::Range::clear(), moab::OrientedBox::closest_location_in_box(), moab::GeomUtil::closest_location_on_polygon(), moab::GeomUtil::closest_location_on_tri(), moab::Range::empty(), moab::Range::end(), moab::OrientedBoxTreeTool::TrvStats::end_traversal(), ErrorCode, moab::Interface::get_child_meshsets(), moab::Interface::get_connectivity(), moab::Interface::get_coords(), moab::Interface::get_entities_by_dimension(), moab::Interface::get_entities_by_type(), moab::OrientedBoxTreeTool::TrvStats::increment(), moab::OrientedBoxTreeTool::TrvStats::increment_leaf(), instance, MB_MULTIPLE_ENTITIES_FOUND, MB_SUCCESS, MBENTITYSET, and moab::Range::size().

{
    ErrorCode rval;
    const CartVect loc( point );
    double smallest_dist_sqr = std::numeric_limits< double >::max();
    double smallest_dist     = smallest_dist_sqr;

    EntityHandle current_set = 0;
    Range sets;
    std::vector< EntityHandle > children( 2 );
    std::vector< double > coords;
    std::vector< OBBTreeCPFrame > stack;
    int max_depth = -1;

    stack.push_back( OBBTreeCPFrame( 0.0, root, current_set, 0 ) );

    while( !stack.empty() )
    {

        // pop from top of stack
        EntityHandle node = stack.back().node;
        double dist_sqr   = stack.back().dist_sqr;
        current_set       = stack.back().mset;
        int current_depth = stack.back().depth;
        stack.pop_back();

        // If current best result is closer than the box, skip this tree node.
        if( dist_sqr > smallest_dist_sqr + tolerance ) continue;

        // increment traversal statistics.
        if( accum )
        {
            accum->increment( current_depth );
            max_depth = std::max( max_depth, current_depth );
        }

        // Check if this node has a set associated with it
        if( sets_out && !current_set )
        {
            sets.clear();
            rval = instance->get_entities_by_type( node, MBENTITYSET, sets );
            if( MB_SUCCESS != rval ) return rval;
            if( !sets.empty() )
            {
                if( sets.size() != 1 ) return MB_MULTIPLE_ENTITIES_FOUND;
                current_set = *sets.begin();
            }
        }

        // Get child boxes
        children.clear();
        rval = instance->get_child_meshsets( node, children );
        if( MB_SUCCESS != rval ) return rval;

        // if not a leaf node
        if( !children.empty() )
        {
            if( children.size() != 2 ) return MB_MULTIPLE_ENTITIES_FOUND;

            // get boxes
            OrientedBox box1, box2;
            rval = box( children[0], box1 );
            if( MB_SUCCESS != rval ) return rval;
            rval = box( children[1], box2 );
            if( MB_SUCCESS != rval ) return rval;

            // get distance from each box
            CartVect pt1, pt2;
            box1.closest_location_in_box( loc, pt1 );
            box2.closest_location_in_box( loc, pt2 );
            pt1 -= loc;
            pt2 -= loc;
            const double dsqr1 = pt1 % pt1;
            const double dsqr2 = pt2 % pt2;

            // push children on tree such that closer one is on top
            if( dsqr1 < dsqr2 )
            {
                stack.push_back( OBBTreeCPFrame( dsqr2, children[1], current_set, current_depth + 1 ) );
                stack.push_back( OBBTreeCPFrame( dsqr1, children[0], current_set, current_depth + 1 ) );
            }
            else
            {
                stack.push_back( OBBTreeCPFrame( dsqr1, children[0], current_set, current_depth + 1 ) );
                stack.push_back( OBBTreeCPFrame( dsqr2, children[1], current_set, current_depth + 1 ) );
            }
        }
        else
        {  // LEAF NODE
            if( accum )
            {
                accum->increment_leaf( current_depth );
            }

            Range facets;
            rval = instance->get_entities_by_dimension( node, 2, facets );
            if( MB_SUCCESS != rval ) return rval;

            const EntityHandle* conn = NULL;
            int len                  = 0;
            CartVect tmp, diff;
            for( Range::iterator i = facets.begin(); i != facets.end(); ++i )
            {
                rval = instance->get_connectivity( *i, conn, len, true );
                if( MB_SUCCESS != rval ) return rval;

                coords.resize( 3 * len );
                rval = instance->get_coords( conn, len, &coords[0] );
                if( MB_SUCCESS != rval ) return rval;

                if( len == 3 )
                    GeomUtil::closest_location_on_tri( loc, (CartVect*)( &coords[0] ), tmp );
                else
                    GeomUtil::closest_location_on_polygon( loc, (CartVect*)( &coords[0] ), len, tmp );

                diff     = tmp - loc;
                dist_sqr = diff % diff;
                if( dist_sqr < smallest_dist_sqr )
                {
                    if( 0.5 * dist_sqr < 0.5 * smallest_dist_sqr + tolerance * ( 0.5 * tolerance - smallest_dist ) )
                    {
                        facets_out.clear();
                        if( sets_out ) sets_out->clear();
                    }
                    smallest_dist_sqr = dist_sqr;
                    smallest_dist     = sqrt( smallest_dist_sqr );
                    // put closest result at start of lists
                    facets_out.push_back( *i );
                    std::swap( facets_out.front(), facets_out.back() );
                    if( sets_out )
                    {
                        sets_out->push_back( current_set );
                        std::swap( sets_out->front(), sets_out->back() );
                    }
                }
                else if( dist_sqr <= smallest_dist_sqr + tolerance * ( tolerance + 2 * smallest_dist ) )
                {
                    facets_out.push_back( *i );
                    if( sets_out ) sets_out->push_back( current_set );
                }
            }
        }  // LEAF NODE
    }

    if( accum )
    {
        accum->end_traversal( max_depth );
    }

    return MB_SUCCESS;
}

ErrorCode moab::OrientedBoxTreeTool::delete_tree ( EntityHandle root_set )

Definition at line 469 of file OrientedBoxTreeTool.cpp.

References children, createdTrees, moab::Interface::delete_entities(), ErrorCode, moab::Interface::get_child_meshsets(), instance, and MB_SUCCESS.

Referenced by build_sets(), build_tree(), and ~OrientedBoxTreeTool().

{
    std::vector< EntityHandle > children;
    ErrorCode rval = instance->get_child_meshsets( set, children, 0 );
    if( MB_SUCCESS != rval ) return rval;

    createdTrees.erase( std::remove( createdTrees.begin(), createdTrees.end(), set ), createdTrees.end() );
    children.insert( children.begin(), set );
    return instance->delete_entities( &children[0], children.size() );
}

ErrorCode moab::OrientedBoxTreeTool::get_close_tris	(	CartVect	int_pt,
		double	tol,
		const EntityHandle *	rootSet,
		const EntityHandle *	geomVol,
		const Tag *	senseTag,
		std::vector< EntityHandle > &	close_tris,
		std::vector< int > &	close_senses
	)

Definition at line 845 of file OrientedBoxTreeTool.cpp.

References moab::CartVect::array(), ErrorCode, get_moab_instance(), MB_SUCCESS, sphere_intersect_triangles(), and moab::Interface::tag_get_data().

Referenced by moab::GQT_IntRegCtxt::register_intersection().

{
    std::vector< EntityHandle > close_surfs;
    ErrorCode rval = sphere_intersect_triangles( int_pt.array(), tol, *rootSet, close_tris, &close_surfs );
    assert( MB_SUCCESS == rval );
    if( MB_SUCCESS != rval ) return rval;

    // for each surface, get the surf sense wrt parent volume
    close_senses.resize( close_surfs.size() );
    for( unsigned i = 0; i < close_surfs.size(); ++i )
    {
        EntityHandle vols[2];
        rval = get_moab_instance()->tag_get_data( *senseTag, &( close_surfs[i] ), 1, vols );
        assert( MB_SUCCESS == rval );
        if( MB_SUCCESS != rval ) return rval;
        if( vols[0] == vols[1] )
        {
            std::cerr << "error: surf has positive and negative sense wrt same volume" << std::endl;
            return MB_FAILURE;
        }
        if( *geomVol == vols[0] )
        {
            close_senses[i] = 1;
        }
        else if( *geomVol == vols[1] )
        {
            close_senses[i] = -1;
        }
        else
        {
            return MB_FAILURE;
        }
    }

    return MB_SUCCESS;
}

Interface* moab::OrientedBoxTreeTool::get_moab_instance ( ) const [inline]

Definition at line 525 of file OrientedBoxTreeTool.hpp.

References instance.

Referenced by get_close_tris(), moab::GQT_IntRegCtxt::register_intersection(), sphere_intersect_triangles(), and moab::GQT_IntRegCtxt::update_orient().

    {
        return instance;
    }

ErrorCode moab::OrientedBoxTreeTool::join_trees	(	const Range &	tree_roots,
		EntityHandle &	root_set_out,
		const Settings *	settings = `0`
	)

Build a tree of sets, where each set contains triangles.

Build a tree of sets. Each set must contain at least one triangle to define its geometry. Each passed set will become a leaf of the OBB tree. Settings controlling tree depth are ignored by this method. The tree will be as deep as it needs to be for each input set to be a leaf.

To build a tree representing the surfaces of a geometric volume, 1) Build an OBB tree for each surface using the 'build' method 2) Add each surface to the contents of the resulting OBB tree root set 3) Build a tree from all the surface OBB tree root sets using this method to get a combined tree for the volume.

Definition at line 123 of file OrientedBoxTreeTool.cpp.

References moab::Range::all_of_type(), moab::Range::begin(), moab::OrientedBoxTreeTool::SetData::box_data, build_sets(), moab::OrientedBox::covariance_data_from_tris(), createdTrees, moab::Range::empty(), moab::Range::end(), ErrorCode, moab::Interface::get_entities_by_dimension(), moab::OrientedBoxTreeTool::SetData::handle, instance, MB_SUCCESS, MB_TYPE_OUT_OF_RANGE, MBENTITYSET, moab::Range::rbegin(), moab::Range::rend(), and moab::OrientedBoxTreeTool::Settings::valid().

Referenced by moab::GeomTopoTool::construct_obb_tree(), and moab::GeomTopoTool::construct_obb_trees().

{
    if( !sets.all_of_type( MBENTITYSET ) ) return MB_TYPE_OUT_OF_RANGE;
    if( settings && !settings->valid() ) return MB_FAILURE;

    // Build initial set data list.
    std::list< SetData > data;
    for( Range::const_iterator i = sets.begin(); i != sets.end(); ++i )
    {
        Range elements;
        ErrorCode rval = instance->get_entities_by_dimension( *i, 2, elements, true );
        if( MB_SUCCESS != rval ) return rval;
        if( elements.empty() ) continue;

        data.push_back( SetData() );
        SetData& set_data = data.back();
        set_data.handle   = *i;
        rval              = OrientedBox::covariance_data_from_tris( set_data.box_data, instance, elements );
        if( MB_SUCCESS != rval ) return rval;
    }

    ErrorCode result = build_sets( data, set_handle_out, 0, settings ? *settings : Settings() );
    if( MB_SUCCESS != result ) return result;

    for( Range::reverse_iterator i = sets.rbegin(); i != sets.rend(); ++i )
    {
        createdTrees.erase( std::remove( createdTrees.begin(), createdTrees.end(), *i ), createdTrees.end() );
    }
    createdTrees.push_back( set_handle_out );
    return MB_SUCCESS;
}

ErrorCode moab::OrientedBoxTreeTool::preorder_traverse	(	EntityHandle	root_set,
		Op &	operation,
		TrvStats *	accum = `0`
	)

Visitor pattern - do operation for each tree node.

Do a preorder traversal of the tree, calling the method in the passed operation instance for each node in the tree. Parent node is visited before either child (pre-order traversal). If operator method passes back the 'descend' argument as false, traversal will not descend to the children of the current node.

Definition at line 498 of file OrientedBoxTreeTool.cpp.

References children, moab::Data::depth, moab::OrientedBoxTreeTool::TrvStats::end_traversal(), ErrorCode, moab::Interface::get_child_meshsets(), moab::OrientedBoxTreeTool::TrvStats::increment(), moab::OrientedBoxTreeTool::TrvStats::increment_leaf(), instance, moab::OrientedBoxTreeTool::Op::leaf(), MB_MULTIPLE_ENTITIES_FOUND, MB_SUCCESS, moab::Data::set, and moab::OrientedBoxTreeTool::Op::visit().

Referenced by obbstat_write(), obbvis_create(), print(), ray_intersect_boxes(), and ray_intersect_sets().

{
    ErrorCode rval;
    std::vector< EntityHandle > children;
    std::vector< Data > the_stack;
    Data data = { set, 0 };
    the_stack.push_back( data );
    int max_depth = -1;

    while( !the_stack.empty() )
    {
        data = the_stack.back();
        the_stack.pop_back();

        // increment traversal statistics
        if( accum )
        {
            accum->increment( data.depth );
            max_depth = std::max( max_depth, data.depth );
        }

        bool descend = true;
        rval         = operation.visit( data.set, data.depth, descend );
        assert( MB_SUCCESS == rval );
        if( MB_SUCCESS != rval ) return rval;

        if( !descend ) continue;

        children.clear();
        rval = instance->get_child_meshsets( data.set, children );
        assert( MB_SUCCESS == rval );
        if( MB_SUCCESS != rval ) return rval;
        if( children.empty() )
        {
            if( accum )
            {
                accum->increment_leaf( data.depth );
            }
            rval = operation.leaf( data.set );
            assert( MB_SUCCESS == rval );
            if( MB_SUCCESS != rval ) return rval;
        }
        else if( children.size() == 2 )
        {
            data.depth++;
            data.set = children[0];
            the_stack.push_back( data );
            data.set = children[1];
            the_stack.push_back( data );
        }
        else
            return MB_MULTIPLE_ENTITIES_FOUND;
    }

    if( accum )
    {
        accum->end_traversal( max_depth );
    }

    return MB_SUCCESS;
}

void moab::OrientedBoxTreeTool::print	(	EntityHandle	tree_root_set,
		std::ostream &	stream,
		bool	list_contents = `false`,
		const char *	id_tag_name = `0`
	)

Print out tree.

Print the tree to an output stream in a human-readable form.

Parameters:

tree_root_set	Entity set representing tree root.
list_contents	If true, list entities in each tree node, If false, just list number of entities.
id_tag_name	If specified, must be the name of an existing integer tag containing an ID for the entities. Not used if list_contents is false.

Definition at line 1653 of file OrientedBoxTreeTool.cpp.

References ErrorCode, instance, MB_SUCCESS, and preorder_traverse().

{
    TreeLayoutPrinter op1( str, instance );
    TreeNodePrinter op2( str, list, true, tag, this );
    ErrorCode r1 = preorder_traverse( set, op1 );
    str << std::endl;
    ErrorCode r2 = preorder_traverse( set, op2 );
    if( r1 != MB_SUCCESS || r2 != MB_SUCCESS )
    {
        std::cerr << "Errors encountered while printing tree\n";
        str << "Errors encountered while printing tree\n";
    }
}

ErrorCode moab::OrientedBoxTreeTool::ray_intersect_boxes	(	Range &	boxes_out,
		EntityHandle	root_set,
		double	tolerance,
		const double	ray_point[3],
		const double	unit_ray_dir[3],
		const double *	ray_length = `0`,
		TrvStats *	accum = `0`
	)

Intersect ray with tree.

Return the tree nodes (as MBENTITYSET handles) for the leaf boxes of the tree intersected by a ray.

Parameters:

boxes_out	The boxes intersected by the ray.
tolerance	The tolerance to use in intersection checks.
ray_point	The base point of the ray.
unit_ray_dir	The ray direction vector (must be unit length)
ray_length	Optional ray length (intersect segment instead of ray.)

Definition at line 815 of file OrientedBoxTreeTool.cpp.

References preorder_traverse().

Referenced by ray_intersect_triangles().

{
    RayIntersector op( this, ray_point, unit_ray_dir, ray_length, tolerance, boxes_out );
    return preorder_traverse( root_set, op, accum );
}

ErrorCode moab::OrientedBoxTreeTool::ray_intersect_sets	(	std::vector< double > &	distances_out,
		std::vector< EntityHandle > &	sets_out,
		std::vector< EntityHandle > &	facets_out,
		EntityHandle	root_set,
		double	tolerance,
		const double	ray_point[3],
		const double	unit_ray_dir[3],
		IntersectSearchWindow &	search_win,
		IntRegCtxt &	register_intersection,
		TrvStats *	accum = `0`
	)

Intersect a ray with the triangles contained within the tree.

Intersect a ray with the triangles contained in the tree and return the distance at which the intersection occured.

Parameters:

distances_out	The output list of intersection points on the ray.
sets_out	The contained set encountered during the tree traversal (see 'set_build'). For the most common use, this is the set corresponding to the geometric surface containing the intersected triangle.
facets_out	Handles of intersected triangles corresponding to distances_out
root_set	The MBENTITYSET representing the root of the tree.
tolerance	The tolerance to use in intersection checks.
ray_point	The base point of the ray.
unit_ray_dir	The ray direction vector (must be unit length)
search_win	An interval that defines the current window in which the an intersection is being sought: (nonnegative, negative)
register_intersection	A context for assessing and registering intersections derived from IntRegCtxt
accum	Optional class for tree traversal statistics.

Definition at line 1036 of file OrientedBoxTreeTool.cpp.

References ErrorCode, moab::OrientedBoxTreeTool::IntRegCtxt::get_facets(), moab::OrientedBoxTreeTool::IntRegCtxt::get_intersections(), moab::OrientedBoxTreeTool::IntRegCtxt::get_sets(), preorder_traverse(), and moab::OrientedBoxTreeTool::TrvStats::ray_tri_tests_count.

Referenced by moab::GeomQueryTool::find_volume(), moab::FBEngine::getPntRayIntsct(), moab::GeomQueryTool::point_in_volume(), moab::GeomQueryTool::ray_fire(), and moab::FBEngine::split_surface_with_direction().

{
    RayIntersectSets op( this, ray_point, unit_ray_dir, tolerance, search_win,
                         accum ? &( accum->ray_tri_tests_count ) : NULL, int_reg_callback );
    ErrorCode rval = preorder_traverse( root_set, op, accum );

    distances_out = int_reg_callback.get_intersections();
    sets_out      = int_reg_callback.get_sets();
    facets_out    = int_reg_callback.get_facets();

    return rval;
}

ErrorCode moab::OrientedBoxTreeTool::ray_intersect_sets	(	std::vector< double > &	distances_out,
		std::vector< EntityHandle > &	sets_out,
		std::vector< EntityHandle > &	facets_out,
		EntityHandle	root_set,
		double	tolerance,
		const double	ray_point[3],
		const double	unit_ray_dir[3],
		const double *	ray_length = `0`,
		TrvStats *	accum = `0`
	)

Definition at line 1059 of file OrientedBoxTreeTool.cpp.

References ErrorCode, moab::OrientedBoxTreeTool::IntRegCtxt::get_facets(), moab::OrientedBoxTreeTool::IntRegCtxt::get_intersections(), moab::OrientedBoxTreeTool::IntRegCtxt::get_sets(), MB_SUCCESS, preorder_traverse(), and moab::OrientedBoxTreeTool::TrvStats::ray_tri_tests_count.

{
    IntRegCtxt int_reg_ctxt;

    OrientedBoxTreeTool::IntersectSearchWindow search_win( ray_length, (double*)0 );

    RayIntersectSets op( this, ray_point, unit_ray_dir, tolerance, search_win,
                         accum ? &( accum->ray_tri_tests_count ) : NULL, int_reg_ctxt );

    ErrorCode rval = preorder_traverse( root_set, op, accum );

    if( MB_SUCCESS != rval )
    {
        return rval;
    }

    distances_out = int_reg_ctxt.get_intersections();
    sets_out      = int_reg_ctxt.get_sets();
    facets_out    = int_reg_ctxt.get_facets();

    return MB_SUCCESS;
}

ErrorCode moab::OrientedBoxTreeTool::ray_intersect_sets	(	EntityHandle	root_set,
		double	tolerance,
		const double	ray_point[3],
		const double	unit_ray_dir[3],
		IntersectSearchWindow &	search_win,
		IntRegCtxt &	register_intersection,
		TrvStats *	accum = `0`
	)

Definition at line 1090 of file OrientedBoxTreeTool.cpp.

References preorder_traverse(), and moab::OrientedBoxTreeTool::TrvStats::ray_tri_tests_count.

{
    RayIntersectSets op( this, ray_point, unit_ray_dir, tolerance, search_win,
                         accum ? &( accum->ray_tri_tests_count ) : NULL, int_reg_callback );
    return preorder_traverse( root_set, op, accum );
}

ErrorCode moab::OrientedBoxTreeTool::ray_intersect_triangles	(	std::vector< double > &	distances_out,
		std::vector< EntityHandle > &	facets_out,
		EntityHandle	root_set,
		double	tolerance,
		const double	ray_point[3],
		const double	unit_ray_dir[3],
		const double *	ray_length = `0`,
		TrvStats *	accum = `0`
	)

Intersect a ray with the triangles contained within the tree.

Intersect a ray with the triangles contained in the tree and return the distance at which the intersection occured.

Parameters:

distances_out	The output list of intersection points on the ray.
facets_out	Handles of intersected triangles corresponding to distances_out
root_set	The MBENTITYSET representing the root of the tree.
tolerance	The tolerance to use in intersection checks.
ray_point	The base point of the ray.
unit_ray_dir	The ray direction vector (must be unit length)
ray_length	Optional ray length (intersect segment instead of ray.)

Definition at line 796 of file OrientedBoxTreeTool.cpp.

References ErrorCode, MB_SUCCESS, ray_intersect_boxes(), and moab::OrientedBoxTreeTool::TrvStats::ray_tri_tests_count.

Referenced by main().

{
    Range boxes;
    ErrorCode rval;

    rval = ray_intersect_boxes( boxes, root_set, tolerance, ray_point, unit_ray_dir, ray_length, accum );
    if( MB_SUCCESS != rval ) return rval;

    return ray_intersect_triangles( intersection_distances_out, intersection_facets_out, boxes, tolerance, ray_point,
                                    unit_ray_dir, ray_length, accum ? &( accum->ray_tri_tests_count ) : NULL );
}

ErrorCode moab::OrientedBoxTreeTool::ray_intersect_triangles	(	std::vector< double > &	intersection_distances_out,
		std::vector< EntityHandle > &	intersection_facets_out,
		const Range &	leaf_boxes_containing_tris,
		double	tolerance,
		const double	ray_point[3],
		const double	unit_ray_dir[3],
		const double *	ray_length = `0`,
		unsigned int *	raytri_test_count = `0`
	)

Intersect ray with triangles contained in passed MBENTITYSETs.

Parameters:

raytri_test_count If non-NULL, count of ray-triangle intersect tests will be added to the value at which this points.

Definition at line 729 of file OrientedBoxTreeTool.cpp.

References moab::Range::begin(), moab::Range::clear(), moab::Range::end(), ErrorCode, moab::Interface::get_connectivity(), moab::Interface::get_coords(), moab::Interface::get_entities_by_handle(), moab::Interface::get_entities_by_type(), instance, MB_SUCCESS, MBTRI, moab::GeomUtil::plucker_ray_tri_intersect(), t, and moab::TYPE_FROM_HANDLE().

{
    ErrorCode rval;
    intersection_distances_out.clear();
#ifdef MB_OBB_USE_VECTOR_QUERIES
    std::vector< EntityHandle > tris;
#endif

    const CartVect point( ray_point );
    const CartVect dir( unit_ray_dir );

    for( Range::iterator b = boxes.begin(); b != boxes.end(); ++b )
    {
#ifndef MB_OBB_USE_VECTOR_QUERIES
        Range tris;
#ifdef MB_OBB_USE_TYPE_QUERIES
        rval = instance->get_entities_by_type( *b, MBTRI, tris );
#else
        rval = instance->get_entities_by_handle( *b, tris );
#endif
#else
        tris.clear();
        rval = instance->get_entities_by_handle( *b, tris );
#endif
        if( MB_SUCCESS != rval ) return rval;
            // dump_fragmentation( tris );

#ifndef MB_OBB_USE_VECTOR_QUERIES
        for( Range::iterator t = tris.begin(); t != tris.end(); ++t )
#else
        for( std::vector< EntityHandle >::iterator t = tris.begin(); t != tris.end(); ++t )
#endif
        {
#ifndef MB_OBB_USE_TYPE_QUERIES
            if( TYPE_FROM_HANDLE( *t ) != MBTRI ) continue;
#endif

            const EntityHandle* conn = NULL;
            int len                  = 0;
            rval                     = instance->get_connectivity( *t, conn, len, true );
            if( MB_SUCCESS != rval ) return rval;

            CartVect coords[3];
            rval = instance->get_coords( conn, 3, coords[0].array() );
            if( MB_SUCCESS != rval ) return rval;

            if( raytri_test_count ) *raytri_test_count += 1;

            double td;
            if( GeomUtil::plucker_ray_tri_intersect( coords, point, dir, td, ray_length ) )
            {
                intersection_distances_out.push_back( td );
                intersection_facets_out.push_back( *t );
            }
        }
    }

    return MB_SUCCESS;
}

ErrorCode moab::OrientedBoxTreeTool::recursive_stats	(	OrientedBoxTreeTool *	tool,
		Interface *	instance,
		EntityHandle	set,
		int	depth,
		StatData &	data,
		unsigned &	count_out,
		CartVect &	dimensions_out
	)		`[private]`

Definition at line 1813 of file OrientedBoxTreeTool.cpp.

References moab::StatData::Ratio::accum(), moab::StatData::Stat< T >::accum(), moab::OrientedBox::area(), moab::StatData::area, box(), children, moab::StatData::count, moab::OrientedBox::dimensions(), moab::StatData::entities, ErrorCode, moab::Interface::get_child_meshsets(), moab::Interface::get_number_entities_by_handle(), moab::OrientedBox::inner_radius(), moab::StatData::leaf_depth, moab::StatData::leaf_ent, MB_MULTIPLE_ENTITIES_FOUND, MB_SUCCESS, moab::measure(), moab::OrientedBox::outer_radius(), moab::StatData::radius, moab::StatData::vol, moab::OrientedBox::volume(), and moab::StatData::volume.

Referenced by stats().

{
    ErrorCode rval;
    OrientedBox tmp_box;
    std::vector< EntityHandle > children( 2 );
    unsigned counts[2];
    bool isleaf;

    ++data.count;

    rval = tool->box( set, tmp_box );
    if( MB_SUCCESS != rval ) return rval;
    children.clear();
    rval = inst->get_child_meshsets( set, children );
    if( MB_SUCCESS != rval ) return rval;
    isleaf = children.empty();
    if( !isleaf && children.size() != 2 ) return MB_MULTIPLE_ENTITIES_FOUND;

    dimensions_out = tmp_box.dimensions();
    data.radius.accum( tmp_box.inner_radius() / tmp_box.outer_radius() );
    data.vol.accum( tmp_box.volume() );
    data.area.accum( tmp_box.area() );

    if( isleaf )
    {
        if( data.leaf_depth.size() <= (unsigned)depth ) data.leaf_depth.resize( depth + 1, 0 );
        ++data.leaf_depth[depth];

        int count;
        rval = inst->get_number_entities_by_handle( set, count );
        if( MB_SUCCESS != rval ) return rval;
        count_out = count;
        data.leaf_ent.accum( count_out );
    }
    else
    {
        for( int i = 0; i < 2; ++i )
        {
            CartVect dims;
            rval = recursive_stats( tool, inst, children[i], depth + 1, data, counts[i], dims );
            if( MB_SUCCESS != rval ) return rval;
            double this_measure, chld_measure;
            int this_dim = measure( dimensions_out, this_measure );
            int chld_dim = measure( dims, chld_measure );
            double ratio;
            if( chld_dim < this_dim )
                ratio = 0;
            else
                ratio = chld_measure / this_measure;

            data.volume.accum( ratio );
        }
        count_out = counts[0] + counts[1];
        data.entities.accum( (double)counts[0] / count_out );
        data.entities.accum( (double)counts[1] / count_out );
    }
    return MB_SUCCESS;
}

ErrorCode moab::OrientedBoxTreeTool::remove_root ( EntityHandle root_set )

Remove obb tree root from the Oriented Box Tree Tool data structure.

Remove obb tree root from the Oriented Box Tree Tool data structure (createdTrees)

Definition at line 480 of file OrientedBoxTreeTool.cpp.

References createdTrees, MB_ENTITY_NOT_FOUND, and MB_SUCCESS.

Referenced by moab::GeomTopoTool::remove_root().

{
    std::vector< EntityHandle >::iterator i = find( createdTrees.begin(), createdTrees.end(), root );
    if( i != createdTrees.end() )
    {
        createdTrees.erase( i );
        return MB_SUCCESS;
    }
    else
        return MB_ENTITY_NOT_FOUND;
}

ErrorCode moab::OrientedBoxTreeTool::sphere_intersect_triangles	(	const double *	center,
		double	radius,
		EntityHandle	tree_root,
		std::vector< EntityHandle > &	facets_out,
		std::vector< EntityHandle > *	sets_out = `0`,
		TrvStats *	accum = `0`
	)

Find facets intersected by a sphere.

Find facets intersected by a spherical volume.

Parameters:

center	Sphere center
radius	Sphere radius
tree_root	Root of OBB tree
facets_out	List of triangles intersecting sphere
sets_out	If not null, sets owning facets are appended to this list in an order corresponding to the entries in facets_out.

Definition at line 570 of file OrientedBoxTreeTool.cpp.

References box(), center(), children, moab::Range::clear(), moab::OrientedBox::closest_location_in_box(), moab::GeomUtil::closest_location_on_tri(), moab::Range::empty(), moab::OrientedBoxTreeTool::TrvStats::end_traversal(), ErrorCode, moab::Range::front(), moab::Interface::get_child_meshsets(), moab::Interface::get_connectivity(), moab::Interface::get_coords(), moab::Interface::get_entities_by_handle(), moab::Interface::get_entities_by_type(), get_moab_instance(), moab::OrientedBoxTreeTool::TrvStats::increment(), moab::OrientedBoxTreeTool::TrvStats::increment_leaf(), instance, MB_SUCCESS, MBENTITYSET, MBTRI, radius, t, and moab::TYPE_FROM_HANDLE().

Referenced by get_close_tris().

{
    const double radsqr = radius * radius;
    OrientedBox b;
    ErrorCode rval;
    Range sets;
    const CartVect center( center_v );
    CartVect closest, coords[3];
    const EntityHandle* conn;
    int num_conn;
#ifndef MB_OBB_USE_VECTOR_QUERIES
    Range tris;
    Range::const_iterator t;
#else
    std::vector< EntityHandle > tris;
    std::vector< EntityHandle >::const_iterator t;
#endif

    std::vector< OBBTreeSITFrame > stack;
    std::vector< EntityHandle > children;
    stack.reserve( 30 );
    stack.push_back( OBBTreeSITFrame( tree_root, 0, 0 ) );

    int max_depth = -1;

    while( !stack.empty() )
    {
        EntityHandle surf = stack.back().surf;
        EntityHandle node = stack.back().node;
        int current_depth = stack.back().depth;
        stack.pop_back();

        // increment traversal statistics.
        if( accum )
        {
            accum->increment( current_depth );
            max_depth = std::max( max_depth, current_depth );
        }

        if( !surf && sets_out )
        {
            rval = get_moab_instance()->get_entities_by_type( node, MBENTITYSET, sets );
            if( !sets.empty() ) surf = sets.front();
            sets.clear();
        }

        // check if sphere intersects box
        rval = box( node, b );
        if( MB_SUCCESS != rval ) return rval;
        b.closest_location_in_box( center, closest );
        closest -= center;
        if( ( closest % closest ) > radsqr ) continue;

        // push child boxes on stack
        children.clear();
        rval = instance->get_child_meshsets( node, children );
        if( MB_SUCCESS != rval ) return rval;
        if( !children.empty() )
        {
            assert( children.size() == 2 );
            stack.push_back( OBBTreeSITFrame( children[0], surf, current_depth + 1 ) );
            stack.push_back( OBBTreeSITFrame( children[1], surf, current_depth + 1 ) );
            continue;
        }

        if( accum )
        {
            accum->increment_leaf( current_depth );
        }
        // if leaf, intersect sphere with triangles
#ifndef MB_OBB_USE_VECTOR_QUERIES
#ifdef MB_OBB_USE_TYPE_QUERIES
        rval = get_moab_instance()->get_entities_by_type( node, MBTRI, tris );
#else
        rval = get_moab_instance()->get_entities_by_handle( node, tris );
#endif
        t = tris.begin();
#else
        rval = get_moab_instance()->get_entities_by_handle( node, tris );
        t    = tris.lower_bound( MBTRI );
#endif
        if( MB_SUCCESS != rval ) return rval;

        for( t = tris.begin(); t != tris.end(); ++t )
        {
#ifndef MB_OBB_USE_VECTOR_QUERIES
            if( TYPE_FROM_HANDLE( *t ) != MBTRI ) break;
#elif !defined( MB_OBB_USE_TYPE_QUERIES )
            if( TYPE_FROM_HANDLE( *t ) != MBTRI ) continue;
#endif
            rval = get_moab_instance()->get_connectivity( *t, conn, num_conn, true );
            if( MB_SUCCESS != rval ) return rval;
            if( num_conn != 3 ) continue;

            rval = get_moab_instance()->get_coords( conn, num_conn, coords[0].array() );
            if( MB_SUCCESS != rval ) return rval;

            GeomUtil::closest_location_on_tri( center, coords, closest );
            closest -= center;
            if( ( closest % closest ) <= radsqr &&
                std::find( facets_out.begin(), facets_out.end(), *t ) == facets_out.end() )
            {
                facets_out.push_back( *t );
                if( sets_out ) sets_out->push_back( surf );
            }
        }
    }

    if( accum )
    {
        accum->end_traversal( max_depth );
    }

    return MB_SUCCESS;
}

ErrorCode moab::OrientedBoxTreeTool::stats	(	EntityHandle	tree_root_set,
		std::ostream &	stream
	)

Print tree statistics.

Print misc. stats. describing tree

Definition at line 1927 of file OrientedBoxTreeTool.cpp.

References moab::StatData::area, moab::StatData::count, moab::StatData::entities, ErrorCode, moab::StatData::Ratio::hist, instance, moab::StatData::leaf_depth, moab::StatData::leaf_ent, moab::StatData::Ratio::max, moab::StatData::Stat< T >::max, MB_SUCCESS, moab::StatData::Ratio::min, moab::StatData::Stat< T >::min, moab::StatData::radius, recursive_stats(), moab::StatData::Ratio::sqr, moab::StatData::Stat< T >::sqr, moab::std_dev(), moab::StatData::Ratio::sum, moab::StatData::Stat< T >::sum, moab::StatData::vol, moab::StatData::volume, WE, and WW.

Referenced by obbstat_write(), and moab::SmoothFace::SmoothFace().

{
    StatData d;
    ErrorCode rval;
    unsigned total_entities, i;
    CartVect total_dim;

    rval = recursive_stats( this, instance, set, 0, d, total_entities, total_dim );
    if( MB_SUCCESS != rval ) return rval;

    unsigned tree_height    = d.leaf_depth.size();
    unsigned min_leaf_depth = tree_height, num_leaves = 0;
    unsigned max_leaf_per_depth = 0;
    double sum_leaf_depth = 0, sqr_leaf_depth = 0;
    for( i = 0; i < d.leaf_depth.size(); ++i )
    {
        unsigned val = d.leaf_depth[i];
        num_leaves += val;
        sum_leaf_depth += (double)val * i;
        sqr_leaf_depth += (double)val * i * i;
        if( val && i < min_leaf_depth ) min_leaf_depth = i;
        if( max_leaf_per_depth < val ) max_leaf_per_depth = val;
    }
    unsigned num_non_leaf = d.count - num_leaves;

    double rv = total_dim[0] * total_dim[1] * total_dim[2];
    s << "entities in tree:  " << total_entities << std::endl
      << "root volume:       " << rv << std::endl
      << "total node volume: " << d.vol.sum << std::endl
      << "total/root volume: " << d.vol.sum / rv << std::endl
      << "tree height:       " << tree_height << std::endl
      << "node count:        " << d.count << std::endl
      << "leaf count:        " << num_leaves << std::endl
      << std::endl;

    double avg_leaf_depth = sum_leaf_depth / num_leaves;
    double rms_leaf_depth = sqrt( sqr_leaf_depth / num_leaves );
    double std_leaf_depth = std_dev( sqr_leaf_depth, sum_leaf_depth, num_leaves );

    double avg_leaf_ent = d.leaf_ent.sum / num_leaves;
    double rms_leaf_ent = sqrt( d.leaf_ent.sqr / num_leaves );
    double std_leaf_ent = std_dev( d.leaf_ent.sqr, d.leaf_ent.sum, num_leaves );

    unsigned num_child = 2 * num_non_leaf;

    double avg_vol_ratio = d.volume.sum / num_child;
    double rms_vol_ratio = sqrt( d.volume.sqr / num_child );
    double std_vol_ratio = std_dev( d.volume.sqr, d.volume.sum, num_child );

    double avg_ent_ratio = d.entities.sum / num_child;
    double rms_ent_ratio = sqrt( d.entities.sqr / num_child );
    double std_ent_ratio = std_dev( d.entities.sqr, d.entities.sum, num_child );

    double avg_rad_ratio = d.radius.sum / d.count;
    double rms_rad_ratio = sqrt( d.radius.sqr / d.count );
    double std_rad_ratio = std_dev( d.radius.sqr, d.radius.sum, d.count );

    double avg_vol = d.vol.sum / d.count;
    double rms_vol = sqrt( d.vol.sqr / d.count );
    double std_vol = std_dev( d.vol.sqr, d.vol.sum, d.count );

    double avg_area = d.area.sum / d.count;
    double rms_area = sqrt( d.area.sqr / d.count );
    double std_area = std_dev( d.area.sqr, d.area.sum, d.count );

    int prec = s.precision();
    s << "                   " WW "Minimum" WW "Average" WW "RMS" WW "Maximum" WW "Std.Dev." << std::endl;
    s << std::setprecision( 1 )
      << "Leaf Depth         " WW min_leaf_depth WW avg_leaf_depth WW rms_leaf_depth WW d.leaf_depth.size() -
             1 WW std_leaf_depth
      << std::endl;
    s << std::setprecision( 0 )
      << "Entities/Leaf      " WW d.leaf_ent.min WW avg_leaf_ent WW rms_leaf_ent WW d.leaf_ent.max WW std_leaf_ent
      << std::endl;
    s << std::setprecision( 3 )
      << "Child Volume Ratio " WW d.volume.min WW avg_vol_ratio WW rms_vol_ratio WW d.volume.max WW std_vol_ratio
      << std::endl;
    s << std::setprecision( 3 )
      << "Child Entity Ratio " WW d.entities.min WW avg_ent_ratio WW rms_ent_ratio WW d.entities.max WW std_ent_ratio
      << std::endl;
    s << std::setprecision( 3 )
      << "Box Radius Ratio   " WW d.radius.min WW avg_rad_ratio WW rms_rad_ratio WW d.radius.max WW std_rad_ratio
      << std::endl;
    s << std::setprecision( 0 ) << "Box volume         " WE d.vol.min WE avg_vol WE rms_vol WE d.vol.max WE std_vol
      << std::endl;
    s << std::setprecision( 0 ) << "Largest side area  " WE d.area.min WE avg_area WE rms_area WE d.area.max WE std_area
      << std::endl;
    s << std::setprecision( prec ) << std::endl;

    s << "Leaf Depth Histogram (Root depth is 0)" << std::endl;
    double f = 60.0 / max_leaf_per_depth;
    for( i = min_leaf_depth; i < d.leaf_depth.size(); ++i )
        s << std::setw( 2 ) << i << " " << std::setw( 5 ) << d.leaf_depth[i] << " |" << std::setfill( '*' )
          << std::setw( (int)floor( f * d.leaf_depth[i] + 0.5 ) ) << "" << std::setfill( ' ' ) << std::endl;
    s << std::endl;

    s << "Child/Parent Volume Ratio Histogram" << std::endl;
    f = 60.0 / *( std::max_element( d.volume.hist, d.volume.hist + 10 ) );
    for( i = 0; i < 10u; ++i )
        s << "0." << i << " " << std::setw( 5 ) << d.volume.hist[i] << " |" << std::setfill( '*' )
          << std::setw( (int)floor( f * d.volume.hist[i] + 0.5 ) ) << "" << std::setfill( ' ' ) << std::endl;
    s << std::endl;

    s << "Child/Parent Entity Count Ratio Histogram" << std::endl;
    f = 60.0 / *( std::max_element( d.entities.hist, d.entities.hist + 10 ) );
    for( i = 0; i < 10u; ++i )
        s << "0." << i << " " << std::setw( 5 ) << d.entities.hist[i] << " |" << std::setfill( '*' )
          << std::setw( (int)floor( f * d.entities.hist[i] + 0.5 ) ) << "" << std::setfill( ' ' ) << std::endl;
    s << std::endl;

    s << "Inner/Outer Radius Ratio Histogram (~0.70 for cube)" << std::endl;
    // max radius ratio for a box is about 0.7071.  Move any boxes
    // in the .7 bucket into .6 and print .0 to .6.
    d.radius.hist[6] += d.radius.hist[7];
    f = 60.0 / *( std::max_element( d.entities.hist, d.entities.hist + 7 ) );
    for( i = 0; i < 7u; ++i )
        s << "0." << i << " " << std::setw( 5 ) << d.entities.hist[i] << " |" << std::setfill( '*' )
          << std::setw( (int)floor( f * d.entities.hist[i] + 0.5 ) ) << "" << std::setfill( ' ' ) << std::endl;
    s << std::endl;

    return MB_SUCCESS;
}

ErrorCode moab::OrientedBoxTreeTool::stats	(	EntityHandle	set,
		unsigned &	entities_in_tree,
		double &	root_volume,
		double &	tot_node_volume,
		double &	tot_to_root_volume,
		unsigned &	tree_height,
		unsigned &	node_count,
		unsigned &	num_leaves
	)

Get tree statistics.

Get summary stats. describing tree

Parameters:

set	Root of tree for which data is requested
total_entities	Entities in tree
root_volume	Total volume of root box
tot_node_volume	Total volume in all nodes
tot_to_root_volume	Ratio of total / root volume
tree_height	Maximum height of tree, from root to leaf
node_count	Number of nodes in tree
num_leaves	Number of leaf nodes in tree

Definition at line 1888 of file OrientedBoxTreeTool.cpp.

References moab::StatData::count, ErrorCode, instance, moab::StatData::leaf_depth, MB_SUCCESS, recursive_stats(), moab::StatData::Stat< T >::sum, and moab::StatData::vol.

{
    StatData d;
    ErrorCode rval;
    unsigned i;
    CartVect total_dim;

    rval = recursive_stats( this, instance, set, 0, d, total_entities, total_dim );
    if( MB_SUCCESS != rval ) return rval;

    tree_height                 = d.leaf_depth.size();
    unsigned min_leaf_depth     = tree_height;
    num_leaves                  = 0;
    unsigned max_leaf_per_depth = 0;
    // double sum_leaf_depth = 0, sqr_leaf_depth = 0;
    for( i = 0; i < d.leaf_depth.size(); ++i )
    {
        unsigned val = d.leaf_depth[i];
        num_leaves += val;
        // sum_leaf_depth += (double)val*i;
        // sqr_leaf_depth += (double)val*i*i;
        if( val && i < min_leaf_depth ) min_leaf_depth = i;
        if( max_leaf_per_depth < val ) max_leaf_per_depth = val;
    }
    rv                 = total_dim[0] * total_dim[1] * total_dim[2];
    tot_node_volume    = d.vol.sum;
    tot_to_root_volume = d.vol.sum / rv;
    node_count         = d.count;

    return MB_SUCCESS;
}