CGAL users discussion list

Thanks for your reply? Might I ask some questions about your comments?

<Here you create an indexation of finite vertices of the triangulation, in
the
<order of iteration using finite_vertices_(begin|end)()... The numeration is
<only for you. Internally the triangulation data structure does not keep such
<an indexation (just handles or iterators).

 1) So how I can retrieve and store the "true" indices? I created my
numeration following the procedure in File_medit.h to write medit files
but I probably misunderstood. I need to store the true indices in order
to solve equations in the corresponding points of my object (domain).

<Here opposite_vertex_index is an integer i, where:
 <0 <= i < 4
<that corresponds to the number of the facet in the cell. It means that the
<vertex that is opposite to the facet (*it) is cell-
<vertex(opposite_vertex_index). Note that this is very different to
<V[opposite_vertex_index]!!

2) So in this way I have not a loop on the vertexes of the cell containing
the intersected facet opposite to the vertex (crossed by the line query)?
The indices of this opposite facet and of the vertex crossed are the
indices I need.

<Here you do not test the equality between opposite_vertex_index and V[vit]
<(anyway, I just said they can only be equal by chance). Instead, you do
<something similar to:
 <opposite_vertex_index = V[vit];
 <if(opposite_vertex_index != 0) { ....}

3) Is this just a syntax error (opposite_vertex_index = V[vit] instead of
opposite_vertex_index == V[vit])? In this case I find less vertexes
fulfilling my condition but not yet just two intersection points. I think
the problem is again the numeration I created.

<"setting"? Did you mean "testing"? c3t3.is_in_complex(cell) returns true if
<the cell is in the volume (3D complex). The facets of a cell that is in the
<volume can either be internal to the volume or be on the surface (2D
complex).

4) So, in order to select only the cells belonging the complex and not all
the cells in the triangulation I have to use the condition
c3t3.is_in_complex(cell) ==true?

<I do not know. I had a look at your code but I do not know how you fill the
<vector myvector. Maybe the error is here. Actually I found other errors (see
<below)...

I attach the procedure to fill the vector myvector.

        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;
if(opposite_vertex_index==V[vit]) {
if(c3t3.is_in_complex(cell)) {
                                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;
}
}
}
}

Thanks again!