From dbc48009bc7e74d4df6ee9bd81b11094b8f83c5d Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Mael=20Rouxel-Labb=C3=A9?= Date: Fri, 12 Jan 2018 16:44:23 +0100 Subject: [PATCH] Various minor to natural/regular coordinates --- .../CGAL/Interpolation/internal/helpers.h | 4 +- .../include/CGAL/interpolation_functions.h | 33 +++--- .../CGAL/natural_neighbor_coordinates_2.h | 107 +++++++++++------- 3 files changed, 86 insertions(+), 58 deletions(-) diff --git a/Interpolation/include/CGAL/Interpolation/internal/helpers.h b/Interpolation/include/CGAL/Interpolation/internal/helpers.h index 1cee56c1947..79780b35520 100644 --- a/Interpolation/include/CGAL/Interpolation/internal/helpers.h +++ b/Interpolation/include/CGAL/Interpolation/internal/helpers.h @@ -23,6 +23,8 @@ #include +#include + #include namespace CGAL { @@ -40,7 +42,7 @@ struct Extract_bare_point typedef typename Traits::Point_d Point; typedef typename Traits::Weighted_point_d Weighted_point; - Extract_bare_point(const Traits& traits) : traits(traits) {} + Extract_bare_point(const Traits& traits = Traits()) : traits(traits) {} const Point& operator()(const Point& p) const { return p; } diff --git a/Interpolation/include/CGAL/interpolation_functions.h b/Interpolation/include/CGAL/interpolation_functions.h index 4d700b10df1..d7455f67ae4 100644 --- a/Interpolation/include/CGAL/interpolation_functions.h +++ b/Interpolation/include/CGAL/interpolation_functions.h @@ -73,14 +73,13 @@ linear_interpolation(ForwardIterator first, ForwardIterator beyond, { val = value_function(first->first); CGAL_assertion(val.second); - result += (first->second/norm) * val.first; + result += (first->second / norm) * val.first; } return result; } -template < class ForwardIterator, class ValueFunctor, class GradFunctor, - class Traits, class Point > +template < class ForwardIterator, class ValueFunctor, class GradFunctor, class Traits, class Point > std::pair< typename ValueFunctor::result_type::first_type, bool> quadratic_interpolation(ForwardIterator first, ForwardIterator beyond, const typename std::iterator_traits::value_type::second_type& norm, @@ -101,21 +100,23 @@ quadratic_interpolation(ForwardIterator first, ForwardIterator beyond, { f = value_function(first->first); grad = gradient_function(first->first); + //test if value and gradient are correctly retrieved: CGAL_assertion(f.second); if(!grad.second) return std::make_pair(Value_type(0), false); + result += (first->second/norm) * (f.first + grad.first * traits.construct_scaled_vector_d_object() - (traits.construct_vector_d_object()(cp(first->first), p),0.5)); + (traits.construct_vector_d_object()(cp(first->first), cp(p)), 0.5)); } + return std::make_pair(result, true); } -template < class ForwardIterator, class ValueFunctor, class GradFunctor, - class Traits, class Point > -std::pair< typename ValueFunctor::result_type::first_type, bool> +template < class ForwardIterator, class ValueFunctor, class GradFunctor, class Traits, class Point > +std::pair< typename ValueFunctor::result_type::first_type, bool > sibson_c1_interpolation(ForwardIterator first, ForwardIterator beyond, const typename std::iterator_traits::value_type::second_type& norm, const Point& p, @@ -127,23 +128,27 @@ sibson_c1_interpolation(ForwardIterator first, ForwardIterator beyond, typedef typename ValueFunctor::result_type::first_type Value_type; typedef typename Traits::FT Coord_type; + typedef typename Traits::Point_d Bare_point; Interpolation::internal::Extract_bare_point cp(traits); + const Bare_point& bp = cp(p); + Coord_type term1(0), term2(term1), term3(term1), term4(term1); Value_type linear_int(0), gradient_int(0); typename ValueFunctor::result_type f; typename GradFunctor::result_type grad; - for(; first !=beyond; ++first) + for(; first!=beyond; ++first) { f = value_function(first->first); grad = gradient_function(first->first); + CGAL_assertion(f.second); if(!grad.second) return std::make_pair(Value_type(0), false); //the values are not correct Coord_type coeff = first->second/norm; - Coord_type squared_dist = traits.compute_squared_distance_d_object()(cp(first->first), p); + Coord_type squared_dist = traits.compute_squared_distance_d_object()(cp(first->first), bp); Coord_type dist = CGAL_NTS sqrt(squared_dist); if(squared_dist == 0) @@ -163,7 +168,7 @@ sibson_c1_interpolation(ForwardIterator first, ForwardIterator beyond, linear_int += coeff * f.first; gradient_int += (coeff/dist) * (f.first + grad.first * - traits.construct_vector_d_object()(cp(first->first), p)); + traits.construct_vector_d_object()(cp(first->first), bp)); } term4 = term3 / term1; @@ -184,8 +189,7 @@ sibson_c1_interpolation(ForwardIterator first, ForwardIterator beyond, // term3 += coeff*(squared_dist/inv_weight); // gradient_int += (coeff/inv_weight) * (vh->get_value()+ vh->get_gradient() * (p - vh->point())); -template < class ForwardIterator, class ValueFunctor, class GradFunctor, - class Traits, class Point > +template < class ForwardIterator, class ValueFunctor, class GradFunctor, class Traits, class Point > std::pair< typename ValueFunctor::result_type::first_type, bool > sibson_c1_interpolation_square(ForwardIterator first, ForwardIterator beyond, const typename std::iterator_traits::value_type::second_type& norm, @@ -198,8 +202,11 @@ sibson_c1_interpolation_square(ForwardIterator first, ForwardIterator beyond, typedef typename ValueFunctor::result_type::first_type Value_type; typedef typename Traits::FT Coord_type; + typedef typename Traits::Point_d Bare_point; Interpolation::internal::Extract_bare_point cp(traits); + const Bare_point& bp = cp(p); + Coord_type term1(0), term2(term1), term3(term1), term4(term1); Value_type linear_int(0), gradient_int(0); typename ValueFunctor::result_type f; @@ -215,7 +222,7 @@ sibson_c1_interpolation_square(ForwardIterator first, ForwardIterator beyond, return std::make_pair(Value_type(0), false); // the gradient is not known Coord_type coeff = first->second/norm; - Coord_type squared_dist = traits.compute_squared_distance_d_object()(cp(first->first), cp(p)); + Coord_type squared_dist = traits.compute_squared_distance_d_object()(cp(first->first), bp); if(squared_dist ==0) { diff --git a/Interpolation/include/CGAL/natural_neighbor_coordinates_2.h b/Interpolation/include/CGAL/natural_neighbor_coordinates_2.h index d8f29af5d55..a0a566ff172 100644 --- a/Interpolation/include/CGAL/natural_neighbor_coordinates_2.h +++ b/Interpolation/include/CGAL/natural_neighbor_coordinates_2.h @@ -25,11 +25,15 @@ #include +#include +#include #include #include #include #include -#include + +#include +#include #include #include @@ -195,15 +199,63 @@ natural_neighbors_2(const Dt& dt, } +/////////////////////////////////////////////////////////////////////////////// +/////////////////////////////////////////////////////////////////////////////// + + +// This function is called when the conflict zone is known. +// OutputIterator has value type `pair< Dt::Geom_traits::Point_2, Dt::Geom_traits::FT>` +template < class Dt, class OutputIterator, class OutputFunctor, class EdgeIterator > +Triple< OutputIterator, typename Dt::Geom_traits::FT, bool > +natural_neighbor_coordinates_2(const Dt& dt, + const typename Dt::Geom_traits::Point_2& p, + OutputIterator out, + OutputFunctor fct, + EdgeIterator hole_begin, EdgeIterator hole_end) +{ + CGAL_precondition(dt.dimension() == 2); + + typedef Interpolation::internal:: + Project_vertex_output_iterator OutputIteratorWithFunctor; + typedef typename Dt::Geom_traits::FT FT; + + OutputIteratorWithFunctor op(out, fct); + Triple result = natural_neighbors_2(dt, p, op, hole_begin, hole_end); + + return make_triple(result.first.base(), result.second, result.third); +} + + +// Same as above but without OutputFunctor. Default to extracting the point, for backward compatibility. +template < class Dt, class OutputIterator, class EdgeIterator > +Triple< OutputIterator, typename Dt::Geom_traits::FT, bool > +natural_neighbor_coordinates_2(const Dt& dt, + const typename Dt::Geom_traits::Point_2& p, + OutputIterator out, + EdgeIterator hole_begin, EdgeIterator hole_end) +{ + typedef typename Dt::Geom_traits::FT FT; + typedef Interpolation::internal::Extract_point_in_pair OutputFunctor; + + return natural_neighbor_coordinates_2(dt, p, out, OutputFunctor(), hole_begin, hole_end); +} + + +/////////////////////////////////////////////////////////////////////////////// +/////////////////////////////////////////////////////////////////////////////// + + +// Function with a point. Can take a default starting face. template < class Dt, class OutputIterator, class OutputFunctor > Triple< OutputIterator, typename Dt::Geom_traits::FT, bool > natural_neighbor_coordinates_2(const Dt& dt, const typename Dt::Geom_traits::Point_2& p, OutputIterator out, OutputFunctor fct, - typename Dt::Face_handle start = - CGAL_TYPENAME_DEFAULT_ARG Dt::Face_handle() ) - + typename Dt::Face_handle start = CGAL_TYPENAME_DEFAULT_ARG Dt::Face_handle(), + typename boost::disable_if_c< + is_iterator::value + >::type* = 0) { CGAL_precondition(dt.dimension() == 2); @@ -233,48 +285,11 @@ natural_neighbor_coordinates_2(const Dt& dt, } -// This function is called when the conflict zone is known. -// OutputIterator has value type `pair< Dt::Geom_traits::Point_2, Dt::Geom_traits::FT>` -template < class Dt, class OutputIterator, class OutputFunctor, class EdgeIterator > -Triple< OutputIterator, typename Dt::Geom_traits::FT, bool > -natural_neighbor_coordinates_2(const Dt& dt, - const typename Dt::Geom_traits::Point_2& p, - OutputIterator out, OutputFunctor fct, - EdgeIterator hole_begin, EdgeIterator hole_end) -{ - CGAL_precondition(dt.dimension() == 2); - - typedef Interpolation::internal:: - Project_vertex_output_iterator OutputIteratorWithFunctor; - typedef typename Dt::Geom_traits::FT FT; - - OutputIteratorWithFunctor op(out, fct); - Triple result = natural_neighbors_2(dt, p, op, hole_begin, hole_end); - - return make_triple(result.first.base(), result.second, result.third); -} +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// -// Same as above but without OutputFunctor. Default to extracting the point, for backward compatibility. -template < class Dt, class OutputIterator, class EdgeIterator > -Triple< OutputIterator, typename Dt::Geom_traits::FT, bool > -natural_neighbor_coordinates_2(const Dt& dt, - const typename Dt::Geom_traits::Point_2& p, - OutputIterator out, - EdgeIterator hole_begin, EdgeIterator hole_end) -{ - /// @todo TEST ME - typedef typename Dt::Geom_traits::FT FT; - typedef Interpolation::internal::Extract_point_in_pair OutputFunctor; - - return natural_neighbor_coordinates_2(dt,p, out, OutputFunctor(), hole_begin, hole_end); -} - - -/**********************************************************/ -// Compute the coordinates for a vertex of the triangulation -// with respect to the other points in the triangulation. -// OutputIterator has value type `pair< Dt::Geom_traits::Point_2, Dt::Geom_traits::FT>` +// Function with a Vertex_handle template Triple< OutputIterator, typename Dt::Geom_traits::FT, bool > natural_neighbor_coordinates_2(const Dt& dt, @@ -328,6 +343,10 @@ natural_neighbor_coordinates_2(const Dt& dt, } +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// + + // OutputIterator has value type `std::pair< Dt::Geom_traits::Point_2, Dt::Geom_traits::FT>` template < class Dt, class OutputIterator, class OutputFunctor > class natural_neighbor_coordinates_2_object