From d3dc9be69c28c4efb56ef8b7f25fefd142f671cc Mon Sep 17 00:00:00 2001 From: Aymeric PELLE Date: Wed, 16 Apr 2014 20:05:53 +0200 Subject: [PATCH] Adds two compatibilty headers. CGAL/Mesh_3/Labeled_mesh_domain_3.h CGAL/Mesh_3/Implicit_to_labeled_function_wrapper.h So, some branches in development using theses features won't be broken. --- .../Implicit_to_labeled_function_wrapper.h | 153 +++++ .../CGAL/Mesh_3/Labeled_mesh_domain_3.h | 600 ++++++++++++++++++ 2 files changed, 753 insertions(+) create mode 100644 Mesh_3/include/CGAL/Mesh_3/Implicit_to_labeled_function_wrapper.h create mode 100644 Mesh_3/include/CGAL/Mesh_3/Labeled_mesh_domain_3.h diff --git a/Mesh_3/include/CGAL/Mesh_3/Implicit_to_labeled_function_wrapper.h b/Mesh_3/include/CGAL/Mesh_3/Implicit_to_labeled_function_wrapper.h new file mode 100644 index 00000000000..ccfd64789bd --- /dev/null +++ b/Mesh_3/include/CGAL/Mesh_3/Implicit_to_labeled_function_wrapper.h @@ -0,0 +1,153 @@ +// Copyright (c) 2009 INRIA Sophia-Antipolis (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) : Stéphane Tayeb +// +//****************************************************************************** +// File Description : +// Implicit_to_labeling_function_wrapper and +// Implicit_vector_to_labeling_function_wrapper class declaration +// and implementation. +// +// See classes description to have more information. +//****************************************************************************** + +#ifndef CGAL_MESH_3_IMPLICIT_TO_LABELED_FUNCTION_WRAPPER_H +#define CGAL_MESH_3_IMPLICIT_TO_LABELED_FUNCTION_WRAPPER_H + +#if defined(BOOST_MSVC) +# pragma warning(push) +# pragma warning(disable:4180) // qualifier applied to function type has no meaning; ignored +#endif + +#warning deprecated : use instead. + +#include + +namespace CGAL { + +namespace Mesh_3 { + +#include + +/** + * @class Implicit_to_labeled_function_wrapper + * + * This class is designed to wrap an implicit function which describes a domain + * by [p is inside if f(p)<0] to a function which takes its values into {0,1}. + * f(p)=0 means that p is outside the domain. + */ +template +class Implicit_to_labeled_function_wrapper +{ +public: + // Types + typedef int return_type; + typedef typename BGT::Point_3 Point_3; + + /// Constructor + Implicit_to_labeled_function_wrapper(const Function_& f) + : r_f_(f) {} + + // Default copy constructor and assignment operator are ok + + /// Destructor + ~Implicit_to_labeled_function_wrapper() {} + + /// Operator () + return_type operator()(const Point_3& p, const bool = true) const + { + return ( (r_f_(p)<0) ? 1 : 0 ); + } + +private: + /// Function to wrap + const Function_& r_f_; + +}; // end class Implicit_to_labeled_function_wrapper + + + +/** + * \deprecated + * + * @class Implicit_vector_to_labeled_function_wrapper + * + * Wraps a set of implicit function [f1,f2,...] to one function F which + * takes its values into N. + * + * Let p be a point. + * F(p) = 0b000000(f2(p)<0)(f1(p)<0) + * + * It can handle at most 8 functions. + */ +template +class Implicit_vector_to_labeled_function_wrapper +{ +public: + // Types + typedef int return_type; + typedef std::vector Function_vector; + typedef typename BGT::Point_3 Point_3; + + /// Constructor + Implicit_vector_to_labeled_function_wrapper(const std::vector& v) + : function_vector_(v) {} + + // Default copy constructor and assignment operator are ok + + /// Destructor + ~Implicit_vector_to_labeled_function_wrapper() {} + + /// Operator () + return_type operator()(const Point_3& p, const bool = true) const + { + int nb_func = function_vector_.size(); + if ( nb_func > 8 ) + { + CGAL_error_msg("We support at most 8 functions !"); + } + + char bits = 0; + for ( int i = 0 ; i < nb_func ; ++i ) + { + // Insert value into bits : we compute fi(p) and insert result at + // bit i of bits + bits |= ( ((*function_vector_[i])(p) < 0) << i ); + } + + return ( static_cast(bits) ); + } + +private: + /// Functions to wrap + const Function_vector function_vector_; + +}; // end class Implicit_to_labeled_function_wrapper + +} // end namespace Mesh_3 + +} // end namespace CGAL + + + +#if defined(BOOST_MSVC) +# pragma warning(pop) +#endif + +#endif // CGAL_MESH_3_IMPLICIT_TO_LABELED_FUNCTION_WRAPPER_H diff --git a/Mesh_3/include/CGAL/Mesh_3/Labeled_mesh_domain_3.h b/Mesh_3/include/CGAL/Mesh_3/Labeled_mesh_domain_3.h new file mode 100644 index 00000000000..7391ac85b98 --- /dev/null +++ b/Mesh_3/include/CGAL/Mesh_3/Labeled_mesh_domain_3.h @@ -0,0 +1,600 @@ +// Copyright (c) 2009 INRIA Sophia-Antipolis (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) : Stéphane Tayeb +// +//****************************************************************************** +// File Description : +// class Labeled_mesh_domain_3. See class description. +//****************************************************************************** + +#ifndef CGAL_MESH_3_LABELED_MESH_DOMAIN_3_H +#define CGAL_MESH_3_LABELED_MESH_DOMAIN_3_H + +#warning deprecated : use instead. + +#include + +#include + +#include +#include + +#include +#include +#include +#include +#include + +namespace CGAL { + +namespace Mesh_3 { + +/** + * \class Labeled_mesh_domain_3 + * + * Function f must take his values into N. + * Let p be a Point. + * - f(p)=0 means that p is outside domain. + * - f(p)=a, a!=0 means that p is inside subdomain a. + * + * Any boundary facet is labelled , a, where b!=0. + */ +template +class Labeled_mesh_domain_3 +{ +public: + /// Geometric object types + typedef typename BGT::Point_3 Point_3; + typedef typename BGT::Segment_3 Segment_3; + typedef typename BGT::Ray_3 Ray_3; + typedef typename BGT::Line_3 Line_3; + typedef typename BGT::Vector_3 Vector_3; + typedef typename BGT::Sphere_3 Sphere_3; + typedef CGAL::Bbox_3 Bbox_3; + +protected: + typedef typename BGT::Iso_cuboid_3 Iso_cuboid_3; + +public: + // Kernel_traits compatibility + typedef BGT R; + + //------------------------------------------------------- + // Index Types + //------------------------------------------------------- + /// Type of indexes for cells of the input complex + typedef typename Function::return_type Subdomain_index; + typedef boost::optional Subdomain; + /// Type of indexes for surface patch of the input complex + typedef std::pair Surface_patch_index; + typedef boost::optional Surface_patch; + /// Type of indexes to characterize the lowest dimensional face of the input + /// complex on which a vertex lie + typedef boost::variant Index; + typedef CGAL::cpp11::tuple Intersection; + + + typedef typename BGT::FT FT; + typedef BGT Geom_traits; + + /** + * @brief Constructor + */ + Labeled_mesh_domain_3(const Function& f, + const Sphere_3& bounding_sphere, + const FT& error_bound = FT(1e-3)); + + Labeled_mesh_domain_3(const Function& f, + const Bbox_3& bbox, + const FT& error_bound = FT(1e-3)); + + /// Destructor + virtual ~Labeled_mesh_domain_3() {} + + /** + * Constructs a set of \ccc{n} points on the surface, and output them to + * the output iterator \ccc{pts} whose value type is required to be + * \ccc{std::pair}. + */ + struct Construct_initial_points + { + Construct_initial_points(const Labeled_mesh_domain_3& domain) + : r_domain_(domain) {} + + template + OutputIterator operator()(OutputIterator pts, const int n = 12) const; + + private: + const Labeled_mesh_domain_3& r_domain_; + }; + + /// Returns Construct_initial_points object + Construct_initial_points construct_initial_points_object() const + { + return Construct_initial_points(*this); + } + + /** + * Returns true if point~\ccc{p} is in the domain. If \ccc{p} is in the + * domain, the parameter index is set to the index of the subdomain + * including $p$. It is set to the default value otherwise. + */ + struct Is_in_domain + { + Is_in_domain(const Labeled_mesh_domain_3& domain) : r_domain_(domain) {} + + Subdomain operator()(const Point_3& p) const + { + // f(p)==0 means p is outside the domain + Subdomain_index index = (r_domain_.function_)(p); + if ( Subdomain_index() == index ) + return Subdomain(); + else + return Subdomain(index); + } + private: + const Labeled_mesh_domain_3& r_domain_; + }; + + /// Returns Is_in_domain object + Is_in_domain is_in_domain_object() const { return Is_in_domain(*this); } + + /** + * Returns true is the element \ccc{type} intersect properly any of the + * surface patches describing the either the domain boundary or some + * subdomain boundary. + * \ccc{Type} is either \ccc{Segment_3}, \ccc{Ray_3} or \ccc{Line_3}. + * Parameter index is set to the index of the intersected surface patch + * if \ccc{true} is returned and to the default \ccc{Surface_patch_index} + * value otherwise. + */ + struct Do_intersect_surface + { + Do_intersect_surface(const Labeled_mesh_domain_3& domain) + : r_domain_(domain) {} + + Surface_patch operator()(const Segment_3& s) const + { + return this->operator()(s.source(), s.target()); + } + + Surface_patch operator()(const Ray_3& r) const + { + return clip_to_segment(r); + } + + Surface_patch operator()(const Line_3& l) const + { + return clip_to_segment(l); + } + + private: + /// Returns true if points \c a & \c b do not belong to the same subdomain + /// \c index is set to the surface index of subdomains f(a), f(b) + Surface_patch operator()(const Point_3& a, const Point_3& b) const + { + // If f(a) != f(b), then [a,b] intersects some surface. Here we consider + // [a,b] intersects surface_patch labelled (or ). + // It may be false, further rafinement will improve precision + const Subdomain_index value_a = r_domain_.function_(a); + const Subdomain_index value_b = r_domain_.function_(b); + + if ( value_a != value_b ) + return Surface_patch(r_domain_.make_surface_index(value_a, value_b)); + else + return Surface_patch(); + } + + /** + * Clips \c query to a segment \c s, and call operator()(s) + */ + template + Surface_patch clip_to_segment(const Query& query) const + { + typename cpp11::result_of::type + clipped = CGAL::intersection(query, r_domain_.bbox_); + + if(clipped) +#if CGAL_INTERSECTION_VERSION > 1 + if(const Segment_3* s = boost::get(&*clipped)) + return this->operator()(*s); +#else + if(const Segment_3* s = object_cast(&clipped)) + return this->operator()(*s); +#endif + + return Surface_patch(); + } + + private: + const Labeled_mesh_domain_3& r_domain_; + }; + + /// Returns Do_intersect_surface object + Do_intersect_surface do_intersect_surface_object() const + { + return Do_intersect_surface(*this); + } + + /** + * Returns a point in the intersection of the primitive \ccc{type} + * with some boundary surface. + * \ccc{Type1} is either \ccc{Segment_3}, \ccc{Ray_3} or \ccc{Line_3}. + * The integer \ccc{dimension} is set to the dimension of the lowest + * dimensional face in the input complex containing the returned point, and + * \ccc{index} is set to the index to be stored at a mesh vertex lying + * on this face. + */ + struct Construct_intersection + { + Construct_intersection(const Labeled_mesh_domain_3& domain) + : r_domain_(domain) {} + + Intersection operator()(const Segment_3& s) const + { +#ifndef CGAL_MESH_3_NO_LONGER_CALLS_DO_INTERSECT_3 + CGAL_precondition(r_domain_.do_intersect_surface_object()(s)); +#endif // NOT CGAL_MESH_3_NO_LONGER_CALLS_DO_INTERSECT_3 + return this->operator()(s.source(),s.target()); + } + + Intersection operator()(const Ray_3& r) const + { + return clip_to_segment(r); + } + + Intersection operator()(const Line_3& l) const + { + return clip_to_segment(l); + } + + private: + /** + * Returns a point in the intersection of [a,b] with the surface + * \c a must be the source point, and \c b the out point. It's important + * because it drives bisection cuts. + * Indeed, the returned point is the first intersection from \c [a,b] + * with a subdomain surface. + */ + Intersection operator()(const Point_3& a, const Point_3& b) const + { + // Functors + typename BGT::Compute_squared_distance_3 squared_distance = + BGT().compute_squared_distance_3_object(); + typename BGT::Construct_midpoint_3 midpoint = + BGT().construct_midpoint_3_object(); + + // Non const points + Point_3 p1 = a; + Point_3 p2 = b; + Point_3 mid = midpoint(p1, p2); + + // Cannot be const: those values are modified below. + Subdomain_index value_at_p1 = r_domain_.function_(p1); + Subdomain_index value_at_p2 = r_domain_.function_(p2); + Subdomain_index value_at_mid = r_domain_.function_(mid,true); + + // If both extremities are in the same subdomain, + // there is no intersection. + // This should not happen... + if( value_at_p1 == value_at_p2 ) + { + return Intersection(); + } + + // Construct the surface patch index and index from the values at 'a' + // and 'b'. Even if the bissection find out a different pair of + // values, the reported index will be constructed from the initial + // values. + const Surface_patch_index sp_index = + r_domain_.make_surface_index(value_at_p1, value_at_p2); + const Index index = r_domain_.index_from_surface_patch_index(sp_index); + + // Else lets find a point (by bisection) + // Bisection ends when the point is near than error bound from surface + while(true) + { + // If the two points are enough close, then we return midpoint + if ( squared_distance(p1, p2) < r_domain_.squared_error_bound_ ) + { + CGAL_assertion(value_at_p1 != value_at_p2); + return Intersection(mid, index, 2); + } + + // Else we must go on + // Here we consider that p1(a) is the source point. Thus, we keep p1 and + // change p2 if f(p1)!=f(p2). + // That allows us to find the first intersection from a of [a,b] with + // a surface. + if ( value_at_p1 != value_at_mid ) + { + p2 = mid; + value_at_p2 = value_at_mid; + } + else + { + p1 = mid; + value_at_p1 = value_at_mid; + } + + mid = midpoint(p1, p2); + value_at_mid = r_domain_.function_(mid,true); + } + } + + /// Clips \c query to a segment \c s, and call operator()(s) + template + Intersection clip_to_segment(const Query& query) const + { + typename cpp11::result_of::type + clipped = CGAL::intersection(query, r_domain_.bbox_); + + if(clipped) +#if CGAL_INTERSECTION_VERSION > 1 + if(const Segment_3* s = boost::get(&*clipped)) + return this->operator()(*s); +#else + if(const Segment_3* s = object_cast(&clipped)) + return this->operator()(*s); +#endif + + return Intersection(); + } + + private: + const Labeled_mesh_domain_3& r_domain_; + }; + + /// Returns Construct_intersection object + Construct_intersection construct_intersection_object() const + { + return Construct_intersection(*this); + } + + /** + * Returns the index to be stored in a vertex lying on the surface identified + * by \c index. + */ + Index index_from_surface_patch_index(const Surface_patch_index& index) const + { return Index(index); } + + /** + * Returns the index to be stored in a vertex lying in the subdomain + * identified by \c index. + */ + Index index_from_subdomain_index(const Subdomain_index& index) const + { return Index(index); } + + /** + * Returns the \c Surface_patch_index of the surface patch + * where lies a vertex with dimension 2 and index \c index. + */ + Surface_patch_index surface_patch_index(const Index& index) const + { return boost::get(index); } + + /** + * Returns the index of the subdomain containing a vertex + * with dimension 3 and index \c index. + */ + Subdomain_index subdomain_index(const Index& index) const + { return boost::get(index); } + + // ----------------------------------- + // Backward Compatibility + // ----------------------------------- +#ifndef CGAL_MESH_3_NO_DEPRECATED_SURFACE_INDEX + typedef Surface_patch_index Surface_index; + + Index index_from_surface_index(const Surface_index& index) const + { return index_from_surface_patch_index(index); } + + Surface_index surface_index(const Index& index) const + { return surface_patch_index(index); } +#endif // CGAL_MESH_3_NO_DEPRECATED_SURFACE_INDEX + // ----------------------------------- + // End backward Compatibility + // ----------------------------------- + + +private: + /// Returns Surface_patch_index from \c i and \c j + Surface_patch_index make_surface_index(const Subdomain_index i, + const Subdomain_index j) const + { + if ( i < j ) return Surface_patch_index(i,j); + else return Surface_patch_index(j,i); + } + + /// Returns squared error bound from \c bbox and \c error + FT squared_error_bound(const Iso_cuboid_3& bbox, const FT& error) const + { + typename BGT::Compute_squared_distance_3 squared_distance = + BGT().compute_squared_distance_3_object(); + return squared_distance((bbox.min)(), (bbox.max)())*error*error/4; + } + + /// Returns squared error bound from \c sphere and \c error + FT squared_error_bound(const Sphere_3& sphere, const FT& error) const + { + typename BGT::Compute_squared_radius_3 squared_radius = + BGT().compute_squared_radius_3_object(); + return squared_radius(sphere)*error*error; + } + + /// Returns the bounding sphere of an Iso_cuboid_3 + Sphere_3 bounding_sphere(const Iso_cuboid_3& bbox) const + { + typename BGT::Construct_sphere_3 sphere = BGT().construct_sphere_3_object(); + return sphere((bbox.min)(), (bbox.max)()); + } + + /// Returns and Iso_cuboid_3 from a Bbox_3 + Iso_cuboid_3 iso_cuboid(const Bbox_3& bbox) + { + const Point_3 p_min(bbox.xmin(), bbox.ymin(), bbox.zmin()); + const Point_3 p_max(bbox.xmax(), bbox.ymax(), bbox.zmax()); + + return Iso_cuboid_3(p_min,p_max); + } + +<<<<<<< HEAD +public: + const Iso_cuboid_3& periodic_cuboid() const { return bounding_box(); } + +======= +>>>>>>> Labeled_mesh_domain_3 has a new constructor taking an Iso_cuboid. +protected: + /// Returns bounding box + const Iso_cuboid_3& bounding_box() const { return bbox_; } + +private: + /// The function which answers subdomain queries + const Function function_; + /// The bounding box + const Iso_cuboid_3 bbox_; + /// Error bound relative to sphere radius + FT squared_error_bound_; + + +private: + // Disabled copy constructor & assignment operator + typedef Labeled_mesh_domain_3 Self; + Labeled_mesh_domain_3(const Self& src); + Self& operator=(const Self& src); + +}; // end class Labeled_mesh_domain_3 + + + + +//------------------------------------------------------- +// Method implementation +//------------------------------------------------------- +template +Labeled_mesh_domain_3::Labeled_mesh_domain_3( + const F& f, + const Sphere_3& bounding_sphere, + const FT& error_bound ) +: function_(f) +, bbox_(iso_cuboid(bounding_sphere.bbox())) +, squared_error_bound_(squared_error_bound(bounding_sphere,error_bound)) +{ + // TODO : CGAL_ASSERT(0 < f(bounding_sphere.get_center()) ) ? +} + +template +Labeled_mesh_domain_3::Labeled_mesh_domain_3( + const F& f, + const Bbox_3& bbox, + const FT& error_bound ) +: function_(f) +, bbox_(iso_cuboid(bbox)) +, squared_error_bound_(squared_error_bound(bbox_,error_bound)) +{ + // TODO : CGAL_ASSERT(0 < f(bounding_sphere.get_center()) ) ? +} + + +template +template +OutputIterator +Labeled_mesh_domain_3::Construct_initial_points::operator()( + OutputIterator pts, + const int nb_points) const +{ + // Create point_iterator on and in bounding_sphere + typedef Random_points_on_sphere_3 Random_points_on_sphere_3; + typedef Random_points_in_sphere_3 Random_points_in_sphere_3; + + + const FT squared_radius = BGT().compute_squared_radius_3_object()( + r_domain_.bounding_sphere(r_domain_.bbox_)); + + const double radius = std::sqrt(CGAL::to_double(squared_radius)); + + Random_points_on_sphere_3 random_point_on_sphere(radius); + Random_points_in_sphere_3 random_point_in_sphere(radius); + + // Get some functors + typename BGT::Construct_segment_3 segment_3 = + BGT().construct_segment_3_object(); + typename BGT::Construct_vector_3 vector_3 = + BGT().construct_vector_3_object(); + typename BGT::Construct_translated_point_3 translate = + BGT().construct_translated_point_3_object(); + typename BGT::Construct_center_3 center = BGT().construct_center_3_object(); + + // Get translation from origin to sphere center + Point_3 center_pt = center(r_domain_.bounding_sphere(r_domain_.bbox_)); + const Vector_3 sphere_translation = vector_3(CGAL::ORIGIN, center_pt); + + // Create nb_point points + int n = nb_points; +#ifdef CGAL_MESH_3_VERBOSE + std::cerr << "construct initial points:\n"; +#endif + while ( 0 != n ) + { + // Get a random segment + const Point_3 random_point = translate(*random_point_on_sphere, + sphere_translation); + const Segment_3 random_seg = segment_3(center_pt, random_point); + + // Add the intersection to the output if it exists + Surface_patch surface = r_domain_.do_intersect_surface_object()(random_seg); + if ( surface ) + { + const Point_3 intersect_pt = CGAL::cpp11::get<0>( + r_domain_.construct_intersection_object()(random_seg)); + *pts++ = std::make_pair(intersect_pt, + r_domain_.index_from_surface_patch_index(*surface)); + --n; + +#ifdef CGAL_MESH_3_VERBOSE + std::cerr << boost::format("\r \r" + "%1%/%2% initial point(s) found...") + % (nb_points - n) + % nb_points; +#endif + } + else + { + // Get a new random point into sphere as center of object + // It may be necessary if the center of the domain is empty, e.g. torus + // In general case, it is good for input point dispersion + ++random_point_in_sphere; + center_pt = translate(*random_point_in_sphere, sphere_translation); + } + ++random_point_on_sphere; + } + +#ifdef CGAL_MESH_3_VERBOSE + std::cerr << "\n"; +#endif + return pts; +} + + +} // end namespace Mesh_3 + +} // end namespace CGAL + +#endif // LABELLED_MESH_TRAITS_3_H_