CGAL users discussion list

["Sebastien Loriot (GeometryFactory)" <>]

I think the fastest solution is to do your intersection tests on the set 
of segments that is incident to a facet in the complex. That way you
don't compute too much.

It might need to be benchmarked but if you pre-compute in parallel the 
position each vertex incident to that set of segments with reference
to the plane, then plane segment intersection predicate could be written
based on that precomputed orientation and no duplicated computation will
be done.

Best,

Sebastien.


On 4/27/23 22:58, Maria Eduarda Veras Martins ( via 
cgal-discuss Mailing List) wrote:
> Hi!
> I need to intersect several planes with a volumetric mesh. I'm making an 
> aabb tree with the faces of the tetrahedron but I don't know how to 
> efficiently intersect the planes.
> Here is my code so far (only one plane):
>
> #include<list>
> #include<CGAL/AABB_tree.h>
> #include<CGAL/AABB_traits.h>
> #include<CGAL/AABB_triangle_primitive.h>
> #include<CGAL/Exact_predicates_inexact_constructions_kernel.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/iterator.h>
> #include<iostream>
> #include<fstream>
> #include<map>
> #include<vector>
> #include<string>
> #include<CGAL/utility.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
> 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;
> typedef Mesh_criteria::Facet_criteria Facet_criteria;
> typedef Mesh_criteria::Cell_criteria Cell_criteria;
>
> typedef C3t3::Triangulation Tr;
> typedef C3t3::Cell_handle Cell_handle;
> typedef Tr::Facet Facet;
> typedef Tr::Finite_facets_iterator Finite_facets_iterator;
> typedef Tr::Finite_vertices_iterator Finite_vertices_iterator;
> typedef Tr::Finite_cells_iterator Finite_cells_iterator;
> typedef Tr::Vertex_handle Vertex_handle;
>
>
> typedef K::Plane_3 Plane;
> typedef K::Line_3 Line;
> typedef K::Direction_3 Direction;
> typedef K::Triangle_3 Triangle;
> typedef K::Segment_3 Segment;
> typedef std::list<Triangle>::iterator Iterator;
> typedef CGAL::AABB_triangle_primitive<K,Iterator> Primitive;
> typedef CGAL::AABB_traits<K, Primitive> AABB_triangle_traits;
> typedef CGAL::AABB_tree<AABB_triangle_traits> Tree;
> typedef Tree::Object_and_primitive_id Object_and_primitive_id;
> typedef Tree::Primitive_id Primitive_id;
>
> // Function
> FT sphere_function (const Point& p)
> {
> const FT x2=p.x()*p.x();
> const FT y2=p.y()*p.y();
> const FT z2=p.z()*p.z();
>
> return x2+y2+z2-1;
> }
>
> int main() {
> // Domain
> Mesh_domain domain(sphere_function, K::Sphere_3(CGAL::ORIGIN, 1.));
>
> // Criteria
> Facet_criteria facet_criteria(30, 0.1, 0.01);// angle, size, approximation
> Cell_criteria cell_criteria(3, 0.1);// radius-edge ratio, size
> Mesh_criteria criteria(facet_criteria, cell_criteria);
>
> // Mesh generation
> C3t3 c3t3 = CGAL::make_mesh_3<C3t3>(domain, criteria);
>
> // Output
> std::ofstream medit_file("outmatrix.mesh");
> c3t3.output_to_medit(medit_file);
>
> //std::cout<< "Number of verticesDT"<<" "<< 
> c3t3.triangulation().number_of_vertices() << std::endl;
> //std::cout<< "Number of facetsDT"<<" "<< 
> c3t3.triangulation().number_of_finite_facets() << std::endl;
> //std::cout<< "Number of facetsC3T3"<<" "<<c3t3.number_of_facets() << 
> std::endl;
> //std::cout<< "Number of cellsC3T3"<<" "<<c3t3.number_of_cells() << 
> std::endl;
> //std::cout<< "Number of cellsDT"<<" 
> "<<c3t3.triangulation().number_of_finite_cells() << std::endl;
>
>
> // for each facet build a triangle
> std::list<Triangle> triangles;
> for(Finite_facets_iterator fi = 
> c3t3.triangulation().finite_facets_begin(); fi != 
> c3t3.triangulation().finite_facets_end(); fi++) {
>
> Cell_handle cell = fi->first;
> int opposite_vertex_index = fi->second;
> if(c3t3.is_in_complex(cell)) {
> Triangle t=c3t3.triangulation().triangle(cell,opposite_vertex_index);
> triangles.push_back(t);
> }
> else if (c3t3.is_in_complex(*fi)) {
> Triangle t=c3t3.triangulation().triangle(*fi);
> triangles.push_back(t);
> }
>
> }
>
> Tree tree(triangles.begin(),triangles.end());
>
> // set plane and count intersections
> Point p(0, 0, 0);
> Plane plane(p, K::Vector_3(1, 0, 0));
>
> std::list <Object_and_primitive_id> intersections;
> tree.all_intersections(plane, std::back_inserter(intersections));
>
> //how to do it with many planes?
> //I'm thinking of looping the plans, but I wanted something more efficient
>
>
> return 0;
> }
>
> --
> You are currently subscribed to cgal-discuss.
> To unsubscribe or access the archives, go to
> https://sympa.inria.fr/sympa/info/cgal-discuss