From 9a8f26ced7cf757e85fad70fbb1b40ca9313fc41 Mon Sep 17 00:00:00 2001 From: Maxime Gimeno Date: Tue, 22 Oct 2019 14:00:46 +0200 Subject: [PATCH] More changes after review --- .../CGAL/Intersections_3/Bbox_3_Bbox_3.h | 15 +- .../Iso_cuboid_3_Plane_3_intersection.h | 215 +++++++++--------- .../Tetrahedron_3_Line_3_intersection.h | 9 +- .../Tetrahedron_3_Plane_3_intersection.h | 2 +- .../Intersections_3/test_intersections_3.cpp | 1 - 5 files changed, 119 insertions(+), 123 deletions(-) diff --git a/Intersections_3/include/CGAL/Intersections_3/Bbox_3_Bbox_3.h b/Intersections_3/include/CGAL/Intersections_3/Bbox_3_Bbox_3.h index fae33a52e50..8ff63ab2883 100644 --- a/Intersections_3/include/CGAL/Intersections_3/Bbox_3_Bbox_3.h +++ b/Intersections_3/include/CGAL/Intersections_3/Bbox_3_Bbox_3.h @@ -47,15 +47,12 @@ intersection(const CGAL::Bbox_3& a, if(!do_intersect(a,b)) return Result_type(); - double xmin, xmax, ymin, ymax, zmin, zmax; - xmin = (std::max)(a.xmin(), b.xmin()); - xmax = (std::min)(a.xmax(), b.xmax()); - - ymin = (std::max)(a.ymin(), b.ymin()); - ymax = (std::min)(a.ymax(), b.ymax()); - - zmin = (std::max)(a.zmin(), b.zmin()); - zmax = (std::min)(a.zmax(), b.zmax()); + double xmin = (std::max)(a.xmin(), b.xmin()); + double xmax = (std::min)(a.xmax(), b.xmax()); + double ymin = (std::max)(a.ymin(), b.ymin()); + double ymax = (std::min)(a.ymax(), b.ymax()); + double zmin = (std::max)(a.zmin(), b.zmin()); + double zmax = (std::min)(a.zmax(), b.zmax()); return Result_type(std::forward(Bbox_3(xmin, ymin, zmin, xmax, ymax, zmax))); } diff --git a/Intersections_3/include/CGAL/Intersections_3/internal/Iso_cuboid_3_Plane_3_intersection.h b/Intersections_3/include/CGAL/Intersections_3/internal/Iso_cuboid_3_Plane_3_intersection.h index c684578035f..d4d11fd50cd 100644 --- a/Intersections_3/include/CGAL/Intersections_3/internal/Iso_cuboid_3_Plane_3_intersection.h +++ b/Intersections_3/include/CGAL/Intersections_3/internal/Iso_cuboid_3_Plane_3_intersection.h @@ -33,20 +33,17 @@ #include namespace CGAL { - namespace Intersections { - namespace internal { template -void filter_points(const std::vector& input, +void filter_points(std::vector& input, std::vector& output) { - std::set tmp; - for( auto p : input) - tmp.insert(p); - for(auto p: tmp) - output.push_back(p); + std::sort(input.begin(), input.end()); + auto last = std::unique(input.begin(), input.end()); + input.erase(last, input.end()); + output = std::move(input); } //Tetrahedron_3 Line_3 @@ -55,7 +52,7 @@ typename Intersection_traits:: intersection( const typename K::Iso_cuboid_3 &cub, const typename K::Plane_3 &pl, - const K&) + const K& k) { typedef typename K::Point_3 Point_3; typedef typename K::Segment_3 Segment_3; @@ -74,11 +71,11 @@ intersection( //get all edges of cub for(int i=0; i< 4; ++i) - edges.push_back(Segment_3(cub.vertex(i), cub.vertex((i+1)%4))); + edges.emplace_back(cub.vertex(i), cub.vertex((i+1)%4)); for(int i=0; i < 4; ++i) - edges.push_back(Segment_3(cub.vertex(i+4), cub.vertex((i+1)%4+4))); + edges.emplace_back(cub.vertex(i+4), cub.vertex((i+1)%4+4)); for(int i=0; i < 4; ++i) - edges.push_back(Segment_3(cub.vertex(i), cub.vertex((i+1)%4+4))); + edges.emplace_back(cub.vertex(i), cub.vertex((i+1)%4+4)); //get all intersections between pl and cub edges std::vector segments; @@ -98,8 +95,14 @@ intersection( } } + if(segments.empty() && points.empty()) + return result_type(); + switch(segments.size()) { + case 0: + //points dealt with later + break; case 1: { //adj to an edge @@ -107,48 +110,48 @@ intersection( { return result_type(std::forward(segments.front())); } - //through an edge not on diagonal + //plane intersecting through an edge (not 2) else { Poly res(4); - Segment_3 front(segments.front()); + const Segment_3& entry_seg(segments.front()); Point_3 p1, p2; - bool has_p1(false), - has_p2(false); + bool p1_found(false), + p2_found(false); - for(auto p : points) + for(const Point_3& p : points) { - if(!front.has_on(p)) + if(!k.equal_3_object()(entry_seg.source(), p) + && ! k.equal_3_object()(entry_seg.target(), p)) { - if(!has_p1) + if(!p1_found) { p1 = p; - has_p1 = true; + p1_found = true; } else { p2 = p; - has_p2 = true; + p2_found = true; break; } } } - CGAL_assertion(has_p1 && has_p2); - Segment_3 back(p1,p2); - res[0] = front.target(); + CGAL_assertion(p1_found && p2_found); + res[0] = entry_seg.target(); - if((front.target() - front.source()) - * (back.target() - back.source()) > 0) + if((entry_seg.target() - p1) + * (p2 - p1) > 0) { - res[1] = back.target(); - res[2] = back.source(); + res[1] = p2; + res[2] = p1; } else { - res[1] = back.source(); - res[2] = back.target(); + res[1] = p1; + res[2] = p2; } - res[3] = front.source(); + res[3] = entry_seg.source(); return result_type(std::forward(res)); } @@ -198,96 +201,92 @@ intersection( Poly filtered_points; filter_points(points, filtered_points); - if(filtered_points.empty()) - return result_type(); //adjacent to a vertex if(filtered_points.size() == 1) { return result_type(std::forward(filtered_points.front())); } - else + + //get intersections between pl and each face -> line. Foreach line, creates segment with points. Then use helper_function to recover polygon. + typedef typename Intersection_traits, + CGAL::Plane_3 >::result_type Pl_pl_type; + + std::vector plane_intersections; + Pl_pl_type pl_inter = CGAL::intersection(pl, Plane_3(cub.vertex(0), + cub.vertex(1), + cub.vertex(5))); + if(const Line_3* line = boost::get(&*pl_inter)){ + plane_intersections.push_back(*line); + } + pl_inter = CGAL::intersection(pl, Plane_3(cub.vertex(0), + cub.vertex(3), + cub.vertex(4))); + if(const Line_3* line = boost::get(&*pl_inter)){ + plane_intersections.push_back(*line); + } + pl_inter = CGAL::intersection(pl, Plane_3(cub.vertex(0), + cub.vertex(1), + cub.vertex(3))); + if(const Line_3* line = boost::get(&*pl_inter)){ + plane_intersections.push_back(*line); + } + pl_inter = CGAL::intersection(pl, Plane_3(cub.vertex(7), + cub.vertex(6), + cub.vertex(1))); + if(const Line_3* line = boost::get(&*pl_inter)){ + plane_intersections.push_back(*line); + } + pl_inter = CGAL::intersection(pl, Plane_3(cub.vertex(7), + cub.vertex(4), + cub.vertex(3))); + if(const Line_3* line = boost::get(&*pl_inter)){ + plane_intersections.push_back(*line); + } + pl_inter = CGAL::intersection(pl, Plane_3(cub.vertex(7), + cub.vertex(6), + cub.vertex(4))); + if(const Line_3* line = boost::get(&*pl_inter)){ + plane_intersections.push_back(*line); + } + + std::list tmp_segs; + for(auto line : plane_intersections) { - //get intersections between pl and each face -> line. Foreach line, creates segment with points. Then use helper_function to recover polygon. - typedef typename Intersection_traits, - CGAL::Plane_3 >::result_type Pl_pl_type; - - std::vector plane_intersections; - Pl_pl_type pl_inter = CGAL::intersection(pl, Plane_3(cub.vertex(0), - cub.vertex(1), - cub.vertex(5))); - if(const Line_3* line = boost::get(&*pl_inter)){ - plane_intersections.push_back(*line); - } - pl_inter = CGAL::intersection(pl, Plane_3(cub.vertex(0), - cub.vertex(3), - cub.vertex(4))); - if(const Line_3* line = boost::get(&*pl_inter)){ - plane_intersections.push_back(*line); - } - pl_inter = CGAL::intersection(pl, Plane_3(cub.vertex(0), - cub.vertex(1), - cub.vertex(3))); - if(const Line_3* line = boost::get(&*pl_inter)){ - plane_intersections.push_back(*line); - } - pl_inter = CGAL::intersection(pl, Plane_3(cub.vertex(7), - cub.vertex(6), - cub.vertex(1))); - if(const Line_3* line = boost::get(&*pl_inter)){ - plane_intersections.push_back(*line); - } - pl_inter = CGAL::intersection(pl, Plane_3(cub.vertex(7), - cub.vertex(4), - cub.vertex(3))); - if(const Line_3* line = boost::get(&*pl_inter)){ - plane_intersections.push_back(*line); - } - pl_inter = CGAL::intersection(pl, Plane_3(cub.vertex(7), - cub.vertex(6), - cub.vertex(4))); - if(const Line_3* line = boost::get(&*pl_inter)){ - plane_intersections.push_back(*line); - } - - std::list tmp_segs; - for(auto line : plane_intersections) + bool first_found = false; + Point_3 first_p; + for(auto p : filtered_points) { - bool first_found = false; - Point_3 first_p; - for(auto p : filtered_points) + if(line.has_on(p)) { - if(line.has_on(p)) + if(!first_found) { - if(!first_found) - { - first_found = true; - first_p = p; - } - else - { - tmp_segs.push_back(Segment_3(first_p, p)); - break; - } + first_found = true; + first_p = p; + } + else + { + tmp_segs.push_back(Segment_3(first_p, p)); + break; } } } - if(tmp_segs.size() < 3) - return result_type(); - std::list tmp_pts; - fill_points_list(tmp_segs,tmp_pts); - Poly res; - for(auto p : tmp_pts) - res.push_back(p); - if(res.size() == 3){ - typename K::Triangle_3 tr(res[0], res[1], res[2]); - return result_type(std::forward(tr)); - } - else - { - return result_type(std::forward(res)); - } + } + if(tmp_segs.size() < 3) + return result_type(); + std::list tmp_pts; + fill_points_list(tmp_segs,tmp_pts); + Poly res; + for(auto p : tmp_pts) + res.push_back(p); + if(res.size() == 3){ + typename K::Triangle_3 tr(res[0], res[1], res[2]); + return result_type(std::forward(tr)); + } + else + { + return result_type(std::forward(res)); } } diff --git a/Intersections_3/include/CGAL/Intersections_3/internal/Tetrahedron_3_Line_3_intersection.h b/Intersections_3/include/CGAL/Intersections_3/internal/Tetrahedron_3_Line_3_intersection.h index 51fc6c4a46e..f9a8e8f7895 100644 --- a/Intersections_3/include/CGAL/Intersections_3/internal/Tetrahedron_3_Line_3_intersection.h +++ b/Intersections_3/include/CGAL/Intersections_3/internal/Tetrahedron_3_Line_3_intersection.h @@ -37,12 +37,13 @@ template struct Tetrahedron_Line_intersection_3 :public Tetrahedron_lines_intersection_3_base > { - typedef typename K::Line_3 O; typedef Tetrahedron_lines_intersection_3_base > Base; typedef typename Base::Result_type Result_type; - Tetrahedron_Line_intersection_3(const typename K::Tetrahedron_3& tet, - const O& o):Base(tet,o) {} + Tetrahedron_Line_intersection_3( + const typename K::Tetrahedron_3& tet, + const typename K::Line_3& o):Base(tet,o) + {} bool all_inside_test() { @@ -64,7 +65,7 @@ intersection( const typename K::Line_3 &line, const K&) { - Tetrahedron_lines_intersection_3_base > solver(tet, line); + Tetrahedron_Line_intersection_3 solver(tet, line); solver.do_procede(); return solver.output; } diff --git a/Intersections_3/include/CGAL/Intersections_3/internal/Tetrahedron_3_Plane_3_intersection.h b/Intersections_3/include/CGAL/Intersections_3/internal/Tetrahedron_3_Plane_3_intersection.h index 8dcb4cd085d..9206ddcd4e5 100644 --- a/Intersections_3/include/CGAL/Intersections_3/internal/Tetrahedron_3_Plane_3_intersection.h +++ b/Intersections_3/include/CGAL/Intersections_3/internal/Tetrahedron_3_Plane_3_intersection.h @@ -36,7 +36,7 @@ namespace Intersections { namespace internal { -//Tetrahedron_3 Segment_3 +//Tetrahedron_3 Plane_3 template typename Intersection_traits::result_type intersection( diff --git a/Intersections_3/test/Intersections_3/test_intersections_3.cpp b/Intersections_3/test/Intersections_3/test_intersections_3.cpp index 0eca13bfb97..7292f661203 100644 --- a/Intersections_3/test/Intersections_3/test_intersections_3.cpp +++ b/Intersections_3/test/Intersections_3/test_intersections_3.cpp @@ -1215,6 +1215,5 @@ int main() Test< CGAL::Homogeneous >().run(); Test< CGAL::Epeck >().run(true); Test< CGAL::Homogeneous >().run(true); - // TODO : test more kernels. }