diff --git a/Operations_on_polyhedra/include/CGAL/internal/Polyhedron_plane_clipping_3.h b/Operations_on_polyhedra/include/CGAL/internal/Polyhedron_plane_clipping_3.h deleted file mode 100644 index 04698711d42..00000000000 --- a/Operations_on_polyhedra/include/CGAL/internal/Polyhedron_plane_clipping_3.h +++ /dev/null @@ -1,437 +0,0 @@ -// Copyright (c) 2015 GeometryFactory (France). -// All rights reserved. -// -// This file is part of CGAL (www.cgal.org). -// You can redistribute it and/or modify it under the terms of the GNU -// General Public License as published by the Free Software Foundation, -// either version 3 of the License, or (at your option) any later version. -// -// Licensees holding a valid commercial license may use this file in -// accordance with the commercial license agreement provided with the software. -// -// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE -// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. -// -// $URL$ -// $Id$ -// -// -// Author(s) : Sebastien Loriot - -#ifndef CGAL_INTERNAL_POLYHEDRON_PLANE_CLIPPING_3_H -#define CGAL_INTERNAL_POLYHEDRON_PLANE_CLIPPING_3_H - -#include -#include -#include -#include -#include -#include -#include -#include -#include - -namespace CGAL{ -namespace corefinement{ - -namespace internal{ - template - class Builder_from_T_2 : public CGAL::Modifier_base { - typedef std::map Vertex_map; - - const T& t; - template - static unsigned get_vertex_index( Vertex_map& vertex_map, - typename T::Vertex_handle vh, - Builder& builder, - unsigned& vindex) - { - std::pair - res=vertex_map.insert(std::make_pair(vh,vindex)); - if (res.second){ - builder.add_vertex(vh->point()); - ++vindex; - } - return res.first->second; - } - - public: - Builder_from_T_2(const T& t_):t(t_) - { - CGAL_assertion(t.dimension()==2); - } - void operator()( HDS& hds) { - // Postcondition: `hds' is a valid polyhedral surface. - CGAL::Polyhedron_incremental_builder_3 B( hds, true); - Vertex_map vertex_map; - //start the surface - B.begin_surface( t.number_of_vertices(), t.number_of_faces()); - unsigned vindex=0; - for (typename T::Finite_faces_iterator it=t.finite_faces_begin(); - it!=t.finite_faces_end();++it) - { - unsigned i0=get_vertex_index(vertex_map,it->vertex(0),B,vindex); - unsigned i1=get_vertex_index(vertex_map,it->vertex(1),B,vindex); - unsigned i2=get_vertex_index(vertex_map,it->vertex(2),B,vindex); - B.begin_facet(); - B.add_vertex_to_facet( i0 ); - B.add_vertex_to_facet( i1 ); - B.add_vertex_to_facet( i2 ); - B.end_facet(); - } - B.end_surface(); - } - }; -} // end of namespace internal - -template -Polyhedron clip_to_bbox(const Bbox_3& bbox, const Plane_3& plane) -{ - typedef typename Polyhedron::Traits::Kernel Kernel; - typedef typename Kernel::Point_3 Point_3; - typedef typename Kernel::Segment_3 Segment_3; - cpp11::array corners= {{ - Point_3(bbox.xmin(),bbox.ymin(),bbox.zmin()), - Point_3(bbox.xmin(),bbox.ymax(),bbox.zmin()), - Point_3(bbox.xmax(),bbox.ymax(),bbox.zmin()), - Point_3(bbox.xmax(),bbox.ymin(),bbox.zmin()), - Point_3(bbox.xmin(),bbox.ymin(),bbox.zmax()), - Point_3(bbox.xmin(),bbox.ymax(),bbox.zmax()), - Point_3(bbox.xmax(),bbox.ymax(),bbox.zmax()), - Point_3(bbox.xmax(),bbox.ymin(),bbox.zmax()) - }}; - - cpp11::array orientations = {{ - plane.oriented_side(corners[0]), - plane.oriented_side(corners[1]), - plane.oriented_side(corners[2]), - plane.oriented_side(corners[3]), - plane.oriented_side(corners[4]), - plane.oriented_side(corners[5]), - plane.oriented_side(corners[6]), - plane.oriented_side(corners[7]) - }}; - - std::vector intersection_points; - // first look for intersections at corners - for (int i=0; i<8; ++i) - if (orientations[i]==ON_ORIENTED_BOUNDARY) - intersection_points.push_back(corners[i]); - // second look for intersections on edges - cpp11::array edge_indices = {{ // 2 *12 edges - 0,1, 1,2, 2,3, 3,0, // bottom face edges - 4,5, 5,6, 6,7, 7,4, // top face edges - 0,4, 1,5, 2,6, 3,7 - }}; - - for (int i=0; i<12; ++i) - { - int i1=edge_indices[2*i], i2=edge_indices[2*i+1]; - if (orientations[i1]==ON_ORIENTED_BOUNDARY) continue; - if (orientations[i2]==ON_ORIENTED_BOUNDARY) continue; - if (orientations[i1]!=orientations[i2]) - intersection_points.push_back( - boost::get( - *CGAL::intersection(plane, Segment_3(corners[i1], corners[i2]) ) - ) - ); - } - - Polyhedron P; - - //if less that 3 points there will be nothing to clipping. - if (intersection_points.size()<3) return P; - - //triangulate the set of intersection points (I know it's overkill) - typedef CGAL::Triangulation_2_projection_traits_3 P_traits; - typedef CGAL::Delaunay_triangulation_2 DT; - DT dt(P_traits(plane.orthogonal_vector())); - dt.insert(intersection_points.begin(), - intersection_points.end()); - - // tangency with the bbox -> no intersection - if (dt.dimension()!=2) return P; - - //now create the polyhedron from the triangulation - internal::Builder_from_T_2< typename Polyhedron::HalfedgeDS,DT > builder(dt); - P.delegate(builder); - - return P; -} - -template -Polyhedron clip_bbox(const Bbox_3& bbox, const Plane_3& plane) -{ - typedef typename Polyhedron::Traits::Kernel Kernel; - typedef typename Kernel::Point_3 Point_3; - typedef typename Kernel::Segment_3 Segment_3; - cpp11::array corners= {{ - Point_3(bbox.xmin(),bbox.ymin(),bbox.zmin()), - Point_3(bbox.xmin(),bbox.ymax(),bbox.zmin()), - Point_3(bbox.xmax(),bbox.ymax(),bbox.zmin()), - Point_3(bbox.xmax(),bbox.ymin(),bbox.zmin()), - Point_3(bbox.xmin(),bbox.ymin(),bbox.zmax()), - Point_3(bbox.xmin(),bbox.ymax(),bbox.zmax()), - Point_3(bbox.xmax(),bbox.ymax(),bbox.zmax()), - Point_3(bbox.xmax(),bbox.ymin(),bbox.zmax()) - }}; - - cpp11::array orientations = {{ - plane.oriented_side(corners[0]), - plane.oriented_side(corners[1]), - plane.oriented_side(corners[2]), - plane.oriented_side(corners[3]), - plane.oriented_side(corners[4]), - plane.oriented_side(corners[5]), - plane.oriented_side(corners[6]), - plane.oriented_side(corners[7]) - }}; - - std::vector points; - // first look for intersections at corners - for (int i=0; i<8; ++i) - if (orientations[i]==ON_ORIENTED_BOUNDARY) - points.push_back(corners[i]); - // second look for intersections on edges - cpp11::array edge_indices = {{ // 2 *12 edges - 0,1, 1,2, 2,3, 3,0, // bottom face edges - 4,5, 5,6, 6,7, 7,4, // top face edges - 0,4, 1,5, 2,6, 3,7 - }}; - - for (int i=0; i<12; ++i) - { - int i1=edge_indices[2*i], i2=edge_indices[2*i+1]; - if (orientations[i1]==ON_ORIENTED_BOUNDARY) continue; - if (orientations[i2]==ON_ORIENTED_BOUNDARY) continue; - if (orientations[i1]!=orientations[i2]) - points.push_back( - boost::get( - *CGAL::intersection(plane, Segment_3(corners[i1], corners[i2]) ) - ) - ); - } - - Polyhedron P; - - //if less that 3 points there will be nothing to clipping. - if (points.size()<3) return P; - - for (int i=0; i<8; ++i) - if (orientations[i]==ON_NEGATIVE_SIDE) - points.push_back(corners[i]); - - // take the convex hull of the points on the negative side+intersection points - // overkill... - CGAL::convex_hull_3(points.begin(), points.end(), P); - - return P; -} - -template -Polyhedron* clip_polyhedron(Polyhedron& P, Polyhedron& clipping_polyhedron) -{ - std::pair result; - typedef CGAL::Polyhedron_corefinement Corefinement; - Corefinement coref; - CGAL::Emptyset_iterator emptyset_iterator; - coref(P, clipping_polyhedron, emptyset_iterator, - &result, Corefinement::Intersection_tag); - - return result.first; -} - -template -Polyhedron* clip_polyhedron(const Polyhedron& P, const Plane_3& p) -{ - if(P.empty()) return new Polyhedron(); - CGAL::Bbox_3 bbox( CGAL::bbox_3(P.points_begin(), P.points_end()) ); - //extend the bbox a bit to avoid border cases - double xd=(bbox.xmax()-bbox.xmin())/100; - double yd=(bbox.ymax()-bbox.ymin())/100; - double zd=(bbox.zmax()-bbox.zmin())/100; - bbox=CGAL::Bbox_3(bbox.xmin()-xd, bbox.ymin()-yd, bbox.zmin()-zd, - bbox.xmax()+xd, bbox.ymax()+yd, bbox.zmax()+zd); - Polyhedron clipping_polyhedron=clip_bbox(bbox, p); - - if (clipping_polyhedron.empty()) //no intersection, result is all or nothing - { - if (p.oriented_side(*P.points_begin())==ON_POSITIVE_SIDE) - return new Polyhedron(); - else - return new Polyhedron(P); - } - Polyhedron copy(P); - return clip_polyhedron(copy, clipping_polyhedron); -} - -template -std::pair split_polyhedron(const Polyhedron& P, const Plane_3& p) -{ - std::pair res; - if(P.empty()) {res.first=res.second= new Polyhedron(); return res;} - CGAL::Bbox_3 bbox( CGAL::bbox_3(P.points_begin(), P.points_end()) ); - //extend the bbox a bit to avoid border cases - double xd=(bbox.xmax()-bbox.xmin())/100; - double yd=(bbox.ymax()-bbox.ymin())/100; - double zd=(bbox.zmax()-bbox.zmin())/100; - bbox=CGAL::Bbox_3(bbox.xmin()-xd, bbox.ymin()-yd, bbox.zmin()-zd, - bbox.xmax()+xd, bbox.ymax()+yd, bbox.zmax()+zd); - //First Polyhedron - Polyhedron clipping_polyhedron=clip_bbox(bbox, p); - - - if (clipping_polyhedron.empty()) //no intersection, result is all or nothing - { - if (p.oriented_side(*P.points_begin())==ON_POSITIVE_SIDE) - { - res.first = new Polyhedron(); res.second = new Polyhedron(P); - } - else - { - res.first = new Polyhedron(P); res.second = new Polyhedron(); - } - return res; - } - - Polyhedron copy(P); - res.first = clip_polyhedron(copy, clipping_polyhedron); - //Second Polyhedron - Plane_3 p_op = p.opposite(); - clipping_polyhedron=clip_bbox(bbox, p_op); - - copy = P; - res.second = clip_polyhedron(copy, clipping_polyhedron); - return res; - -} - -namespace internal{ - -template -struct Edge_is_marked4coref{ - std::set& marked_halfedges; - typedef bool value_type; - typedef value_type reference; - typedef std::pair key_type; - typedef boost::read_write_property_map_tag category; - - Edge_is_marked4coref(std::set& mh) - : marked_halfedges(mh) - {} - - friend reference get(Edge_is_marked4coref& map,const key_type& key) { - return map.marked_halfedges.count(key.first)!=0; - } - friend void put(Edge_is_marked4coref& map,key_type key,value_type v) { - if (v) map.marked_halfedges.insert(key.first); - else map.marked_halfedges.erase(key.first); - } -}; - -template -struct Edge_is_marked{ - const std::set* marked_halfedges; - typedef bool value_type; - typedef value_type reference; - typedef typename boost::graph_traits::edge_descriptor key_type; - typedef boost::readable_property_map_tag category; - - Edge_is_marked(){} - Edge_is_marked(const std::set& mh) - : marked_halfedges(&mh) - {} - - friend reference get(const Edge_is_marked& map,const key_type& key) { - return map.marked_halfedges->count(key.halfedge())!=0; - } -}; - -} //end of internal namespace - -template -void inplace_clip_open_polyhedron(Polyhedron& P, const Plane_3& p) -{ - CGAL::Bbox_3 bbox( CGAL::bbox_3(P.points_begin(), P.points_end()) ); - //extend the bbox a bit to avoid border cases - double xd=(bbox.xmax()-bbox.xmin())/100; - double yd=(bbox.ymax()-bbox.ymin())/100; - double zd=(bbox.zmax()-bbox.zmin())/100; - bbox=CGAL::Bbox_3(bbox.xmin()-xd, bbox.ymin()-yd, bbox.zmin()-zd, - bbox.xmax()+xd, bbox.ymax()+yd, bbox.zmax()+zd); - Polyhedron clipping_polyhedron=clip_to_bbox(bbox, p); - - if (clipping_polyhedron.empty()) //no intersection, result is all or nothing - { - if (p.oriented_side(*P.points_begin())==ON_POSITIVE_SIDE) - P.clear(); - return; - } - - // set for marking edges of P intersected by the clipping plane - std::set marked_halfedges; - internal::Edge_is_marked4coref cr_edge_is_marked(marked_halfedges); - typedef typename Polyhedron::Traits::Kernel K; - typedef CGAL::Node_visitor_refine_polyhedra< - Polyhedron,CGAL::Default,K,internal::Edge_is_marked4coref > - Split_visitor; - CGAL::Default_output_builder output_builder; - - Split_visitor visitor(output_builder, - Default_polyhedron_ppmap(), - cr_edge_is_marked); - - CGAL::Intersection_of_Polyhedra_3 - polyline_intersections(visitor); - CGAL::Emptyset_iterator emptyset_iterator; - // corefinement P and clipping_polyhedron - polyline_intersections(P,clipping_polyhedron,emptyset_iterator); - - // extract connected components bounded by marked edges - internal::Edge_is_marked edge_is_marked(marked_halfedges); - namespace PMP=Polygon_mesh_processing; - std::map face_ccs; - std::size_t nb_cc=PMP::connected_components(P, - boost::make_assoc_property_map(face_ccs), - PMP::parameters::edge_is_constrained_map(edge_is_marked) - .face_index_map(get(boost::face_external_index,P)) - ); - - // remove cc on the positive side of the plane - std::vector cc_handled(nb_cc, false); - std::vector ccs_to_remove; - BOOST_FOREACH(typename Polyhedron::Face_handle f, faces(P)) - { - std::size_t cc_id=face_ccs[f]; - if (cc_handled[cc_id]) continue; - - //look for a vertex not on the intersection - typename Polyhedron::Halfedge_handle h=f->halfedge(); - for(int i=0;i<3;++i){ - bool no_marked_edge=true; - BOOST_FOREACH(typename Polyhedron::Halfedge_handle h, halfedges_around_target(h, P)) - if ( marked_halfedges.count(h) ) - no_marked_edge=false; - if (no_marked_edge){ - if ( p.oriented_side(h->vertex()->point())==ON_POSITIVE_SIDE ) - ccs_to_remove.push_back(cc_id); - cc_handled[cc_id]=true; - if (--nb_cc==0) break; - break; - } - h=h->next(); - } - } - - //now remove the faces on the positive side - PMP::remove_connected_components(P, - ccs_to_remove, - boost::make_assoc_property_map(face_ccs), - PMP::parameters::vertex_index_map(get(boost::vertex_external_index,P)) - ); -} - -} } // CGAL::corefinement - - -#endif // CGAL_INTERNAL_POLYHEDRON_PLANE_CLIPPING_3_H diff --git a/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/internal/clip.h b/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/internal/clip.h new file mode 100644 index 00000000000..80922c497e7 --- /dev/null +++ b/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/internal/clip.h @@ -0,0 +1,367 @@ +// Copyright (c) 2016 GeometryFactory (France). +// All rights reserved. +// +// This file is part of CGAL (www.cgal.org). +// You can redistribute it and/or modify it under the terms of the GNU +// General Public License as published by the Free Software Foundation, +// either version 3 of the License, or (at your option) any later version. +// +// Licensees holding a valid commercial license may use this file in +// accordance with the commercial license agreement provided with the software. +// +// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE +// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. +// +// $URL$ +// $Id$ +// +// +// Author(s) : Sebastien Loriot + +#ifndef CGAL_POLYGON_MESH_PROCESSING_CLIP_H +#define CGAL_POLYGON_MESH_PROCESSING_CLIP_H + +#include +#include +#include + +#include + +namespace CGAL{ +namespace Polygon_mesh_processing { + +namespace internal +{ +template +bool +clip_open_impl( TriangleMesh& tm, + TriangleMesh& clipper, + Ecm ecm, + const NamedParameters1& np_tm, + const NamedParameters2& np_c) +{ + // first corefine the meshes + corefine(tm, clipper, np_tm, np_c); + + typedef typename GetVertexPointMap::type Vpm; + typedef typename GetFaceIndexMap::type Fid_map; + typedef typename GetVertexIndexMap::type Vid_map; + + typedef boost::graph_traits GT; + typedef typename GetGeomTraits::type GeomTraits; + typedef typename GT::halfedge_descriptor halfedge_descriptor; + typedef typename GT::face_descriptor face_descriptor; + + Fid_map fid_map = boost::choose_param(get_param(np_tm, face_index), + get_property_map(face_index, tm)); + Vid_map vid_map = boost::choose_param(get_param(np_tm, boost::vertex_index), + get_property_map(boost::vertex_index, tm)); + Vpm vpm1 = boost::choose_param(get_param(np_tm, vertex_point), + get_property_map(vertex_point, tm)); + Vpm vpm_c = boost::choose_param(get_param(np_c, vertex_point), + get_property_map(vertex_point, clipper)); + + // init indices if needed + Corefinement::init_face_indices(tm, fid_map); + Corefinement::init_vertex_indices(tm, vid_map); + + // set the connected component id of each face + std::vector face_cc(num_faces(tm), std::size_t(-1)); + std::size_t nb_cc = + connected_components(tm, + Corefinement::bind_maps(fid_map, make_property_map(face_cc)), + parameters::face_index_map(fid_map). + edge_is_constrained_map(ecm)); + + + boost::dynamic_bitset<> cc_not_handled(nb_cc); + cc_not_handled.set(); + std::vector ccs_to_remove; + /// \todo clipper has been modified, this is not robust if inexact constructions are used + CGAL::Side_of_triangle_mesh + side_of(clipper, vpm_c); + BOOST_FOREACH(face_descriptor f, faces(tm)) + { + std::size_t cc_id = face_cc[ get(fid_map, f) ]; + if ( !cc_not_handled.test(cc_id) ) continue; + + halfedge_descriptor h=halfedge(f, tm); + for(int i=0;i<3;++i) + { + bool no_marked_edge=true; + BOOST_FOREACH(halfedge_descriptor h2, halfedges_around_target(h, tm)) + if ( get(ecm, edge(h2, tm)) ){ + no_marked_edge=false; + break; + } + if (no_marked_edge){ + if ( side_of( get(vpm1, target(h, tm) ) ) == ON_UNBOUNDED_SIDE ) + ccs_to_remove.push_back(cc_id); + cc_not_handled.reset(cc_id); + break; + } + h=next(h, tm); + } + if (!cc_not_handled.any()) break; + } + + if (cc_not_handled.any()) + { + ///\todo handle using barycenters? won't work for coplanar faces + } + //now remove the cc + remove_connected_components(tm, + ccs_to_remove, + Corefinement::bind_maps(fid_map, make_property_map(face_cc)), + np_tm); + + return true; +} + +/// \TODO move this to property_map.h? +template +struct Constrained_edge_map +{ + typedef boost::read_write_property_map_tag category; + typedef bool value_type; + typedef bool reference; + typedef typename Set::value_type key_type; + + Constrained_edge_map() + : edge_set(NULL) + {} + + Constrained_edge_map(Set& set) + : edge_set(&set) + {} + + friend bool get(const Constrained_edge_map& map, key_type k) + { + return map.edge_set->count(k)!=0; + } + + friend void put(Constrained_edge_map& map, key_type k, bool b) + { + if (b) + map.edge_set->insert(k); + else + map.edge_set->erase(k); + } +private: + Set* edge_set; +}; + +template +bool +clip_open_impl( TriangleMesh& tm, + TriangleMesh& clipper, + boost::param_not_found, + const NamedParameters1& np_tm, + const NamedParameters2& np_c) +{ + typedef typename boost::graph_traits + ::edge_descriptor edge_descriptor; + boost::unordered_set constrained_edges; + Constrained_edge_map > + cst_map(constrained_edges); + + return clip_open_impl(tm, clipper, + cst_map, + np_tm.edge_is_constrained_map(cst_map), + np_c); +} + +} // end of internal namespace + +#ifndef DOXYGEN_RUNNING + +///\todo clipper const! +/// requires face_index_map, vertex_index_map for np_tm +/// requires face_index_map for np_c +/// if edge_is_constrained_map is not provided in np_tm a default one is +/// provided using boost::unordered_set +template +bool +clip( TriangleMesh& tm, + /*const*/ TriangleMesh& clipper, + bool close, + const NamedParameters1& np_tm, + const NamedParameters2& np_c) +{ + if (close && is_closed(tm)) + return corefine_and_compute_intersection(tm, clipper, tm, np_tm, np_c); + + return internal::clip_open_impl(tm, clipper, + get_param(np_tm, edge_is_constrained), np_tm, np_c); +} + +/// \todo document me +template +Oriented_side +clip_to_bbox(const Plane_3& plane, + const Bbox_3& bbox, + TriangleMesh& tm_out, + const NamedParameters& np ) +{ + typedef typename GetGeomTraits::type Geom_traits; + typedef typename Geom_traits::Point_3 Point_3; + typedef typename Geom_traits::Segment_3 Segment_3; + typedef typename GetVertexPointMap::type Vpm; + + Vpm vpm_out = boost::choose_param(get_param(np, vertex_point), + get_property_map(vertex_point, tm_out)); + + + cpp11::array corners= {{ + Point_3(bbox.xmin(),bbox.ymin(),bbox.zmin()), + Point_3(bbox.xmin(),bbox.ymax(),bbox.zmin()), + Point_3(bbox.xmax(),bbox.ymax(),bbox.zmin()), + Point_3(bbox.xmax(),bbox.ymin(),bbox.zmin()), + Point_3(bbox.xmin(),bbox.ymin(),bbox.zmax()), + Point_3(bbox.xmin(),bbox.ymax(),bbox.zmax()), + Point_3(bbox.xmax(),bbox.ymax(),bbox.zmax()), + Point_3(bbox.xmax(),bbox.ymin(),bbox.zmax()) + }}; + + cpp11::array orientations = {{ + plane.oriented_side(corners[0]), + plane.oriented_side(corners[1]), + plane.oriented_side(corners[2]), + plane.oriented_side(corners[3]), + plane.oriented_side(corners[4]), + plane.oriented_side(corners[5]), + plane.oriented_side(corners[6]), + plane.oriented_side(corners[7]) + }}; + + std::vector points; + + // look for intersections on edges + cpp11::array edge_indices = {{ // 2 *12 edges + 0,1, 1,2, 2,3, 3,0, // bottom face edges + 4,5, 5,6, 6,7, 7,4, // top face edges + 0,4, 1,5, 2,6, 3,7 + }}; + + for (int i=0; i<12; ++i) + { + int i1=edge_indices[2*i], i2=edge_indices[2*i+1]; + if (orientations[i1]==ON_ORIENTED_BOUNDARY) continue; + if (orientations[i2]==ON_ORIENTED_BOUNDARY) continue; + if (orientations[i1]!=orientations[i2]) + points.push_back( + boost::get( + *intersection(plane, Segment_3(corners[i1], corners[i2]) ) + ) + ); + } + + + Oriented_side last_os = ON_ORIENTED_BOUNDARY; + for (int i=0; i<8; ++i) + if (orientations[i]!=ON_ORIENTED_BOUNDARY) + { + if (last_os==ON_ORIENTED_BOUNDARY) + last_os=orientations[i]; + else + { + if(last_os!=orientations[i]) + { + last_os=ON_ORIENTED_BOUNDARY; + break; + } + } + } + + // the intersection is the full bbox + if (last_os!=ON_ORIENTED_BOUNDARY) + return last_os; + + //add points on negative side and on the plane + for (int i=0; i<8; ++i) + if (orientations[i]!=ON_POSITIVE_SIDE) + points.push_back(corners[i]); + + // take the convex hull of the points on the negative side+intersection points + // overkill... + Polyhedron_3 P; + CGAL::convex_hull_3(points.begin(), points.end(), P); + copy_face_graph(P, tm_out, + Emptyset_iterator(), Emptyset_iterator(), Emptyset_iterator(), + get(vertex_point, P), vpm_out); + return ON_ORIENTED_BOUNDARY; +} + + +// convenience overload +template +bool +clip( TriangleMesh& tm, + /*const*/ TriangleMesh& clipper, + bool close, + const NamedParameters1& np_tm) +{ + return clip(tm, clipper, close, np_tm, parameters::all_default()); +} + +// convenience overload +template +bool +clip( TriangleMesh& tm, + /*const*/ TriangleMesh& clipper, + bool close) +{ + return clip(tm, clipper, close, parameters::all_default()); +} + +// works only with the default point map, for more complex use cases, use +// clip_to_bbox first and the other overload of clip with two meshes +/// \todo document me +template +void clip( TriangleMesh& tm, + const Plane_3& plane, + bool close) +{ + if( boost::begin(faces(tm))==boost::end(faces(tm)) ) return; + CGAL::Bbox_3 bbox = bbox_3(tm); + //extend the bbox a bit to avoid border cases + double xd=(bbox.xmax()-bbox.xmin())/100; + double yd=(bbox.ymax()-bbox.ymin())/100; + double zd=(bbox.zmax()-bbox.zmin())/100; + bbox=CGAL::Bbox_3(bbox.xmin()-xd, bbox.ymin()-yd, bbox.zmin()-zd, + bbox.xmax()+xd, bbox.ymax()+yd, bbox.zmax()+zd); + TriangleMesh clipper; + Oriented_side os = clip_to_bbox(plane, bbox, clipper, parameters::all_default()); + + switch(os) + { + case ON_NEGATIVE_SIDE: + return; // nothing to clip, the full mesh is on the negative side + case ON_POSITIVE_SIDE: + clear(tm); // clear the mesh that is fully on the positive side + return; + default: + clip(tm, clipper, close); + } +} + +#endif // !DOXYGEN_RUNNING + +} } //end of namespace CGAL::Polygon_mesh_processing + +#endif // CGAL_POLYGON_MESH_PROCESSING_CLIP_H diff --git a/Polygon_mesh_processing/test/Polygon_mesh_processing/CMakeLists.txt b/Polygon_mesh_processing/test/Polygon_mesh_processing/CMakeLists.txt index 6d50031616e..97cc3a3cacf 100644 --- a/Polygon_mesh_processing/test/Polygon_mesh_processing/CMakeLists.txt +++ b/Polygon_mesh_processing/test/Polygon_mesh_processing/CMakeLists.txt @@ -81,7 +81,7 @@ if (EIGEN3_FOUND) create_single_source_cgal_program("test_corefinement_bool_op.cpp") create_single_source_cgal_program("test_corefine.cpp") create_single_source_cgal_program("test_does_bound_a_volume.cpp") - create_single_source_cgal_program("polyhedron_clipping.cpp") + create_single_source_cgal_program("test_pmp_clip.cpp") if(NOT (${EIGEN3_VERSION} VERSION_LESS 3.2.0)) create_single_source_cgal_program("fairing_test.cpp") diff --git a/Polygon_mesh_processing/test/Polygon_mesh_processing/polyhedron_clipping.cpp b/Polygon_mesh_processing/test/Polygon_mesh_processing/polyhedron_clipping.cpp deleted file mode 100644 index 85cbab8f474..00000000000 --- a/Polygon_mesh_processing/test/Polygon_mesh_processing/polyhedron_clipping.cpp +++ /dev/null @@ -1,67 +0,0 @@ -#include -#include -#include -#include - -#include -#include - -typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel; -typedef CGAL::Polyhedron_3 Polyhedron; - -int main(int argc, char** argv) -{ - const char* fname = argc>1?argv[1]:"data/elephant.off"; - std::ifstream input(fname); - - Polyhedron P; - input >> P; - std::cout << "Polyhedron with " << P.size_of_vertices() << " vertices\n"; -// test clipping by a plane -{ - Polyhedron* P_clipped= - CGAL::corefinement::clip_polyhedron(P, Kernel::Plane_3(1,1,0,0)); - std::ofstream output("clipped1.off"); - output << *P_clipped; - delete P_clipped; -} -// test clipping with the opposite plane -{ - Polyhedron* P_clipped= - CGAL::corefinement::clip_polyhedron(P, Kernel::Plane_3(-1,-1,0,0)); - std::ofstream output("clipped2.off"); - output << *P_clipped; - delete P_clipped; -} -// test clipping with a plane with no intersection having the elephant -// on its positive side, result should be empty -{ - Polyhedron* P_clipped= - CGAL::corefinement::clip_polyhedron(P, Kernel::Plane_3(1,0,0,1)); - if (!P_clipped->empty()){ - std::cerr << "Error: Polyhedron should be empty!\n"; - return 1; - } - delete P_clipped; -} -// test clipping with a plane with no intersection having the elephant -// on its negative side, result should be empty -{ - Polyhedron* P_clipped= - CGAL::corefinement::clip_polyhedron(P, Kernel::Plane_3(-1,0,0,-1)); - if (P_clipped->size_of_vertices()!=P.size_of_vertices()){ - std::cerr << "Error: Polyhedron should be full!\n"; - return 1; - } - delete P_clipped; -} - -// clip with a plane but do not close the polyhedron -{ - CGAL::corefinement::inplace_clip_open_polyhedron(P, Kernel::Plane_3(1,1,0,0)); - std::ofstream output("open_clipped.off"); - output << P; -} - - return 0; -} diff --git a/Polygon_mesh_processing/test/Polygon_mesh_processing/test_pmp_clip.cpp b/Polygon_mesh_processing/test/Polygon_mesh_processing/test_pmp_clip.cpp new file mode 100644 index 00000000000..9a9d24e4ee7 --- /dev/null +++ b/Polygon_mesh_processing/test/Polygon_mesh_processing/test_pmp_clip.cpp @@ -0,0 +1,65 @@ +#include +#include +#include +#include + +#include +#include + +namespace PMP = CGAL::Polygon_mesh_processing; +namespace params = PMP::parameters; + +typedef CGAL::Exact_predicates_inexact_constructions_kernel K; +typedef CGAL::Surface_mesh Triangle_mesh; +typedef CGAL::Polyhedron_3 Polyhedron; +typedef Triangle_mesh::Property_map Constrained_edge_map; + + +int main() +{ + { + // test open clipping with Surface_mesh + Triangle_mesh tm1, tm2; + std::ifstream input("data-coref/elephant.off"); + input >> tm1; + input.close(); + input.open("data-coref/sphere.off"); + input >> tm2; + input.close(); + + Constrained_edge_map ecm1 = + tm1.add_property_map("e:cst", false).first; + + PMP::clip(tm1, tm2, false, params::edge_is_constrained_map(ecm1)); + std::ofstream output("clipped_opened.off"); + output << tm1; + + // test open clipping with Polyhedron + Polyhedron P, Q; + input.open("data-coref/elephant.off"); + input >> P; + input.close(); + input.open("data-coref/sphere.off"); + input >> Q; + + PMP::clip(P, Q, false, + params::face_index_map(get(CGAL::face_external_index, P)). + vertex_index_map(get(CGAL::vertex_external_index, P)), + params::face_index_map(get(CGAL::face_external_index, Q))); + assert(P.size_of_vertices() == tm1.number_of_vertices()); + } + { + Triangle_mesh tm1, tm2; + std::ifstream input("data-coref/elephant.off"); + input >> tm1; + input.close(); + input.open("data-coref/sphere.off"); + input >> tm2; + + PMP::clip(tm1, tm2, true); + std::ofstream output("clipped_closed.off"); + output << tm1; + } + + return 0; +} diff --git a/Polyhedron/demo/Polyhedron/Plugins/Operations_on_polyhedra/Clip_polyhedron_plugin.cpp b/Polyhedron/demo/Polyhedron/Plugins/Operations_on_polyhedra/Clip_polyhedron_plugin.cpp index 060bf93b19c..3f87c1c71f8 100644 --- a/Polyhedron/demo/Polyhedron/Plugins/Operations_on_polyhedra/Clip_polyhedron_plugin.cpp +++ b/Polyhedron/demo/Polyhedron/Plugins/Operations_on_polyhedra/Clip_polyhedron_plugin.cpp @@ -7,7 +7,8 @@ #include "Scene_plane_item.h" #include #include -#include +#include + #include "ui_Clip_polyhedron_plugin.h" #include "Viewer.h" @@ -187,100 +188,43 @@ public Q_SLOTS: //apply the clipping function Q_FOREACH(Scene_polyhedron_item* poly, polyhedra) { - - if(ui_widget.close_checkBox->isChecked() && poly->polyhedron()->is_closed()) + if (ui_widget.clip_radioButton->isChecked()) { - std::pairpolyhedron; - if(ui_widget.clip_radioButton->isChecked()) - { - polyhedron.first = CGAL::corefinement::clip_polyhedron(*(poly->polyhedron()),plane->plane()); - polyhedron.second = NULL; - } - else - { - polyhedron = CGAL::corefinement::split_polyhedron(*(poly->polyhedron()),plane->plane()); - } - if(polyhedron.first != NULL) - { - Scene_polyhedron_item* new_item = new Scene_polyhedron_item(polyhedron.first); - if(polyhedron.second != NULL) - new_item->setName(QString("%1 %2").arg(poly->name()).arg("1")); - else - new_item->setName(poly->name()); - new_item->setColor(poly->color()); - new_item->setRenderingMode(poly->renderingMode()); - new_item->setVisible(poly->visible()); - new_item->invalidateOpenGLBuffers(); - new_item->setProperty("source filename", poly->property("source filename")); - scene->replaceItem(scene->item_id(poly),new_item); - new_item->invalidateOpenGLBuffers(); - viewer->updateGL(); - if(ui_widget.clip_radioButton->isChecked()) - messages->information(QString("%1 clipped").arg(new_item->name())); - } - else - { - messages->information(QString("Could not clip %1 : returned polyhedron is null.").arg(poly->name())); - delete polyhedron.first; - } - if(polyhedron.second!= NULL) - { - Scene_polyhedron_item* new_item = new Scene_polyhedron_item(polyhedron.second); - new_item->setName(QString("%1 %2").arg(poly->name()).arg("2")); - new_item->setColor(poly->color()); - new_item->setRenderingMode(poly->renderingMode()); - new_item->setVisible(poly->visible()); - new_item->invalidateOpenGLBuffers(); - new_item->setProperty("source filename", poly->property("source filename")); - scene->addItem(new_item); - new_item->invalidateOpenGLBuffers(); - viewer->updateGL(); - if(!ui_widget.clip_radioButton->isChecked()) - messages->information(QString("%1 splitted").arg(poly->name())); - } - if(polyhedron.first != NULL) - delete poly; + CGAL::Polygon_mesh_processing::clip(*(poly->polyhedron()), + plane->plane(), + ui_widget.close_checkBox->isChecked()); + poly->invalidateOpenGLBuffers(); + viewer->updateGL(); } - else { - Scene_polyhedron_item* poly1 = new Scene_polyhedron_item(*(poly->polyhedron())); - poly1->setProperty("source filename", poly->property("source filename")); - Scene_polyhedron_item* poly2 = NULL; - if(!ui_widget.clip_radioButton->isChecked()) - { - poly2 = new Scene_polyhedron_item(*(poly->polyhedron())); - poly2->setProperty("source filename", poly->property("source filename")); - } + //part on the negative side + Polyhedron* neg_side = new Polyhedron(*poly->polyhedron()); + CGAL::Polygon_mesh_processing::clip(*neg_side, + plane->plane(), + ui_widget.close_checkBox->isChecked()); + Scene_polyhedron_item* new_item = new Scene_polyhedron_item(neg_side); + new_item->setName(QString("%1 on %2").arg(poly->name()).arg("negative side")); + new_item->setColor(poly->color()); + new_item->setRenderingMode(poly->renderingMode()); + new_item->setVisible(poly->visible()); + scene->addItem(new_item); + new_item->invalidateOpenGLBuffers(); + // part on the positive side + Polyhedron* pos_side = new Polyhedron(*poly->polyhedron()); + CGAL::Polygon_mesh_processing::clip(*pos_side, + plane->plane().opposite(), + ui_widget.close_checkBox->isChecked()); + new_item = new Scene_polyhedron_item(pos_side); + new_item->setName(QString("%1 on %2").arg(poly->name()).arg("positive side")); + new_item->setColor(poly->color()); + new_item->setRenderingMode(poly->renderingMode()); + new_item->setVisible(poly->visible()); + scene->addItem(new_item); + new_item->invalidateOpenGLBuffers(); - CGAL::corefinement::inplace_clip_open_polyhedron(*(poly1->polyhedron()),plane->plane()); - if(poly2 != NULL) - poly1->setName(QString("%1 %2").arg(poly->name()).arg("1")); - else - poly1->setName(poly->name()); - poly1->setColor(poly->color()); - poly1->setRenderingMode(poly->renderingMode()); - poly1->setVisible(poly->visible()); - scene->replaceItem(scene->item_id(poly),poly1); - poly1->invalidateOpenGLBuffers(); - viewer->updateGL(); - if(ui_widget.clip_radioButton->isChecked()) - messages->information(QString("%1 clipped").arg(poly->name())); - - if(poly2 != NULL) - { - CGAL::corefinement::inplace_clip_open_polyhedron(*(poly2->polyhedron()),plane->plane().opposite()); - poly2->setName(QString("%1 %2").arg(poly->name()).arg("2")); - poly2->setColor(poly->color()); - poly2->setRenderingMode(poly->renderingMode()); - poly2->setVisible(poly->visible()); - scene->addItem(poly2); - poly2->invalidateOpenGLBuffers(); - viewer->updateGL(); - messages->information(QString("%1 splitted").arg(poly->name())); - } - delete poly; + viewer->updateGL(); } } }