CGAL users discussion list

[Mariette Yvinec <>]

When the line goes through a vertex a the convex hull
you can get
1 when the line is within the cone formed by the triangles
of the convex hull incident to the vertex
and 0 otherwise.

For the direction  and the cube you choose,
its seems that you can get 0 for any vertex on
(or close to)
one of the edges
of the cube.




Le 01/08/11 14:03, 

 a écrit :
> Hi,
> I mesh a cube. Then I consider a line crossing each mesh vertex and I look
> for the intersections between the line and the triangle of each incident
> cell opposite to the vertex.  At the end of each iteration on the incident
> cells I expect one intersection (with a finite cell) if one of the
> intersected cells is on the convex hull (vertex on boundary), two
> intersection if both the intersected cells are finite (vertex in the
> volume). Otherwise I found also the case of intersection number=0. For the
> case of finite cells incident to the vertex I guess this is due to the
> convexity of the small neighborhood of the vertex and it should occur only
> for the direction tangent to the meshed object . This is not the case in
> the minimal example I attach. (When both finite and  infinite incident
> cells are considered this
> should corresponds to have two intersections with infinite cells.) I
> verify that intersection number=0 occurs only for vertexes incident to the
> infinite vertex.
> Am I experiencing an accuracy problem on the boundary? The result doesn't
> change refining the facet parameters.
> Thanks in advance for any suggestion
>
>
> #include<CGAL/Exact_predicates_inexact_constructions_kernel.h>
> #include<CGAL/Exact_predicates_exact_constructions_kernel.h>
> #include<CGAL/Cartesian_converter.h>
> #include<CGAL/Mesh_triangulation_3.h>
> #include<CGAL/Mesh_complex_3_in_triangulation_3.h>
> #include<CGAL/Mesh_criteria_3.h>
> #include<CGAL/Implicit_mesh_domain_3.h>
> #include<CGAL/make_mesh_3.h>
> #include<CGAL/intersections.h>
>
> // Domain
> struct K: public CGAL::Exact_predicates_inexact_constructions_kernel {};
> typedef K::FT FT;
> typedef K::Point_3 Point;
> typedef FT (Function)(const Point&);
> typedef CGAL::Implicit_mesh_domain_3<Function,K>  Mesh_domain;
>
> // Triangulation/domain
> typedef CGAL::Mesh_triangulation_3<Mesh_domain>::type Tr;
> typedef CGAL::Mesh_complex_3_in_triangulation_3<Tr>  C3t3;
>
> // Criteria
> typedef CGAL::Mesh_criteria_3<Tr>  Mesh_criteria;
>
> // Triangulation
> typedef  C3t3::Triangulation Tr;
> typedef  C3t3::Cell_handle Cell_handle;
> typedef  Tr::Finite_vertices_iterator Finite_vertices_iterator;
> typedef  Tr::Vertex_handle Vertex_handle;
>
> //new kernel:exact constructions
> typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel;
> typedef CGAL::Cartesian_converter<K,Kernel>  converter_exact;
> typedef CGAL::Cartesian_converter<Kernel,K>  converter_inexact;
>
> //Geometric objects for both the kernels
> typedef K::Segment_3 Segment;
> typedef K::Triangle_3 Triangle;
> typedef Kernel::Direction_3 Direction;
> typedef Kernel::Line_3 Lineconv;
> typedef Kernel::Point_3 Pointex;
> typedef Kernel::Triangle_3 Triangleconv;
> typedef Kernel::Segment_3 Segmentconv;
>
> // To avoid verbose function and named parameters call
> using namespace CGAL::parameters;
>
> // Function
> FT cube_function (const Point&  p){
>         if( (p.x())<  1&&  -(p.x())<  1&&
>        (p.y())<  1&&  -(p.y())<  1&&
>         (p.z())<  1&&  -(p.z())<  1 )
>                 return -1.;
>         return 1.;
> }
>
> int main()
> {
> // Domain
> Mesh_domain domain(cube_function, K::Sphere_3(CGAL::ORIGIN, 2.));
>
> // Criteria
> Mesh_criteria criteria(facet_angle=30, facet_size=0.3, facet_distance=0.1,
>            cell_radius_edge=2, cell_size=0.3);
>
>         // Mesh generation and optimization
>         C3t3 c3t3 = CGAL::make_mesh_3<C3t3>(domain, criteria,
>                                                                         
> lloyd(time_limit=30),
>                                                                               
>   no_perturb(),
>                                                                               
>   exude(time_limit=30, sliver_bound=10));
>
>         //converter to_exact to evaluate intersection points;
>         converter_exact to_exact;
>         converter_inexact to_inexact;
>         Point picos(0.2,0.5,0.1);
>         Point porig(0.0,0.0,0.0);
>         Pointex picos_ex=to_exact(picos);
>         Pointex porig_ex=to_exact(porig);
>         Segmentconv seg(porig_ex,picos_ex);
>         Direction dir(seg.direction());
>
>         int count_vert=0;
>         for(Finite_vertices_iterator vit =
> c3t3.triangulation().finite_vertices_begin();
>                 vit != c3t3.triangulation().finite_vertices_end();
>                 ++vit)
>         {
>                 std::cout<<"SCAN VERTEX"<<++count_vert<<"th 
> vertex"<<std::endl;
>                 Point vertexpoint = vit->point();
>                 Pointex vertexpoint_ex=to_exact(vertexpoint);
>                 Lineconv line_vertex(vertexpoint_ex,dir);
>                 std::vector<Cell_handle>  incident_cells;
>                 bool is_on_chull=false;
>                 int nb_inter=0;
>                 
> c3t3.triangulation().incident_cells(vit,std::back_inserter(incident_cells));
>                         for (std::vector<Cell_handle>::iterator 
> cit=incident_cells.begin();
>                                 cit!=incident_cells.end();
>                                 ++cit)
>                         {
>                                 if (c3t3.triangulation().is_infinite(*cit)){
>                                         is_on_chull=true;
>                                         continue;
>                                 }
>                                 int index=(*cit)->index(vit);
>                                 K::Triangle_3 
> triangle((*cit)->vertex((index+1)%4)->point(),
>                                                                            
> (*cit)->vertex((index+2)%4)->point(),
>                                                                            
> (*cit)->vertex((index+3)%4)->point());
>                                 Kernel::Triangle_3 tr2=to_exact(triangle);
>                                 if( CGAL::do_intersect(line_vertex,tr2) ) 
> ++nb_inter;
>                         }
>                     //or CGAL_assertion
>                         if( ((is_on_chull&&  nb_inter==1) || (!is_on_chull&&
> nb_inter==2))==true ){
>                         std::cout<<"OK"<<std::endl;
>
>                         }
>                         else {
>                                 std::cout<<"ERRORINTERSEC"<<std::endl;
>                                 }
>                         }
>
> return 0;
> }

--
Mariette Yvinec
Geometrica project team
INRIA  Sophia-Antipolis