CGAL users discussion list

Hi,
I'm sorry to ask on the same topic, but I could not interpret my results.
So just a question on my starting hypothesis.
I recall I combined a 3D mesh generation procedure with an AABB tree. My
query object
is a line with fixed direction and crossing a vertex of the mesh (I
consider parallel rays crossing each vertex for different intersection
queries) and I use a kind of AABB_primitive, such as  Id is
Tr::Facet_finite_iterator.
With the condition
f(vit->point()==cell->vertex(opposite_vertex_index)->point()) in the
attached code, for each vertex in the iteration of
Finite_vertices_iterator, I want select the intersected facets,in the list
of intersected primitives, which have this vertex as opposite vertex.
I would expect to find just two cells incident to the vertex and
containing the intersected facets and sort along the direction of the
intersection line query the two intersection points lying in these two
facets opposite to the vertex.
This should allow me to define an intersection point "before" and an
intersection point "after" the vertex of interest for each step.

Conversely I detect almost everywhere just one intersection and not 2 as
expected.Two intersection are very seldom detected. It depends on the mesh
size. If I reduce the size of tetrahedra and facets I detect more double
intersections. Otherwise, when these sizes are increased double
intersections disappear.

1) So, fixed a vertex and a line crossing it there are not necessarily two
incident cells containing two intersected facets which are facets opposite
to the fixed vertex?I'm wondering whether this hypothesis is false, or I
make some mistakes in applying it.

Below the code.....

Thanks for your attention!

int main()
{....MESH A SPHERE
Tree
tree(c3t3.triangulation().finite_facets_begin(),c3t3.triangulation().finite_facets_end());
Point a(0.0, 0.0, 0.0);
Point b(0.0, 1.0, 0.0);
Ray ray_query(a,b);
Direction dir(ray_query);
std::map<Vertex_handle, int> V;
int inum = 0;
for( Finite_vertices_iterator vit =
c3t3.triangulation().finite_vertices_begin();
       vit != c3t3.triangulation().finite_vertices_end();
       ++vit)
{
       V[vit] = inum++;
}

        for( Finite_vertices_iterator vit =
c3t3.triangulation().finite_vertices_begin();
                vit != c3t3.triangulation().finite_vertices_end();
                ++vit)
        {
                Point c = vit->point();
                Line line_vertex(c,dir);
                std::cout << 
tree.number_of_intersected_primitives(line_vertex)
                << " intersected primitives"<<std::endl;
                // computes all intersections with line query
                std::list<Object_and_primitive_id> intersections;
                tree.all_intersections(line_vertex, 
std::back_inserter(intersections));
std::list<Object_and_primitive_id>::iterator it_intersections=
intersections.begin();
                std::list<Primitive_id> primitives;
                std::vector<std::pair<Point,Id> > myvector;
                std::vector<std::pair<Point,Id> >::iterator it;
       tree.all_intersected_primitives(line_vertex,
std::back_inserter(primitives));
                for(; it_intersections != intersections.end() ; 
++it_intersections)

                {
                        // get intersection object
                        Object_and_primitive_id op = *it_intersections;
                        CGAL::Object object = op.first;
                        Id id=op.second;


                        Point point;
                        Segment segment;
                        if(CGAL::assign(point,object))
                        {
                                it=myvector.begin();
                                myvector.insert(it, make_pair(point,id));
                        }
                        else if(CGAL::assign(segment,object))
                                std::cout << "intersection object is a 
segment" <<std::endl;
                }
                std::sort(myvector.begin(), myvector.end(), sort_pred());

for ( it=myvector.begin() ; it != myvector.end(); it++ ){
                        C3t3::Cell_handle cell    = (*it).second->first;
                        int opposite_vertex_index = (*it).second->second;
const Vertex_handle& vopp=cell->vertex(opposite_vertex_index);
if(vit->point()==cell->vertex(opposite_vertex_index)->point()) {
std::cout<<vit->point()<< " "<<" "<<" "<<"
"<<cell->vertex(opposite_vertex_index)->point()<<std::endl;
if(c3t3.is_in_complex(cell)==true) {
                                std::cout << " intersection/coordinates"<<"
"<<std::setprecision(15)<<(*it).first << " => " << std::endl;
for(int i = 0; i < 4; i++){
const Vertex_handle& vh = (*((*it).second)).first->vertex(i);
                                                                std::cout<<" 
vertex index:"<<" " <<V[vh] << " "<<std::endl;
                                        std::cout <<" vertex coordinates:"<<" 
"<<vh->point().x()<<"
"<<vh->point().y()<<" "<<vh->point().z()<<std::endl;
}
}
}
}