An implementation of an Intersection Registration Context for use GQT ray-firing. More...

Inheritance diagram for moab::GQT_IntRegCtxt:

Collaboration diagram for moab::GQT_IntRegCtxt:

Public Member Functions
	GQT_IntRegCtxt (OrientedBoxTreeTool obbtool, const double ray_point[3], const double ray_dir[3], double tolerance, int min_tolerance_intersections, const EntityHandle root_set, const EntityHandle geom_volume, const Tag sense_tag, const int desired_orient, const std::vector< EntityHandle > prev_facets)
virtual ErrorCode	register_intersection (EntityHandle set, EntityHandle triangle, double distance, OrientedBoxTreeTool::IntersectSearchWindow &, GeomUtil::intersection_type int_type)
virtual ErrorCode	update_orient (EntityHandle set, int *surfTriOrient)
virtual const int *	getDesiredOrient ()
Private Member Functions
void	add_intersection (EntityHandle set, EntityHandle tri, double dist, OrientedBoxTreeTool::IntersectSearchWindow &search_win)
void	append_intersection (EntityHandle set, EntityHandle facet, double dist)
void	set_intersection (int len_idx, EntityHandle set, EntityHandle facet, double dist)
void	add_mode1_intersection (EntityHandle set, EntityHandle facet, double dist, OrientedBoxTreeTool::IntersectSearchWindow &search_win)
bool	edge_node_piercing_intersect (const EntityHandle tri, const CartVect &ray_direction, const GeomUtil::intersection_type int_type, const std::vector< EntityHandle > &close_tris, const std::vector< int > &close_senses, const Interface MBI, std::vector< EntityHandle > neighborhood_tris=0)
	Determine if a ray-edge/node intersection is glancing or piercing. This function avoids asking for upward adjacencies to prevent their creation.
bool	in_prevFacets (const EntityHandle tri)
bool	in_neighborhoods (const EntityHandle tri)
Private Attributes
OrientedBoxTreeTool *	tool
const CartVect	ray_origin
const CartVect	ray_direction
const double	tol
const int	minTolInt
const EntityHandle *	rootSet
const EntityHandle *	geomVol
const Tag *	senseTag
const int *	desiredOrient
const std::vector< EntityHandle > *	prevFacets
std::vector< std::vector < EntityHandle > >	neighborhoods
std::vector< EntityHandle >	neighborhood

Detailed Description

An implementation of an Intersection Registration Context for use GQT ray-firing.

This context uses a variety of tests and conditions to confirm whether or not to accumulate an intersection, to ensure robustness for ray firing.

This context only accumulates intersections that are oriented parallel to the 'desiredOrient', if provided, with respect to 'geomVol', using information in the in 'senseTag'.

This context only accumulates a single intersection out of a set of multiple intersections that fall in the same 'neighborhood', where a 'neighborhood' is defined as facets that share edges or vertices.

This context only accumulates piercing intersections. This is relevant for intersections that are found to be on an edge or vertex by the Plucker test. Such intersections are piercing if the ray has the same orientation w.r.t. to all fecets that share that edge or vertex.

This context tests intersections against a list of 'prevFacets' to prevent a ray from crossing the same facet more than once. The user is responsible for ensuring that this list is reset when appropriate.

This context accumulates all intersections within 'tol' of the start of the ray and if the number of intersections within the 'tol' of the ray start point is less than 'minTolInt', the next closest intersection. If the desired result is only the closest intersection, 'minTolInt' should be 0. This function will return all intersections, regardless of distance from the start of the ray, if 'minTolInt' is negative.

Definition at line 123 of file GeomQueryTool.cpp.

Constructor & Destructor Documentation

moab::GQT_IntRegCtxt::GQT_IntRegCtxt	(	OrientedBoxTreeTool *	obbtool,
		const double	ray_point[3],
		const double	ray_dir[3],
		double	tolerance,
		int	min_tolerance_intersections,
		const EntityHandle *	root_set,
		const EntityHandle *	geom_volume,
		const Tag *	sense_tag,
		const int *	desired_orient,
		const std::vector< EntityHandle > *	prev_facets
	)		`[inline]`

Definition at line 176 of file GeomQueryTool.cpp.

        : tool( obbtool ), ray_origin( ray_point ), ray_direction( ray_dir ), tol( tolerance ),
          minTolInt( min_tolerance_intersections ), rootSet( root_set ), geomVol( geom_volume ), senseTag( sense_tag ),
          desiredOrient( desired_orient ), prevFacets( prev_facets ){

                                           };

Member Function Documentation

void moab::GQT_IntRegCtxt::add_intersection	(	EntityHandle	set,
		EntityHandle	tri,
		double	dist,
		OrientedBoxTreeTool::IntersectSearchWindow &	search_win
	)		`[private]`

Definition at line 532 of file GeomQueryTool.cpp.

References add_mode1_intersection(), append_intersection(), moab::OrientedBoxTreeTool::IntRegCtxt::intersections, minTolInt, neighborhood, neighborhoods, set_intersection(), and tol.

Referenced by register_intersection().

{

    // Mode 1, detected by non-null neg_ray_len pointer
    // keep the closest nonneg intersection and one negative intersection, if closer
    if( search_win.second && search_win.first )
    {
        return add_mode1_intersection( set, facet, dist, search_win );
    }

    // ---------------------------------------------------------------------------
    /*   Mode 2: Used if neg_ray_len is not specified
         variables used:     min_tol_int, tol, search_win.first
         variables not used: neg_ray_len
         1) if(min_tol_int<0) return all intersections
         2) otherwise return all inside tolerance and unless there are >min_tol_int
         inside of tolerance, return the closest outside of tolerance */
    // Mode 2
    // If minTolInt is less than zero, return all intersections
    if( minTolInt < 0 && dist > -tol )
    {
        append_intersection( set, facet, dist );
        neighborhoods.push_back( neighborhood );
        return;
    }

    // Check if the 'len' pointer is pointing into the intersection
    // list.  If this is the case, then the list contains, at that
    // location, an intersection greater than the tolerance away from
    // the base point of the ray.
    int len_idx = -1;
    if( search_win.first && search_win.first >= &intersections[0] &&
        search_win.first < &intersections[0] + intersections.size() )
        len_idx = search_win.first - &intersections[0];

    // If the intersection is within tol of the ray base point, we
    // always add it to the list.
    if( dist <= tol )
    {
        // If the list contains an intersection outside the tolerance...
        if( len_idx >= 0 )
        {
            // If we no longer want an intersection outside the tolerance,
            // remove it.
            if( (int)intersections.size() >= minTolInt )
            {
                set_intersection( len_idx, set, facet, dist );
                // From now on, we want only intersections within the tolerance,
                // so update length accordingly
                search_win.first = &tol;
            }
            // Otherwise appended to the list and update pointer
            else
            {
                append_intersection( set, facet, dist );
                search_win.first = &intersections[len_idx];
            }
        }
        // Otherwise just append it
        else
        {
            append_intersection( set, facet, dist );
            // If we have all the intersections we want, set
            // length such that we will only find further intersections
            // within the tolerance
            if( (int)intersections.size() >= minTolInt ) search_win.first = &tol;
        }
    }
    // Otherwise the intersection is outside the tolerance
    // If we already have an intersection outside the tolerance and
    // this one is closer, replace the existing one with this one.
    else if( len_idx >= 0 )
    {
        if( dist <= *search_win.first )
        {
            set_intersection( len_idx, set, facet, dist );
        }
    }
    // Otherwise if we want an intersection outside the tolerance
    // and don'thave one yet, add it.
    else if( (int)intersections.size() < minTolInt )
    {
        append_intersection( set, facet, dist );
        // update length.  this is currently the closest intersection, so
        // only want further intersections that are closer than this one.
        search_win.first = &intersections.back();
    }
}

void moab::GQT_IntRegCtxt::add_mode1_intersection	(	EntityHandle	set,
		EntityHandle	facet,
		double	dist,
		OrientedBoxTreeTool::IntersectSearchWindow &	search_win
	)		`[private]`

Definition at line 494 of file GeomQueryTool.cpp.

References moab::OrientedBoxTreeTool::IntRegCtxt::facets, moab::OrientedBoxTreeTool::IntRegCtxt::intersections, set_intersection(), and moab::OrientedBoxTreeTool::IntRegCtxt::sets.

Referenced by add_intersection().

{
    if( 2 != intersections.size() )
    {
        intersections.resize( 2, 0 );
        sets.resize( 2, 0 );
        facets.resize( 2, 0 );
        // must initialize this for comparison below
        intersections[0] = -std::numeric_limits< double >::max();
    }

    // negative case
    if( 0.0 > dist )
    {
        set_intersection( 0, set, facet, dist );
        search_win.second = &intersections[0];
        // nonnegative case
    }
    else
    {
        set_intersection( 1, set, facet, dist );
        search_win.first = &intersections[1];
        // if the intersection is closer than the negative one, remove the negative one
        if( dist < -*( search_win.second ) )
        {
            set_intersection( 0, 0, 0, -intersections[1] );
            search_win.second = &intersections[0];
        }
    }
    //    std::cout << "add_intersection: dist = " << dist << " search_win.second=" <<
    //    *search_win.second
    //          << " search_win.first=" << *search_win.first << std::endl;
    return;
}

void moab::GQT_IntRegCtxt::append_intersection	(	EntityHandle	set,
		EntityHandle	facet,
		double	dist
	)		`[private]`

Definition at line 472 of file GeomQueryTool.cpp.

References moab::OrientedBoxTreeTool::IntRegCtxt::facets, moab::OrientedBoxTreeTool::IntRegCtxt::intersections, neighborhood, neighborhoods, and moab::OrientedBoxTreeTool::IntRegCtxt::sets.

Referenced by add_intersection().

{
    intersections.push_back( dist );
    sets.push_back( set );
    facets.push_back( facet );
    neighborhoods.push_back( neighborhood );
    return;
}

bool moab::GQT_IntRegCtxt::edge_node_piercing_intersect	(	const EntityHandle	tri,
		const CartVect &	ray_dir,
		const GeomUtil::intersection_type	int_type,
		const std::vector< EntityHandle > &	close_tris,
		const std::vector< int > &	close_senses,
		const Interface *	MBI,
		std::vector< EntityHandle > *	neighborhood_tris = `0`
	)		`[private]`

Determine if a ray-edge/node intersection is glancing or piercing. This function avoids asking for upward adjacencies to prevent their creation.

Parameters:

tri	The intersected triangle
ray_dir	The direction of the ray
int_type	The type of intersection (EDGE0, EDGE1, NODE2, ...)
close_tris	Vector of triangles in the proximity of the intersection
close_senses	Vector of surface senses for tris in the proximity of the intersection
neighborhood	Vector of triangles in the topological neighborhood of the intersection

Returns:: True if piercing, false otherwise.

Definition at line 278 of file GeomQueryTool.cpp.

References moab::GeomUtil::EDGE0, moab::GeomUtil::EDGE1, moab::GeomUtil::EDGE2, ErrorCode, moab::Interface::get_connectivity(), moab::Interface::get_coords(), moab::Interface::list_entities(), MB_CHK_ERR_RET_VAL, moab::GeomUtil::NODE0, moab::GeomUtil::NODE1, and moab::GeomUtil::NODE2.

Referenced by register_intersection().

{

    // get the node of the triangle
    const EntityHandle* conn = NULL;
    int len                  = 0;
    ErrorCode rval           = MBI->get_connectivity( tri, conn, len );MB_CHK_ERR_RET_VAL( rval, false );
    // NOTE: original code is next line, but return type was wrong; returning true to get same net effect
    if( 3 != len ) return false;

    // get adjacent tris (and keep their corresponding senses)
    std::vector< EntityHandle > adj_tris;
    std::vector< int > adj_senses;

    // node intersection
    if( GeomUtil::NODE0 == int_type || GeomUtil::NODE1 == int_type || GeomUtil::NODE2 == int_type )
    {

        // get the intersected node
        EntityHandle node;
        if( GeomUtil::NODE0 == int_type )
            node = conn[0];
        else if( GeomUtil::NODE1 == int_type )
            node = conn[1];
        else
            node = conn[2];

        // get tris adjacent to node
        for( unsigned i = 0; i < close_tris.size(); ++i )
        {
            const EntityHandle* con = NULL;
            rval                    = MBI->get_connectivity( close_tris[i], con, len );MB_CHK_ERR_RET_VAL( rval, false );
            if( 3 != len ) return false;

            if( node == con[0] || node == con[1] || node == con[2] )
            {
                adj_tris.push_back( close_tris[i] );
                adj_senses.push_back( close_senses[i] );
            }
        }
        if( adj_tris.empty() )
        {
            std::cerr << "error: no tris are adjacent to the node" << std::endl;
            return false;
        }
        // edge intersection
    }
    else if( GeomUtil::EDGE0 == int_type || GeomUtil::EDGE1 == int_type || GeomUtil::EDGE2 == int_type )
    {

        // get the endpoints of the edge
        EntityHandle endpts[2];
        if( GeomUtil::EDGE0 == int_type )
        {
            endpts[0] = conn[0];
            endpts[1] = conn[1];
        }
        else if( GeomUtil::EDGE1 == int_type )
        {
            endpts[0] = conn[1];
            endpts[1] = conn[2];
        }
        else
        {
            endpts[0] = conn[2];
            endpts[1] = conn[0];
        }

        // get tris adjacent to edge
        for( unsigned i = 0; i < close_tris.size(); ++i )
        {
            const EntityHandle* con = NULL;
            rval                    = MBI->get_connectivity( close_tris[i], con, len );MB_CHK_ERR_RET_VAL( rval, false );
            if( 3 != len ) return false;

            // check both orientations in case close_tris are not on the same surface
            if( ( endpts[0] == con[0] && endpts[1] == con[1] ) || ( endpts[0] == con[1] && endpts[1] == con[0] ) ||
                ( endpts[0] == con[1] && endpts[1] == con[2] ) || ( endpts[0] == con[2] && endpts[1] == con[1] ) ||
                ( endpts[0] == con[2] && endpts[1] == con[0] ) || ( endpts[0] == con[0] && endpts[1] == con[2] ) )
            {
                adj_tris.push_back( close_tris[i] );
                adj_senses.push_back( close_senses[i] );
            }
        }
        // In a 2-manifold each edge is adjacent to exactly 2 tris
        if( 2 != adj_tris.size() )
        {
            std::cerr << "error: edge of a manifold must be topologically adjacent to exactly 2 tris" << std::endl;
            MBI->list_entities( endpts, 2 );
            return true;
        }
    }
    else
    {
        std::cerr << "error: special case not an node/edge intersection" << std::endl;
        return false;
    }

    // The close tris were in proximity to the intersection. The adj_tris are
    // topologically adjacent to the intersection (the neighborhood).
    if( neighborhood_tris ) ( *neighborhood_tris ).assign( adj_tris.begin(), adj_tris.end() );

    // determine glancing/piercing
    // If a desired_orientation was used in this call to ray_intersect_sets,
    // the plucker_ray_tri_intersect will have already used it. For a piercing
    // intersection, the normal of all tris must have the same orientation.
    int sign = 0;
    for( unsigned i = 0; i < adj_tris.size(); ++i )
    {
        const EntityHandle* con = NULL;
        rval                    = MBI->get_connectivity( adj_tris[i], con, len );MB_CHK_ERR_RET_VAL( rval, false );
        if( 3 != len ) return false;

        CartVect coords[3];
        rval = MBI->get_coords( con, len, coords[0].array() );MB_CHK_ERR_RET_VAL( rval, false );

        // get normal of triangle
        CartVect v0     = coords[1] - coords[0];
        CartVect v1     = coords[2] - coords[0];
        CartVect norm   = adj_senses[i] * ( v0 * v1 );
        double dot_prod = norm % ray_dir;

        // if the sign has not yet been decided, choose it
        if( 0 == sign && 0 != dot_prod )
        {
            if( 0 < dot_prod )
                sign = 1;
            else
                sign = -1;
        }

        // intersection is glancing if tri and ray do not point in same direction
        // for every triangle
        if( 0 != sign && 0 > sign * dot_prod ) return false;
    }
    return true;
}

virtual const int* moab::GQT_IntRegCtxt::getDesiredOrient ( ) [inline, virtual]

Reimplemented from moab::OrientedBoxTreeTool::IntRegCtxt.

Definition at line 199 of file GeomQueryTool.cpp.

References desiredOrient.

    {
        return desiredOrient;
    };

bool moab::GQT_IntRegCtxt::in_neighborhoods ( const EntityHandle tri ) [private]

Definition at line 252 of file GeomQueryTool.cpp.

References neighborhoods.

Referenced by register_intersection().

{
    bool same_neighborhood = false;
    for( unsigned i = 0; i < neighborhoods.size(); ++i )
    {
        if( neighborhoods[i].end() != find( neighborhoods[i].begin(), neighborhoods[i].end(), tri ) )
        {
            same_neighborhood = true;
            continue;
        }
    }
    return same_neighborhood;
}

bool moab::GQT_IntRegCtxt::in_prevFacets ( const EntityHandle tri ) [private]

Definition at line 247 of file GeomQueryTool.cpp.

References prevFacets.

Referenced by register_intersection().

{
    return ( prevFacets && ( ( *prevFacets ).end() != find( ( *prevFacets ).begin(), ( *prevFacets ).end(), tri ) ) );
}

ErrorCode moab::GQT_IntRegCtxt::register_intersection	(	EntityHandle	set,
		EntityHandle	triangle,
		double	distance,
		OrientedBoxTreeTool::IntersectSearchWindow &	search_win,
		GeomUtil::intersection_type	int_type
	)		`[virtual]`

Reimplemented from moab::OrientedBoxTreeTool::IntRegCtxt.

Definition at line 422 of file GeomQueryTool.cpp.

References add_intersection(), edge_node_piercing_intersect(), ErrorCode, geomVol, moab::OrientedBoxTreeTool::get_close_tris(), moab::OrientedBoxTreeTool::get_moab_instance(), in_neighborhoods(), in_prevFacets(), moab::GeomUtil::INTERIOR, MB_SUCCESS, neighborhood, ray_direction, ray_origin, rootSet, senseTag, tol, and tool.

{
    ErrorCode rval;

    // Do not accept intersections if they are in the vector of previously intersected
    // facets.
    if( in_prevFacets( t ) ) return MB_SUCCESS;

    // Do not accept intersections if they are in the neighborhood of previous
    // intersections.
    if( in_neighborhoods( t ) ) return MB_SUCCESS;

    neighborhood.clear();

    // Handle special case of edge/node intersection. Accept piercing
    // intersections and reject glancing intersections.
    // The edge_node_intersection function needs to know surface sense wrt volume.
    // A less-robust implementation could work without sense information.
    // Would it ever be useful to accept a glancing intersection?
    if( GeomUtil::INTERIOR != int_type && rootSet && geomVol && senseTag )
    {
        // get triangles in the proximity of the intersection
        std::vector< EntityHandle > close_tris;
        std::vector< int > close_senses;
        rval = tool->get_close_tris( ray_origin + int_dist * ray_direction, tol, rootSet, geomVol, senseTag, close_tris,
                                     close_senses );

        if( MB_SUCCESS != rval ) return rval;

        if( !edge_node_piercing_intersect( t, ray_direction, int_type, close_tris, close_senses,
                                           tool->get_moab_instance(), &neighborhood ) )
            return MB_SUCCESS;
    }
    else
    {
        neighborhood.push_back( t );
    }

    // NOTE: add_intersection may modify the 'neg_ray_len' and 'nonneg_ray_len'
    //       members, which will affect subsequent calls to ray_tri_intersect
    //       in this loop.
    add_intersection( set, t, int_dist, search_win );

    return MB_SUCCESS;
}

void moab::GQT_IntRegCtxt::set_intersection	(	int	len_idx,
		EntityHandle	set,
		EntityHandle	facet,
		double	dist
	)		`[private]`

Definition at line 481 of file GeomQueryTool.cpp.

References moab::OrientedBoxTreeTool::IntRegCtxt::facets, moab::OrientedBoxTreeTool::IntRegCtxt::intersections, and moab::OrientedBoxTreeTool::IntRegCtxt::sets.

Referenced by add_intersection(), and add_mode1_intersection().

{
    intersections[len_idx] = dist;
    sets[len_idx]          = set;
    facets[len_idx]        = facet;
    return;
}

ErrorCode moab::GQT_IntRegCtxt::update_orient	(	EntityHandle	set,
		int *	surfTriOrient
	)		`[virtual]`

Reimplemented from moab::OrientedBoxTreeTool::IntRegCtxt.

Definition at line 205 of file GeomQueryTool.cpp.

References desiredOrient, ErrorCode, geomVol, moab::OrientedBoxTreeTool::get_moab_instance(), MB_SUCCESS, senseTag, moab::Interface::tag_get_data(), and tool.

{

    ErrorCode rval;

    // Get desired orientation of surface wrt volume. Use this to return only
    // exit or entrance intersections.
    if( geomVol && senseTag && desiredOrient && surfTriOrient )
    {
        if( 1 != *desiredOrient && -1 != *desiredOrient )
        {
            std::cerr << "error: desired orientation must be 1 (forward) or -1 (reverse)" << std::endl;
        }
        EntityHandle vols[2];
        rval = tool->get_moab_instance()->tag_get_data( *senseTag, &set, 1, vols );
        assert( MB_SUCCESS == rval );
        if( MB_SUCCESS != rval ) return rval;
        if( vols[0] == vols[1] )
        {
            std::cerr << "error: surface has positive and negative sense wrt same volume" << std::endl;
            return MB_FAILURE;
        }
        // surfTriOrient will be used by plucker_ray_tri_intersect to avoid
        // intersections with wrong orientation.
        if( *geomVol == vols[0] )
        {
            *surfTriOrient = *desiredOrient * 1;
        }
        else if( *geomVol == vols[1] )
        {
            *surfTriOrient = *desiredOrient * ( -1 );
        }
        else
        {
            assert( false );
            return MB_FAILURE;
        }
    }

    return MB_SUCCESS;
}