Readability changes in Interpolation/test

This commit is contained in:
Mael Rouxel-Labbé 2018-01-11 17:20:26 +01:00
parent 91620fbec6
commit c7b2757405
5 changed files with 528 additions and 618 deletions

View File

@ -15,20 +15,10 @@
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0+
//
//
//
// Author(s) : Julia Floetotto
#include <iostream>
#include <cassert>
#include <utility>
#include <CGAL/double.h>
#include <CGAL/Random.h>
#include <CGAL/Origin.h>
#include <CGAL/point_generators_2.h>
#include <CGAL/algorithm.h>
#include <CGAL/natural_neighbor_coordinates_2.h>
#include <CGAL/Interpolation_traits_2.h>
@ -36,35 +26,40 @@
#include <CGAL/Interpolation_gradient_fitting_traits_2.h>
#include <CGAL/sibson_gradient_fitting.h>
#include <CGAL/double.h>
#include <CGAL/Random.h>
#include <CGAL/Origin.h>
#include <CGAL/point_generators_2.h>
#include <CGAL/algorithm.h>
#include <iostream>
#include <cassert>
#include <utility>
template < class ForwardIterator >
bool test_norm(ForwardIterator first, ForwardIterator beyond,
typename std::iterator_traits<ForwardIterator>::
value_type::second_type norm)
typename std::iterator_traits<ForwardIterator>::value_type::second_type norm)
{
typename std::iterator_traits<ForwardIterator>::value_type::second_type
sum(0);
for(; first !=beyond; first++)
sum+= first->second;
typename std::iterator_traits<ForwardIterator>::value_type::second_type sum(0);
for(; first !=beyond; first++)
sum += first->second;
return norm==sum;
return norm == sum;
}
template < class ForwardIterator >
bool test_barycenter(ForwardIterator first, ForwardIterator beyond,
typename std::iterator_traits<ForwardIterator>::
value_type::second_type norm,
const typename std::iterator_traits<ForwardIterator>::
value_type::first_type& p, const typename
std::iterator_traits<ForwardIterator>::
value_type::second_type& tolerance)
typename std::iterator_traits<ForwardIterator>::value_type::second_type norm,
const typename std::iterator_traits<ForwardIterator>::value_type::first_type& p,
const typename std::iterator_traits<ForwardIterator>::value_type::second_type& tolerance)
{
typedef typename
std::iterator_traits<ForwardIterator>::value_type::first_type Point;
typedef typename std::iterator_traits<ForwardIterator>::value_type::first_type Point;
Point b(CGAL::ORIGIN);
for(; first !=beyond; first++)
b = b+ (first->second/norm) * (first->first - CGAL::ORIGIN);
return CGAL::squared_distance(p,b) <= tolerance;
b = b + (first->second/norm) * (first->first - CGAL::ORIGIN);
return CGAL::squared_distance(p, b) <= tolerance;
}
/////////////////////////////////////////////////////////////////
@ -73,141 +68,112 @@ bool test_barycenter(ForwardIterator first, ForwardIterator beyond,
template < class ForwardIterator, class Functor, class GradFunctor,
class Gt>
bool
_test_sibson_c1_interpolation_sqrt(ForwardIterator ,
ForwardIterator ,
const typename
std::iterator_traits<ForwardIterator>
::value_type::second_type&
/* norm */,
const typename
std::iterator_traits<ForwardIterator>
::value_type::first_type& /* p */,
Functor /* f */,
_test_sibson_c1_interpolation_sqrt(ForwardIterator , ForwardIterator ,
const typename std::iterator_traits<ForwardIterator>::value_type::second_type& /* norm */,
const typename std::iterator_traits<ForwardIterator>::value_type::first_type& /* p */,
Functor /* f */,
GradFunctor /* grad_f */,
const Gt& /*geom_traits */,
const typename
std::iterator_traits<ForwardIterator>
::value_type::second_type& /* tolerance */,
const typename
std::iterator_traits<ForwardIterator>
::value_type::second_type& /* exact_value */,
const Gt& /*geom_traits */,
const typename std::iterator_traits<ForwardIterator>::value_type::second_type& /* tolerance */,
const typename std::iterator_traits<ForwardIterator>::value_type::second_type& /* exact_value */,
CGAL::Integral_domain_without_division_tag)
{
// bool UNTESTED_STUFF_BECAUSE_SQRT_IS_NOT_SUPPORTED;
// bool UNTESTED_STUFF_BECAUSE_SQRT_IS_NOT_SUPPORTED;
std::cout << std::endl
<< "NOTE: FT doesn't support sqrt(),"
" hence sibson_c1_interpolation is not tested." << std::endl;
" hence sibson_c1_interpolation is not tested." << std::endl;
return true;
}
template < class ForwardIterator, class Functor, class GradFunctor,
class Gt>
bool
_test_sibson_c1_interpolation_sqrt(ForwardIterator first,
ForwardIterator beyond,
const typename
std::iterator_traits<ForwardIterator>
::value_type::second_type&
norm, const typename
std::iterator_traits<ForwardIterator>
::value_type
::first_type& p,
Functor f, GradFunctor grad_f,
const Gt& geom_traits,
const typename
std::iterator_traits<ForwardIterator>
::value_type::second_type&
tolerance,
const typename
std::iterator_traits<ForwardIterator>
::value_type::second_type& exact_value,
CGAL::Field_with_sqrt_tag)
template < class ForwardIterator, class Functor, class GradFunctor, class Gt>
bool _test_sibson_c1_interpolation_sqrt(ForwardIterator first, ForwardIterator beyond,
const typename std::iterator_traits<ForwardIterator>::value_type::second_type& norm,
const typename std::iterator_traits<ForwardIterator>::value_type::first_type& p,
Functor f, GradFunctor grad_f,
const Gt& geom_traits,
const typename std::iterator_traits<ForwardIterator>::value_type::second_type& tolerance,
const typename std::iterator_traits<ForwardIterator>::value_type::second_type& exact_value,
CGAL::Field_with_sqrt_tag)
{
typename Functor::result_type res =
CGAL::sibson_c1_interpolation(first, beyond,
norm,p,f,
grad_f, geom_traits);
typename Functor::result_type res = CGAL::sibson_c1_interpolation(first, beyond,
norm, p, f,
grad_f, geom_traits);
return res.second && (CGAL_NTS abs(res.first-exact_value)<= tolerance);
}
template < class ForwardIterator, class Functor, class GradFunctor,
class Gt>
bool test_interpolation(ForwardIterator first, ForwardIterator beyond,
const typename
std::iterator_traits<ForwardIterator>
::value_type::second_type&
norm, const typename
std::iterator_traits<ForwardIterator>::value_type
::first_type& p,
Functor f, GradFunctor grad_f,
const Gt& geom_traits,
const int& i,
const typename std::iterator_traits<ForwardIterator>
::value_type::second_type& tolerance)
const typename std::iterator_traits<ForwardIterator>::value_type::second_type& norm,
const typename std::iterator_traits<ForwardIterator>::value_type::first_type& p,
Functor f, GradFunctor grad_f,
const Gt& geom_traits,
const int& i,
const typename std::iterator_traits<ForwardIterator>::value_type::second_type& tolerance)
{
typedef typename Functor::result_type::first_type Value_type;
typedef typename Functor::result_type::first_type Value_type;
assert(f(p).second);
Value_type exact_value = f(p).first;
if(i==0){
Value_type val = CGAL::linear_interpolation(first, beyond, norm,f);
if(i == 0)
{
Value_type val = CGAL::linear_interpolation(first, beyond, norm, f);
assert(CGAL_NTS abs(val-exact_value)<= tolerance);
}
typename Functor::result_type
res = CGAL::quadratic_interpolation(first, beyond, norm,p,f,
grad_f, geom_traits);
assert(res.second && (CGAL_NTS abs(res.first-exact_value)<=tolerance));
typename Functor::result_type res = CGAL::quadratic_interpolation(first, beyond, norm, p, f,
grad_f, geom_traits);
assert(res.second && (CGAL_NTS abs(res.first - exact_value) <= tolerance));
if(i<2){
if(i<2)
{
//without sqrt:
res = CGAL::sibson_c1_interpolation_square(first, beyond,
norm,p,f,
grad_f, geom_traits);
assert(res.second && (CGAL_NTS abs(res.first-exact_value)<= tolerance));
res = CGAL::sibson_c1_interpolation_square(first, beyond,
norm,p,f,
grad_f, geom_traits);
assert(res.second && (CGAL_NTS abs(res.first - exact_value) <= tolerance));
//with sqrt:
typedef CGAL::Algebraic_structure_traits<Value_type> AST;
assert(_test_sibson_c1_interpolation_sqrt(
first, beyond, norm,p,f, grad_f,
geom_traits, tolerance, exact_value,
typename AST::Algebraic_category()));
assert(_test_sibson_c1_interpolation_sqrt(first, beyond, norm, p, f, grad_f,
geom_traits, tolerance, exact_value,
typename AST::Algebraic_category()));
}
res = CGAL::farin_c1_interpolation(first, beyond, norm,p,f,grad_f,
geom_traits);
res = CGAL::farin_c1_interpolation(first, beyond, norm, p, f, grad_f, geom_traits);
assert(res.second);
assert(CGAL_NTS abs(res.first-exact_value)<= tolerance);
assert(CGAL_NTS abs(res.first - exact_value) <= tolerance);
return true;
}
template <class Triangul>
void
_test_interpolation_functions_2_delaunay( const Triangul &,
const typename Triangul::Geom_traits::FT& tolerance)
template <class Tr>
void _test_interpolation_functions_2_delaunay(const Tr &, const typename Tr::Geom_traits::FT& tolerance)
{
CGAL::Set_ieee_double_precision pfr;
Triangul T;
Tr T;
int n=20, m=20;
double r = 3;
double max_value = 5;
typedef typename Triangul::Geom_traits Gt;
typedef typename Tr::Geom_traits Gt;
typedef CGAL::Interpolation_traits_2<Gt> Traits;
typedef typename Triangul::Face_handle Face_handle;
typedef typename Tr::Face_handle Face_handle;
typedef typename Gt::Point_2 Point;
typedef typename Gt::FT Coord_type;
typedef typename Gt::Vector_2 Vector;
typedef std::map<Point, Coord_type, typename Gt::Less_xy_2> Point_value_map ;
typedef std::map<Point, Vector, typename Gt::Less_xy_2 > Point_vector_map;
typedef std::map<Point, Vector, typename Gt::Less_xy_2> Point_vector_map;
typedef std::vector< std::pair<Point, Coord_type> > Point_coordinate_vector;
typedef std::vector<std::pair<Point, Coord_type> > Point_coordinate_vector;
std::cout << "NN2: Testing random points." << std::endl;
//test random points in a square of length r:
std::vector<Point> points;
points.reserve(n+m);
@ -220,7 +186,7 @@ _test_interpolation_functions_2_delaunay( const Triangul &,
// Create n+m-4 points within a disc of radius 2
CGAL::Random_points_in_square_2<Point> g(r);
CGAL::cpp11::copy_n( g, n+m, std::back_inserter(points));
CGAL::cpp11::copy_n(g, n+m, std::back_inserter(points));
CGAL::Random random;
@ -228,167 +194,150 @@ _test_interpolation_functions_2_delaunay( const Triangul &,
Point_vector_map gradients[3];
Coord_type alpha = Coord_type(random.get_double(-max_value,max_value)),
beta1 = Coord_type(random.get_double(-max_value,max_value)),
beta2 = Coord_type(random.get_double(-max_value,max_value)),
gamma1 = Coord_type(random.get_double(-max_value,max_value)),
gamma2 = Coord_type(random.get_double(-max_value,max_value)),
gamma3 = Coord_type(random.get_double(-max_value,max_value));
beta1 = Coord_type(random.get_double(-max_value,max_value)),
beta2 = Coord_type(random.get_double(-max_value,max_value)),
gamma1 = Coord_type(random.get_double(-max_value,max_value)),
gamma2 = Coord_type(random.get_double(-max_value,max_value)),
gamma3 = Coord_type(random.get_double(-max_value,max_value));
//INSERTION + DET. of GRADIENT for n DATA POINTS :
for(int j=0; j<n; j++){
for(int j=0; j<n; j++)
{
T.insert(points[j]);
gradients[0].insert(std::make_pair(points[j], Vector( beta1, beta2)));
gradients[1].insert(std::make_pair
(points[j],
Vector( beta1+Coord_type(2)*gamma1*points[j].x(),
beta2+Coord_type(2)*gamma1* points[j].y())));
gradients[2].insert
(std::make_pair(points[j],Vector(beta1+Coord_type(2)*gamma1
*points[j].x()+gamma3* points[j].y(),
beta2 +
Coord_type(2) * gamma2* points[j].y()+
gamma3* points[j].x())));
gradients[1].insert(std::make_pair(points[j],
Vector(beta1 + Coord_type(2)*gamma1*points[j].x(),
beta2 + Coord_type(2)*gamma1*points[j].y())));
gradients[2].insert(std::make_pair(points[j],
Vector(beta1 + Coord_type(2)*gamma1*points[j].x() + gamma3*points[j].y(),
beta2 + Coord_type(2)*gamma2* points[j].y() + gamma3*points[j].x())));
}
//DETERMINE FUNCTION VALUE FOR n DATA POINTS AND m RANDOM TEST POINTS:
for(int j=0; j<n+m; j++){
for(int j=0; j<n+m; j++)
{
//linear function
values[0].insert(std::make_pair(points[j],alpha +
beta1*points[j].x()
+ beta2*points[j].y()));
values[0].insert(std::make_pair(points[j], alpha + beta1*points[j].x() + beta2*points[j].y()));
//spherical function:
values[1].insert(std::make_pair(points[j],alpha +
beta1*points[j].x() +
beta2*points[j].y() +
gamma1*points[j].x()*points[j].x()+
gamma1*points[j].y()*points[j].y()));
values[1].insert(std::make_pair(points[j], alpha + beta1*points[j].x() +
beta2*points[j].y() +
gamma1*points[j].x()*points[j].x()+
gamma1*points[j].y()*points[j].y()));
//quadratic function
values[2].insert(std::make_pair(points[j],alpha +
beta1*points[j].x() +
beta2*points[j].y() +
gamma1*points[j].x()*points[j].x()+
gamma2*points[j].y()*points[j].y()
+ gamma3*points[j].x()*points[j].y()));
values[2].insert(std::make_pair(points[j], alpha + beta1*points[j].x() +
beta2*points[j].y() +
gamma1*points[j].x()*points[j].x() +
gamma2*points[j].y()*points[j].y() +
gamma3*points[j].x()*points[j].y()));
}
//INTERPOLATION OF RANDOM POINTS:
Coord_type norm;
Point_coordinate_vector coords;
for(int j=n;j<n+m;j++){
CGAL::Triple<
std::back_insert_iterator<Point_coordinate_vector>,
Coord_type, bool> coordinate_result =
CGAL::natural_neighbor_coordinates_2(T, points[j],
std::back_inserter(coords));
for(int j=n; j<n+m; j++)
{
CGAL::Triple<std::back_insert_iterator<Point_coordinate_vector>, Coord_type, bool> coordinate_result =
CGAL::natural_neighbor_coordinates_2(T, points[j], std::back_inserter(coords));
assert(coordinate_result.third);
norm = coordinate_result.second;
assert(norm>0);
assert(test_norm( coords.begin(), coords.end(),norm));
assert(test_barycenter( coords.begin(), coords.end(),norm,
points[j], tolerance));
assert(norm > 0);
assert(test_norm(coords.begin(), coords.end(), norm));
assert(test_barycenter( coords.begin(), coords.end(), norm, points[j], tolerance));
for(int i=0; i<3; i++)
assert(test_interpolation(coords.begin(), coords.end(),norm,points[j],
CGAL::Data_access< Point_value_map >
(values[i]),
CGAL::Data_access< Point_vector_map >
(gradients[i]),
Traits(), i, tolerance ));
assert(test_interpolation(coords.begin(), coords.end(), norm, points[j],
CGAL::Data_access< Point_value_map >(values[i]),
CGAL::Data_access< Point_vector_map >(gradients[i]),
Traits(), i, tolerance));
coords.clear();
}
//TESTING THE GRADIENT APPRXIMATION METHOD:
std::cout << "Testing gradient estimation method on random points."
<<std::endl;
std::cout << "Testing gradient estimation method on random points." << std::endl;
typedef CGAL::Interpolation_gradient_fitting_traits_2<Gt> GradTraits;
Point_vector_map approx_gradients[2];
CGAL::sibson_gradient_fitting_nn_2(T,std::inserter(approx_gradients[0],
approx_gradients[0].
begin()),
CGAL::Data_access<Point_value_map>
(values[0]),
GradTraits());
CGAL::sibson_gradient_fitting_nn_2
(T,std::inserter(approx_gradients[1],approx_gradients[1].begin()),
CGAL::Data_access<Point_value_map>(values[1]),GradTraits());
for(int j=0; j<n; j++){
std::pair<Vector, bool> res =
CGAL::Data_access<Point_vector_map>(approx_gradients[0])(points[j]);
if(res.second){
CGAL::sibson_gradient_fitting_nn_2(T,
std::inserter(approx_gradients[0], approx_gradients[0].begin()),
CGAL::Data_access<Point_value_map>(values[0]), GradTraits());
CGAL::sibson_gradient_fitting_nn_2(T,
std::inserter(approx_gradients[1], approx_gradients[1].begin()),
CGAL::Data_access<Point_value_map>(values[1]),GradTraits());
for(int j=0; j<n; j++)
{
std::pair<Vector, bool> res = CGAL::Data_access<Point_vector_map>(approx_gradients[0])(points[j]);
if(res.second)
{
//if it is the exact computation kernel: test the equality:
assert(tolerance > Coord_type(0) || res.first ==
CGAL::Data_access<Point_vector_map>(gradients[0])
(points[j]).first);
res =
CGAL::Data_access<Point_vector_map>(approx_gradients[1])(points[j]);
assert(tolerance > Coord_type(0) ||
res.first == CGAL::Data_access<Point_vector_map>(gradients[0])(points[j]).first);
res = CGAL::Data_access<Point_vector_map>(approx_gradients[1])(points[j]);
//if one exists->the other must also exist
assert(res.second);
assert(tolerance > Coord_type(0) || res.first ==
CGAL::Data_access<Point_vector_map>(gradients[1])
(points[j]).first);
}else
assert(!CGAL::Data_access<Point_vector_map>(approx_gradients[1])
(points[j]).second);
assert(tolerance > Coord_type(0) ||
res.first == CGAL::Data_access<Point_vector_map>(gradients[1])(points[j]).first);
}
else
{
assert(!CGAL::Data_access<Point_vector_map>(approx_gradients[1])(points[j]).second);
}
}
//TESTING A POINT == A DATA POINT:
CGAL::Triple<
std::back_insert_iterator<Point_coordinate_vector>,
Coord_type, bool> coordinate_result =
CGAL::natural_neighbor_coordinates_2(T,points[n/2],std::back_inserter
(coords));
CGAL::Triple<std::back_insert_iterator<Point_coordinate_vector>, Coord_type, bool> coordinate_result =
CGAL::natural_neighbor_coordinates_2(T,points[n/2],std::back_inserter(coords));
assert(coordinate_result.third);
norm = coordinate_result.second;
assert(norm == Coord_type(1));
typename std::vector< std::pair< Point, Coord_type > >::iterator
ci= coords.begin();
typename std::vector< std::pair< Point, Coord_type > >::iterator ci= coords.begin();
assert(ci->first == points[n/2]);
assert(ci->second == Coord_type(1));
ci++;
assert(ci==coords.end());
for(int j=0; j<3; j++) {
assert(test_interpolation(coords.begin(), coords.end(),norm,points[n/2],
CGAL::Data_access< Point_value_map >(values[j]),
CGAL::Data_access< Point_vector_map >
(gradients[j]),
Traits(),j, tolerance));
for(int j=0; j<3; j++)
{
assert(test_interpolation(coords.begin(), coords.end(), norm, points[n/2],
CGAL::Data_access<Point_value_map>(values[j]),
CGAL::Data_access<Point_vector_map>(gradients[j]),
Traits(), j, tolerance));
}
coords.clear();
//FURTHER TESTS FOR NATURAL NEIGHBOR COORDINATES 2:
//TESTING A POINT on an EDGE of the triangulation:
Face_handle fh = T.finite_faces_begin();
int i =0;
int i = 0;
if(T.is_infinite(fh->neighbor(i)))
i++;
assert(!T.is_infinite(fh->neighbor(i)));
Point p = fh->vertex(T.ccw(i))->point() + Coord_type(0.5)*
(fh->vertex(T.cw(i))->point()-fh->vertex(T.ccw(i))->point());
coordinate_result =
CGAL::natural_neighbor_coordinates_2(T, p,
std::back_inserter(coords));
Point p = fh->vertex(T.ccw(i))->point() + Coord_type(0.5)*
(fh->vertex(T.cw(i))->point() - fh->vertex(T.ccw(i))->point());
coordinate_result = CGAL::natural_neighbor_coordinates_2(T, p, std::back_inserter(coords));
assert(coordinate_result.third);
norm = coordinate_result.second;
assert(test_norm( coords.begin(), coords.end(),norm));
assert(test_barycenter( coords.begin(), coords.end(),norm,p, tolerance));
assert(test_barycenter(coords.begin(), coords.end(), norm, p, tolerance));
coords.clear();
//END OF TEST WITH EDGE
//TESTING a GRID POINT SET
std::cout << "NN2: Testing grid points." << std::endl;
Triangul T2;
Tr T2;
//Grid points:
Point p1_2(-2, -2);
Point p2_2(-2,2);
@ -403,33 +352,34 @@ _test_interpolation_functions_2_delaunay( const Triangul &,
Point p34(1,0);
Point p41(0,-1);
T2.insert(p1_2);T2.insert(p2_2);T2.insert(p3_2);T2.insert(p4_2);
T2.insert(p1);T2.insert(p2);T2.insert(p3);T2.insert(p4);
T2.insert(p12);T2.insert(p23);T2.insert(p34);T2.insert(p41);
T2.insert(p1_2); T2.insert(p2_2); T2.insert(p3_2); T2.insert(p4_2);
T2.insert(p1); T2.insert(p2); T2.insert(p3); T2.insert(p4);
T2.insert(p12); T2.insert(p23); T2.insert(p34); T2.insert(p41);
coordinate_result =
CGAL::natural_neighbor_coordinates_2(T2,Point(0,0),
std::back_inserter(coords));
coordinate_result = CGAL::natural_neighbor_coordinates_2(T2,
Point(0,0),
std::back_inserter(coords));
assert(coordinate_result.third);
norm = coordinate_result.second;
assert(norm == Coord_type(1));
ci= coords.begin();
ci = coords.begin();
for(; ci!= coords.end(); ci++)
assert(ci->second == Coord_type(0.25));
assert(test_barycenter( coords.begin(), coords.end(),norm,
Point(0,0), tolerance));
assert(test_barycenter(coords.begin(), coords.end(), norm, Point(0,0), tolerance));
coords.clear();
//add the middle point of the grid
T2.insert(Point(0,0));
//point on a vertex;
coordinate_result =
CGAL::natural_neighbor_coordinates_2(T2,p34,
std::back_inserter(coords));
coordinate_result = CGAL::natural_neighbor_coordinates_2(T2, p34, std::back_inserter(coords));
assert(coordinate_result.third);
norm = coordinate_result.second;
assert(norm == Coord_type(1));
ci= coords.begin();
assert(ci->first == p34);
assert(ci->second == Coord_type(1));
@ -438,25 +388,20 @@ _test_interpolation_functions_2_delaunay( const Triangul &,
coords.clear();
//point on an edge:
p= Point(0,0.5);
coordinate_result =
CGAL::natural_neighbor_coordinates_2(T2,p,
std::back_inserter(coords));
p = Point(0,0.5);
coordinate_result = CGAL::natural_neighbor_coordinates_2(T2, p, std::back_inserter(coords));
assert(coordinate_result.third);
norm = coordinate_result.second;
assert(test_barycenter( coords.begin(), coords.end(),norm, p, tolerance));
assert(test_barycenter(coords.begin(), coords.end(), norm, p, tolerance));
coords.clear();
//Point outside convex hull:
p= Point(3,0.5);
coordinate_result =
CGAL::natural_neighbor_coordinates_2(T2,p,
std::back_inserter(coords));
p = Point(3,0.5);
coordinate_result = CGAL::natural_neighbor_coordinates_2(T2, p, std::back_inserter(coords));
assert(!coordinate_result.third);
//Point on a convex hull edge:
p= Point(2,1);
coordinate_result =
CGAL::natural_neighbor_coordinates_2(T2,p,
std::back_inserter(coords));
p = Point(2,1);
coordinate_result = CGAL::natural_neighbor_coordinates_2(T2, p, std::back_inserter(coords));
assert(coordinate_result.third);
}

View File

@ -15,71 +15,70 @@
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0+
//
//
//
// Author(s) : Naceur MESKINI.
#include <iostream>
#include <cassert>
#include <utility>
#include <CGAL/natural_neighbor_coordinates_2.h>
#include <CGAL/algorithm.h>
#include <CGAL/Random.h>
#include <CGAL/barycenter.h>
#include <CGAL/natural_neighbor_coordinates_2.h>
#include <iostream>
#include <cassert>
#include <utility>
template < class ForwardIterator >
bool test_norm(ForwardIterator first, ForwardIterator beyond,
typename std::iterator_traits<ForwardIterator>
::value_type::second_type norm)
typename std::iterator_traits<ForwardIterator>::value_type::second_type norm)
{
typename
std::iterator_traits<ForwardIterator>::value_type::second_type sum(0);
for(; first !=beyond; first++)
sum+= first->second;
typename
std::iterator_traits<ForwardIterator>::value_type::second_type sum(0);
for(; first !=beyond; first++)
sum += first->second;
return norm==sum;
return norm == sum;
}
template < class ForwardIterator, class Point >
bool test_barycenter(ForwardIterator first, ForwardIterator beyond,
typename std::iterator_traits<ForwardIterator>
::value_type::second_type /*norm*/, const Point& p)
typename std::iterator_traits<ForwardIterator>::value_type::second_type /*norm*/,
const Point& p)
{
return p == CGAL::barycenter(first, beyond);
}
template <class Triangul>
void
_test_natural_neighbors_2( const Triangul & )
template <class Tr>
void _test_natural_neighbors_2( const Tr & )
{
typedef typename Triangul::Geom_traits Gt;
typedef typename Tr::Geom_traits Gt;
typedef typename Gt::Point_2 Point_2;
typedef typename Gt::FT Coord_type;
typedef std::vector< std::pair<Point_2, Coord_type> > Point_coordinate_vector;
//TESTING a GRID POINT SET
// TESTING a GRID POINT SET
std::cout << "NN2: Testing grid points." << std::endl;
Triangul T;
Tr T;
//Grid points:
Point_2 p1_2(-2, -2);
Point_2 p2_2(-2,2);
Point_2 p3_2(2,-2);
Point_2 p4_2(2,2);
Point_2 p2_2(-2, 2);
Point_2 p3_2(2, -2);
Point_2 p4_2(2, 2);
Point_2 p1(-1, -1);
Point_2 p2(-1,1);
Point_2 p3(1,-1);
Point_2 p4(1,1);
Point_2 p2(-1, 1);
Point_2 p3(1, -1);
Point_2 p4(1, 1);
Point_2 p12(-1, 0);
Point_2 p23(0,1);
Point_2 p34(1,0);
Point_2 p41(0,-1);
Point_2 p23(0, 1);
Point_2 p34(1, 0);
Point_2 p41(0, -1);
T.insert(p1_2);
T.insert(p2_2);
@ -95,45 +94,44 @@ _test_natural_neighbors_2( const Triangul & )
T.insert(p41);
T.insert(Point_2(0,0));
Point_coordinate_vector coords;
//point on a vertex;
// point on a vertex;
Point_2 p = Point_2(p34);
CGAL::Triple<std::back_insert_iterator<Point_coordinate_vector>,Coord_type, bool> coordinate_result = CGAL::natural_neighbor_coordinates_2(T,p,
std::back_inserter(coords));
CGAL::Triple<std::back_insert_iterator<Point_coordinate_vector>, Coord_type, bool> coordinate_result =
CGAL::natural_neighbor_coordinates_2(T, p, std::back_inserter(coords));
assert(coordinate_result.third);
Coord_type norm = coordinate_result.second;
assert(norm == Coord_type(1));
typename std::vector< std::pair< Point_2, Coord_type >
>::const_iterator ci= coords.begin();
typename std::vector<std::pair<Point_2, Coord_type> >::const_iterator ci= coords.begin();
assert(ci->first == p);
assert(ci->second == Coord_type(1));
ci++;
assert(ci==coords.end());
coords.clear();
//point on an edge:
// point on an edge:
p = Point_2(0,0.5);
coordinate_result = CGAL::natural_neighbor_coordinates_2
(T,p,std::back_inserter(coords));
coordinate_result = CGAL::natural_neighbor_coordinates_2(T, p, std::back_inserter(coords));
assert(coordinate_result.third);
norm = coordinate_result.second;
assert(test_barycenter( coords.begin(), coords.end(),norm, p));
assert(test_barycenter(coords.begin(), coords.end(), norm, p));
coords.clear();
//outside convex hull:
p= Point_2(3,0.5);
coordinate_result = CGAL::natural_neighbor_coordinates_2
(T,p,std::back_inserter(coords));
// outside convex hull:
p = Point_2(3,0.5);
coordinate_result = CGAL::natural_neighbor_coordinates_2(T, p, std::back_inserter(coords));
assert(!coordinate_result.third);
//on a convex hull edge:
// on a convex hull edge:
coords.clear();
p= Point_2(2,1);
coordinate_result = CGAL::natural_neighbor_coordinates_2
(T,p,std::back_inserter(coords));
p = Point_2(2, 1);
coordinate_result = CGAL::natural_neighbor_coordinates_2(T, p, std::back_inserter(coords));
assert(coordinate_result.third);
norm = coordinate_result.second;
assert(test_barycenter( coords.begin(), coords.end(),norm, p));

View File

@ -15,47 +15,44 @@
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0+
//
//
// Author(s) : Julia Floetotto
#include <iostream>
#include <cassert>
#include <utility>
#include <CGAL/regular_neighbor_coordinates_2.h>
#include <CGAL/Interpolation_gradient_fitting_traits_2.h>
#include <CGAL/sibson_gradient_fitting.h>
#include <CGAL/point_generators_2.h>
#include <CGAL/algorithm.h>
#include <CGAL/Random.h>
#include <CGAL/Origin.h>
#include <CGAL/regular_neighbor_coordinates_2.h>
#include <CGAL/Interpolation_gradient_fitting_traits_2.h>
#include <CGAL/sibson_gradient_fitting.h>
#include <iostream>
#include <cassert>
#include <utility>
template < class ForwardIterator >
bool test_norm(ForwardIterator first, ForwardIterator beyond,
typename std::iterator_traits<ForwardIterator>
::value_type::second_type norm)
typename std::iterator_traits<ForwardIterator>::value_type::second_type norm)
{
typename
std::iterator_traits<ForwardIterator>::value_type::second_type sum(0);
typename std::iterator_traits<ForwardIterator>::value_type::second_type sum(0);
for(; first !=beyond; first++)
sum+= first->second;
sum += first->second;
return norm==sum;
return norm == sum;
}
template < class ForwardIterator, class Point >
bool test_barycenter(ForwardIterator first, ForwardIterator beyond,
typename std::iterator_traits<ForwardIterator>
::value_type::second_type norm, const Point& p)
typename std::iterator_traits<ForwardIterator>::value_type::second_type norm,
const Point& p)
{
typedef typename CGAL::Kernel_traits<Point>::Kernel::Point_2 Bare_point;
Bare_point b = CGAL::ORIGIN;
for(; first != beyond; ++first)
b = b + (first->second/norm) * (first->first.point() - CGAL::ORIGIN);
return p==b;
return p == b;
}
template <class Point, class FT>
@ -67,17 +64,17 @@ struct Functor
result_type operator()(const Point&) { return std::make_pair(0., true); }
};
template <class Triangul>
void test_gradient_fitting(const Triangul& rt)
template <class Tr>
void test_gradient_fitting(const Tr& rt)
{
typedef typename Triangul::Geom_traits::Kernel K;
typedef typename Triangul::Geom_traits::FT FT;
typedef typename Triangul::Geom_traits::Vector_2 Vector_2;
typedef typename Triangul::Geom_traits::Weighted_point_2 Point;
typedef typename Tr::Geom_traits::Kernel K;
typedef typename Tr::Geom_traits::FT FT;
typedef typename Tr::Geom_traits::Vector_2 Vector_2;
typedef typename Tr::Geom_traits::Weighted_point_2 Point;
typedef typename std::back_insert_iterator<
std::vector<
std::pair< Point, Vector_2 > > > OutputIterator;
std::pair<Point, Vector_2> > > OutputIterator;
std::vector<std::pair< Point, Vector_2 > > v;
OutputIterator out = std::back_inserter(v);
@ -88,23 +85,21 @@ void test_gradient_fitting(const Triangul& rt)
sibson_gradient_fitting_rn_2(rt, out, f, traits);
}
template <class Triangul>
void
_test_regular_neighbors_2( const Triangul & )
template <class Tr>
void _test_regular_neighbors_2(const Tr &)
{
Triangul T;
Tr T;
int n=20, m=200;
double r = 3;
double max_weight =1;
typedef typename Triangul::Geom_traits Gt;
typedef typename Gt::Weighted_point_2 Weighted_point;
typedef typename Gt::Point_2 Bare_point;
typedef typename Gt::FT Coord_type;
typedef typename Tr::Geom_traits Gt;
typedef typename Gt::Weighted_point_2 Weighted_point;
typedef typename Gt::Point_2 Bare_point;
typedef typename Gt::FT Coord_type;
typedef std::vector< std::pair< Weighted_point, Coord_type > >
Point_coordinate_vector;
typedef std::vector<std::pair<Weighted_point,Coord_type > > Point_coordinate_vector;
std::cout << "RN2: Testing random points." << std::endl;
@ -120,28 +115,23 @@ _test_regular_neighbors_2( const Triangul & )
// Create n+m-4 points within a disc of radius 2
CGAL::Random_points_in_square_2<Bare_point> g(r);
CGAL::cpp11::copy_n( g, n+m, std::back_inserter(points));
CGAL::cpp11::copy_n(g, n+m, std::back_inserter(points));
CGAL::Random random;
for(int i=0; i<n ; i++)
T.insert(Weighted_point(points[i],random.get_double(-max_weight,
max_weight)));
T.insert(Weighted_point(points[i],random.get_double(-max_weight, max_weight)));
Point_coordinate_vector coords;
for(int i=n;i<n+m;i++){
Weighted_point wp = Weighted_point(points[i],random.get_double
(2*max_weight,3*max_weight));
CGAL::Triple<
std::back_insert_iterator<Point_coordinate_vector>,
Coord_type, bool> coordinate_result =
for(int i=n; i<n+m; i++)
{
Weighted_point wp = Weighted_point(points[i], random.get_double(2*max_weight, 3*max_weight));
CGAL::Triple<std::back_insert_iterator<Point_coordinate_vector>, Coord_type, bool> coordinate_result =
CGAL::regular_neighbor_coordinates_2(T,wp,std::back_inserter(coords));
assert(coordinate_result.third);
Coord_type norm = coordinate_result.second;
assert(norm>0);
assert(test_barycenter( coords.begin(), coords.end(),norm, points[i]));
assert(norm > 0);
assert(test_barycenter(coords.begin(), coords.end(), norm, points[i]));
coords.clear();
}
@ -150,7 +140,7 @@ _test_regular_neighbors_2( const Triangul & )
//TESTING a GRID POINT SET
std::cout << "RN2: Testing grid points." << std::endl;
Triangul T2;
Tr T2;
//Grid points:
Bare_point p1_2(-2, -2);
Bare_point p2_2(-2,2);
@ -166,101 +156,88 @@ _test_regular_neighbors_2( const Triangul & )
Bare_point p41(0,-1);
T2.insert(Weighted_point(p1_2, 0));
T2.insert(Weighted_point(p2_2, 0));
T2.insert(Weighted_point(p2_2, 0));
T2.insert(Weighted_point(p3_2, 0));
T2.insert(Weighted_point(p4_2, 0));
T2.insert(Weighted_point(p1, 0));
T2.insert(Weighted_point(p2, 0));
T2.insert(Weighted_point(p2, 0));
T2.insert(Weighted_point(p3, 0));
T2.insert(Weighted_point(p4, 0));
T2.insert(Weighted_point(p12, 0));
T2.insert(Weighted_point(p23, 0));
T2.insert(Weighted_point(p34, 0));
T2.insert(Weighted_point(p41, 0));
T2.insert(Weighted_point(p12, 0));
T2.insert(Weighted_point(p23, 0));
T2.insert(Weighted_point(p34, 0));
T2.insert(Weighted_point(p41, 0));
//test with 0 weight:
Weighted_point wp(Bare_point(0,0),0);
CGAL::Triple<
std::back_insert_iterator<Point_coordinate_vector>,
Coord_type, bool> coordinate_result =
CGAL::regular_neighbor_coordinates_2(T2,wp,
std::back_inserter(coords));
Weighted_point wp(Bare_point(0,0), 0);
CGAL::Triple<std::back_insert_iterator<Point_coordinate_vector>, Coord_type, bool> coordinate_result =
CGAL::regular_neighbor_coordinates_2(T2, wp, std::back_inserter(coords));
assert(coordinate_result.third);
Coord_type norm = coordinate_result.second;
assert(norm == Coord_type(1));
typename std::vector< std::pair< Weighted_point, Coord_type >
>::const_iterator
ci= coords.begin();
typename std::vector< std::pair<Weighted_point, Coord_type> >::const_iterator ci = coords.begin();
for(; ci!= coords.end(); ci++)
assert(ci->second == Coord_type(0.25));
assert(test_barycenter( coords.begin(), coords.end(),norm, wp));
assert(test_barycenter(coords.begin(), coords.end(), norm, wp));
coords.clear();
//test with hidden_vertices:
wp = Weighted_point(Bare_point(0,0),4);
coordinate_result =
CGAL::regular_neighbor_coordinates_2(T2,wp,
std::back_inserter(coords));
wp = Weighted_point(Bare_point(0,0), 4);
coordinate_result = CGAL::regular_neighbor_coordinates_2(T2, wp, std::back_inserter(coords));
assert(coordinate_result.third);
norm = coordinate_result.second;
assert(test_barycenter( coords.begin(), coords.end(),norm,wp));
assert(test_barycenter(coords.begin(), coords.end(), norm, wp));
coords.clear();
//add the middle point of the grid
T2.insert(Weighted_point(Bare_point(0,0),0));
T2.insert(Weighted_point(Bare_point(0,0), 0));
//point on a vertex;
wp = Weighted_point(p34,0);
coordinate_result =
CGAL::regular_neighbor_coordinates_2(T2,wp,
std::back_inserter(coords));
wp = Weighted_point(p34, 0);
coordinate_result = CGAL::regular_neighbor_coordinates_2(T2, wp, std::back_inserter(coords));
assert(coordinate_result.third);
norm = coordinate_result.second;
assert(norm == Coord_type(1));
ci= coords.begin();
ci = coords.begin();
assert(ci->first == wp);
assert(ci->second == Coord_type(1));
ci++;
assert(ci==coords.end());
assert(ci == coords.end());
coords.clear();
//point on the vertex but creating a hole:
wp = Weighted_point(p34,2);
coordinate_result =
CGAL::regular_neighbor_coordinates_2(T2,wp,
std::back_inserter(coords));
wp = Weighted_point(p34, 2);
coordinate_result = CGAL::regular_neighbor_coordinates_2(T2, wp, std::back_inserter(coords));
assert(coordinate_result.third);
norm = coordinate_result.second;
assert(test_barycenter( coords.begin(), coords.end(),norm, wp));
assert(test_barycenter(coords.begin(), coords.end(), norm, wp));
coords.clear();
//point on an edge:
wp= Weighted_point(Bare_point(0,0.5), 3);
coordinate_result = CGAL::regular_neighbor_coordinates_2
(T2,wp,std::back_inserter(coords));
wp = Weighted_point(Bare_point(0,0.5), 3);
coordinate_result = CGAL::regular_neighbor_coordinates_2(T2, wp, std::back_inserter(coords));
assert(coordinate_result.third);
norm = coordinate_result.second;
assert(test_barycenter( coords.begin(), coords.end(),norm, wp));
assert(test_barycenter(coords.begin(), coords.end(), norm, wp));
coords.clear();
//a vertex v in Reg(P\v->point()):
typename Triangul::Vertex_iterator vit = T2.finite_vertices_end();
coordinate_result = CGAL::regular_neighbor_coordinates_2
(T2, --vit, std::back_inserter(coords));
typename Tr::Vertex_iterator vit = T2.finite_vertices_end();
coordinate_result = CGAL::regular_neighbor_coordinates_2(T2, --vit, std::back_inserter(coords));
assert(coordinate_result.third);
norm = coordinate_result.second;
assert(test_barycenter( coords.begin(), coords.end(),norm,vit->point()));
assert(test_barycenter(coords.begin(), coords.end(), norm,vit->point()));
coords.clear();
//outside convex hull:
wp= Weighted_point(Bare_point(3,0.5), 3);
coordinate_result = CGAL::regular_neighbor_coordinates_2
(T2,wp,std::back_inserter(coords));
wp = Weighted_point(Bare_point(3,0.5), 3);
coordinate_result = CGAL::regular_neighbor_coordinates_2(T2, wp, std::back_inserter(coords));
assert(!coordinate_result.third);
//on a convex hull edge:
wp= Weighted_point(Bare_point(2,1), 3);
coordinate_result = CGAL::regular_neighbor_coordinates_2
(T2,wp,std::back_inserter(coords));
wp = Weighted_point(Bare_point(2,1), 3);
coordinate_result = CGAL::regular_neighbor_coordinates_2(T2, wp, std::back_inserter(coords));
assert(!coordinate_result.third);
}

View File

@ -15,15 +15,10 @@
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0+
//
//
//
// Author(s) : Julia Floetotto
#include <iostream>
#include <cassert>
#include <fstream>
#include <utility>
#include <CGAL/Voronoi_intersection_2_traits_3.h>
#include <CGAL/surface_neighbor_coordinates_3.h>
#include <CGAL/surface_neighbors_3.h>
@ -32,52 +27,50 @@
#include <CGAL/point_generators_2.h>
#include <CGAL/algorithm.h>
#include <CGAL/function_objects.h>
#include <CGAL/Origin.h>
#include <iostream>
#include <cassert>
#include <fstream>
#include <utility>
template < class ForwardIterator >
bool test_norm(ForwardIterator first, ForwardIterator beyond,
typename std::iterator_traits<ForwardIterator>
::value_type::second_type norm)
typename std::iterator_traits<ForwardIterator>::value_type::second_type norm)
{
typename
std::iterator_traits<ForwardIterator>::value_type::second_type sum(0);
for(; first !=beyond; first++)
sum+= first->second;
typename
std::iterator_traits<ForwardIterator>::value_type::second_type sum(0);
for(; first !=beyond; first++)
sum += first->second;
return norm==sum;
return norm == sum;
}
template < class ForwardIterator >
bool test_barycenter(ForwardIterator first, ForwardIterator beyond,
typename std::iterator_traits<ForwardIterator>
::value_type::second_type norm ,
const typename
std::iterator_traits<ForwardIterator>::value_type
::first_type& p,
const typename
std::iterator_traits<ForwardIterator>::value_type
::second_type& tolerance )
typename std::iterator_traits<ForwardIterator>::value_type::second_type norm ,
const typename std::iterator_traits<ForwardIterator>::value_type::first_type& p,
const typename std::iterator_traits<ForwardIterator>::value_type::second_type& tolerance)
{
typedef typename
std::iterator_traits<ForwardIterator>::value_type::first_type Point;
std::iterator_traits<ForwardIterator>::value_type::first_type Point;
Point b(0,0,0);
ForwardIterator it=first;
for(; it !=beyond; it++)
b = b + (it->second/norm)* (it->first - CGAL::ORIGIN);
for(; it!=beyond; it++)
b = b + (it->second/norm) * (it->first - CGAL::ORIGIN);
return CGAL::squared_distance(p,b) <= tolerance;
}
template < class ForwardIteratorCoord, class ForwardIteratorPoint,
class Kernel >
template < class ForwardIteratorCoord, class ForwardIteratorPoint, class Kernel >
bool compare_neighbors(ForwardIteratorCoord first_coord,
ForwardIteratorCoord beyond_coord,
ForwardIteratorPoint first_point,
ForwardIteratorPoint beyond_point,
Kernel )
ForwardIteratorCoord beyond_coord,
ForwardIteratorPoint first_point,
ForwardIteratorPoint beyond_point,
Kernel )
{
typedef std::set<typename std::iterator_traits<ForwardIteratorPoint>::value_type,
typename Kernel::Less_xyz_3 > Point_set;
typename Kernel::Less_xyz_3 > Point_set;
Point_set neighbors, coord_neighbors, diff;
@ -87,16 +80,16 @@ bool compare_neighbors(ForwardIteratorCoord first_coord,
coord_neighbors.insert(first_coord->first);
std::set_difference(neighbors.begin(),
neighbors.end(),
coord_neighbors.begin(),
coord_neighbors.end(),
std::inserter(diff, diff.begin()),
typename Kernel::Less_xyz_3());
neighbors.end(),
coord_neighbors.begin(),
coord_neighbors.end(),
std::inserter(diff, diff.begin()),
typename Kernel::Less_xyz_3());
if(!diff.empty()){
if(!diff.empty())
{
std::cout << "Compare neighbors -- Diff: " << std::endl;
for(typename Point_set::const_iterator it = diff.begin();
it != diff.end(); ++it)
for(typename Point_set::const_iterator it = diff.begin(); it != diff.end(); ++it)
std::cout << " point: " << *it;
std::cout << std::endl;
}
@ -104,214 +97,206 @@ bool compare_neighbors(ForwardIteratorCoord first_coord,
return diff.empty();
}
template < class Triangul, class OutputIterator>
template < class Tr, class OutputIterator>
OutputIterator
test_neighbors(const Triangul& T, const typename
Triangul::Geom_traits::Point_3& p,
const typename Triangul::Geom_traits::Vector_3& n,
test_neighbors(const Tr& T, const typename Tr::Geom_traits::Point_3& p,
const typename Tr::Geom_traits::Vector_3& n,
const int& version,
OutputIterator out)
{
typedef CGAL::Voronoi_intersection_2_traits_3<typename
Triangul::Geom_traits> I_traits;
typedef CGAL::Voronoi_intersection_2_traits_3<typename Tr::Geom_traits> I_traits;
//the result type of the certified version:
typedef std::pair< OutputIterator, bool > NeighborIt_bool_pair;
typename Triangul::Cell_handle start;
typename Tr::Cell_handle start;
//test different function calls
switch(version){
case 0:{
//certified call with Kernel:
NeighborIt_bool_pair
result_pair = CGAL::surface_neighbors_certified_3(T.points_begin(),
T.points_end(),p,n,
out,
T.geom_traits());
assert(result_pair.second);
out =result_pair.first; break;}
case 1: {
//certified call with instantiated traits::
NeighborIt_bool_pair
result_pair = CGAL::surface_neighbors_certified_3(T.points_begin(),
T.points_end(),p,
out,I_traits(p,n));
assert(result_pair.second);
out =result_pair.first; break;}
//both versions with locate:
case 2:{
start = T.locate(p);
//certified call with Kernel and locate:
switch(version)
{
case 0:{
//certified call with Kernel:
NeighborIt_bool_pair
result_pair = CGAL::surface_neighbors_certified_3(T.points_begin(),
T.points_end(), p, n,
out,
T.geom_traits());
assert(result_pair.second);
out = result_pair.first; break;}
case 1: {
//certified call with instantiated traits::
NeighborIt_bool_pair
result_pair = CGAL::surface_neighbors_certified_3(T.points_begin(),
T.points_end(), p,
out, I_traits(p,n));
assert(result_pair.second);
out =result_pair.first; break;}
//both versions with locate:
case 2:{
start = T.locate(p);
//certified call with Kernel and locate:
out =CGAL::surface_neighbors_3(T, p,n,out, start);
break;}
case 3: {
start = T.locate(p);
//with instantiated traits and locate:
break;}
case 3: {
start = T.locate(p);
//with instantiated traits and locate:
out =CGAL::surface_neighbors_3(T,p,out, I_traits(p,n),start);
break;}
//taking all points:
case 4: {
//with instantiated traits and locate:
out =
CGAL::surface_neighbors_3(T,p,out,I_traits(p,n));
break;}
case 5: {
//certified call with Kernel and locate:
out =
CGAL::surface_neighbors_3(T, p,n,out);
break;}
//the last two with certification:
case 6: {
out =
CGAL::surface_neighbors_3(T.points_begin(),
T.points_end(),
p, out,I_traits(p,n));
break;
}
case 7: {
out =
CGAL::surface_neighbors_3(T.points_begin(),
T.points_end(),
p,n,out,T.geom_traits());
break;
}
default:
std::cout << "Switch function calls: Nothing is tested. " <<
std::endl;
break;}
//taking all points:
case 4: {
//with instantiated traits and locate:
out =
CGAL::surface_neighbors_3(T,p,out,I_traits(p,n));
break;}
case 5: {
//certified call with Kernel and locate:
out =
CGAL::surface_neighbors_3(T, p,n,out);
break;}
//the last two with certification:
case 6: {
out =
CGAL::surface_neighbors_3(T.points_begin(),
T.points_end(),
p, out,I_traits(p,n));
break;
}
case 7: {
out =
CGAL::surface_neighbors_3(T.points_begin(),
T.points_end(),
p,n,out,T.geom_traits());
break;
}
default:
std::cout << "Switch function calls: Nothing is tested. " <<
std::endl;
}
return out;
}
template < class Triangul, class OutputIterator>
std::pair< OutputIterator, typename Triangul::Geom_traits::FT>
test_coords(const Triangul& T, const typename
Triangul::Geom_traits::Point_3& p,
const typename Triangul::Geom_traits::Vector_3& n,
const int& version, OutputIterator out)
template < class Tr, class OutputIterator>
std::pair< OutputIterator, typename Tr::Geom_traits::FT>
test_coords(const Tr& T,
const typename Tr::Geom_traits::Point_3& p,
const typename Tr::Geom_traits::Vector_3& n,
const int& version, OutputIterator out)
{
typedef CGAL::Voronoi_intersection_2_traits_3<typename
Triangul::Geom_traits> I_traits;
typedef CGAL::Voronoi_intersection_2_traits_3<typename Tr::Geom_traits> I_traits;
//coordinate computation result types
typedef CGAL::Triple< OutputIterator,
typename Triangul::Geom_traits::FT, bool > Result_triple;
typedef CGAL::Triple< OutputIterator, typename Tr::Geom_traits::FT, bool > Result_triple;
//the result type of the certified version:
typedef CGAL::Quadruple< OutputIterator,
typename Triangul::Geom_traits::FT, bool, bool > Result_quadruple;
typedef CGAL::Quadruple< OutputIterator, typename Tr::Geom_traits::FT, bool, bool > Result_quadruple;
typename Triangul::Cell_handle start;
typename Triangul::Geom_traits::FT norm = 1; // 1 for that default doesn't trigger an assert
typename Tr::Cell_handle start;
typename Tr::Geom_traits::FT norm = 1; // 1 for that default doesn't trigger an assert
//test different function calls
switch(version){
case 0:{
Result_triple result
= CGAL::surface_neighbor_coordinates_3(T, p,n,out);
assert(result.third);
norm = result.second;
break;}
case 1: {
Result_triple result =
CGAL::surface_neighbor_coordinates_3(T, p,out,I_traits(p,n));
assert(result.third);
norm = result.second; break;}
//both versions with locate:
case 2:{
start = T.locate(p);
Result_triple result = CGAL::surface_neighbor_coordinates_3(T, p,n,
out, start);
assert(result.third);
norm = result.second; break;}
case 3: {
start = T.locate(p);
Result_triple result =
CGAL::surface_neighbor_coordinates_3(T, p,out, I_traits(p,n),start);
assert(result.third);
norm = result.second; break;}
//taking all points:
case 4: {
Result_triple result
= CGAL::surface_neighbor_coordinates_3(T.points_begin(),
T.points_end(),p,n,
out,
T.geom_traits());
assert(result.third);
norm = result.second; break;}
case 5: {
Result_triple result
= CGAL::surface_neighbor_coordinates_3(T.points_begin(),
T.points_end(),p,
out ,I_traits(p,n));
assert(result.third);
norm = result.second; break;}
//the last two with certification:
case 6: {
Result_quadruple
result = CGAL::surface_neighbor_coordinates_certified_3
(T.points_begin(), T.points_end(),p,n,
out, T.geom_traits());
assert(result.third && result.fourth);
norm = result.second; break;
case 0:{
Result_triple result
= CGAL::surface_neighbor_coordinates_3(T, p,n,out);
assert(result.third);
norm = result.second;
break;}
case 1: {
Result_triple result =
CGAL::surface_neighbor_coordinates_3(T, p,out,I_traits(p,n));
assert(result.third);
norm = result.second; break;}
//both versions with locate:
case 2:{
start = T.locate(p);
Result_triple result = CGAL::surface_neighbor_coordinates_3(T, p, n, out, start);
assert(result.third);
norm = result.second; break;}
case 3: {
start = T.locate(p);
Result_triple result =
CGAL::surface_neighbor_coordinates_3(T, p, out, I_traits(p,n), start);
assert(result.third);
norm = result.second; break;}
//taking all points:
case 4: {
Result_triple result
= CGAL::surface_neighbor_coordinates_3(T.points_begin(),
T.points_end(), p, n,
out,
T.geom_traits());
assert(result.third);
norm = result.second; break;}
case 5: {
Result_triple result
= CGAL::surface_neighbor_coordinates_3(T.points_begin(),
T.points_end(), p,
out ,I_traits(p,n));
assert(result.third);
norm = result.second; break;}
//the last two with certification:
case 6: {
Result_quadruple
result = CGAL::surface_neighbor_coordinates_certified_3
(T.points_begin(), T.points_end(),p,n,
out, T.geom_traits());
assert(result.third && result.fourth);
norm = result.second; break;
}
case 7: {
Result_quadruple
result = CGAL::surface_neighbor_coordinates_certified_3
(T.points_begin(), T.points_end(),p, out ,I_traits(p,n));
assert(result.third && result.fourth);
norm = result.second;
break;
}
default:
std::cout << "Switch function calls: Nothing is tested. " <<
std::endl;
}
case 7: {
Result_quadruple
result = CGAL::surface_neighbor_coordinates_certified_3
(T.points_begin(), T.points_end(),p, out ,I_traits(p,n));
assert(result.third && result.fourth);
norm = result.second;
break;
}
default:
std::cout << "Switch function calls: Nothing is tested. " <<
std::endl;
}
assert(norm>0);
assert(norm > 0);
return std::make_pair(out,norm);
return std::make_pair(out, norm);
}
//call for test functions:
template < class Triangul>
template < class Tr>
void
test_coords_and_neighbors( const Triangul& T, const typename
Triangul::Geom_traits::Point_3& p,
const typename Triangul::Geom_traits::Vector_3& n,
const typename Triangul::Geom_traits::FT&
tolerance, const int& version )
test_coords_and_neighbors(const Tr& T, const typename
Tr::Geom_traits::Point_3& p,
const typename Tr::Geom_traits::Vector_3& n,
const typename Tr::Geom_traits::FT& tolerance, const int& version)
{
CGAL::Set_ieee_double_precision pfr;
typedef std::pair<typename Triangul::Geom_traits::Point_3,
typename Triangul::Geom_traits::FT > Point_coord_pair;
typedef std::pair<typename Tr::Geom_traits::Point_3,
typename Tr::Geom_traits::FT > Point_coord_pair;
std::vector<Point_coord_pair > coords;
typename Triangul::Geom_traits::FT norm;
std::vector<Point_coord_pair> coords;
typename Tr::Geom_traits::FT norm;
norm = test_coords(T, p, n,version, std::back_inserter(coords)).second;
assert(test_norm( coords.begin(), coords.end(),norm));
assert(test_barycenter(coords.begin(), coords.end(),norm,p,
tolerance));
assert(test_norm( coords.begin(), coords.end(), norm));
assert(test_barycenter(coords.begin(), coords.end(), norm, p, tolerance));
//All function testing surface neighbors are
// grouped together:
std::vector< typename Triangul::Geom_traits::Point_3 > neighbors;
std::vector< typename Tr::Geom_traits::Point_3 > neighbors;
test_neighbors(T, p, n,version, std::back_inserter(neighbors));
assert(compare_neighbors(coords.begin(),
coords.end(),neighbors.begin(),
neighbors.end(), T.geom_traits()));
coords.end(),neighbors.begin(),
neighbors.end(), T.geom_traits()));
//done
}
template <class Triangul>
template <class Tr>
void
_test_surface_neighbors_3_sphere( const Triangul & )
_test_surface_neighbors_3_sphere( const Tr & )
{
Triangul T;
Tr T;
int n=200, m=20;
double r = 3;
typedef typename Triangul::Geom_traits Gt;
typedef typename Gt::Point_3 Point;
typedef typename Tr::Geom_traits Gt;
typedef typename Gt::Point_3 Point;
std::vector<Point> points;
points.reserve(n+m);
@ -325,28 +310,28 @@ _test_surface_neighbors_3_sphere( const Triangul & )
//test with different calls:
int k=0;
for(int i=n;i<n+m;i++){
test_coords_and_neighbors(T, points[i],typename
Gt::Vector_3(points[i]-CGAL::ORIGIN),
typename Gt::FT(0.5),++k % 8);
for(int i=n;i<n+m;i++)
{
test_coords_and_neighbors(T, points[i],
typename Gt::Vector_3(points[i]-CGAL::ORIGIN),
typename Gt::FT(0.5),++k % 8);
}
}
////cube case in the interior of a face
template <class Triangul, class Transformation>
template <class Tr, class Transformation>
void
_test_surface_neighbors_3_cube(const Triangul &, const Transformation&
transform, const int n = 75,
typename Triangul::Geom_traits::FT tolerance
= typename Triangul::Geom_traits::FT(1e-29),
bool grid=true)
_test_surface_neighbors_3_cube(const Tr &, const Transformation&
transform, const int n = 75,
typename Tr::Geom_traits::FT tolerance = typename Tr::Geom_traits::FT(1e-29),
bool grid=true)
{
Triangul T;
Tr T;
int m=10;
double r = 3;
typedef typename Triangul::Geom_traits Gt;
typedef typename Tr::Geom_traits Gt;
typedef typename Gt::FT Coord_type;
typedef typename Gt::Point_3 Point;
typedef typename Gt::Point_2 Point_2;
@ -355,17 +340,19 @@ _test_surface_neighbors_3_cube(const Triangul &, const Transformation&
//data points: generate random points in a square of length r
std::vector<Point_2> points_2_data;
points_2_data.reserve(n);
if(grid)
CGAL::points_on_square_grid_2
(r,n,std::back_inserter(points_2_data),
CGAL::Creator_uniform_2<Coord_type,Point_2>());
else{
CGAL::Random_points_in_square_2<Point_2> g(r);
CGAL::cpp11::copy_n( g, n, std::back_inserter(points_2_data));
{
CGAL::points_on_square_grid_2(r, n, std::back_inserter(points_2_data),
CGAL::Creator_uniform_2<Coord_type,Point_2>());
}
for(int i=0; i < n; i++){
else
{
CGAL::Random_points_in_square_2<Point_2> g(r);
CGAL::cpp11::copy_n(g, n, std::back_inserter(points_2_data));
}
for(int i=0; i < n; i++)
{
T.insert(transform(Point(points_2_data[i].x(),points_2_data[i].y(), -r)));
T.insert(transform(Point(points_2_data[i].x(),points_2_data[i].y(), r)));
T.insert(transform(Point(-r, points_2_data[i].x(),points_2_data[i].y())));
@ -378,39 +365,41 @@ _test_surface_neighbors_3_cube(const Triangul &, const Transformation&
std::vector<Point_2> points_2_test;
points_2_test.reserve(m);
CGAL::Random_points_in_square_2<Point_2> g2(r-1.0);
CGAL::cpp11::copy_n( g2, m, std::back_inserter(points_2_test));
CGAL::cpp11::copy_n(g2, m, std::back_inserter(points_2_test));
int k=0;
for(int i=0;i<m;i++){
for(int i=0;i<m;i++)
{
//test point on z=r plane:
test_coords_and_neighbors(T,transform(Point(points_2_test[i].x(),
points_2_test[i].y(), r)),
transform(Vector(0,0,1)),tolerance, ++k % 8);
points_2_test[i].y(), r)),
transform(Vector(0,0,1)),tolerance, ++k % 8);
//test point on x=-r plane:
test_coords_and_neighbors(T,transform(Point(-r, points_2_test[i].x(),
points_2_test[i].y())),
transform(Vector(-1,0,0)),tolerance, ++k % 8 );
points_2_test[i].y())),
transform(Vector(-1,0,0)),tolerance, ++k % 8 );
//test point on x=r plane:
test_coords_and_neighbors(T,transform(Point(r, points_2_test[i].x(),
points_2_test[i].y())),
transform(Vector(1,0,0)),tolerance,++k % 8 );
points_2_test[i].y())),
transform(Vector(1,0,0)),tolerance,++k % 8 );
//test point on y=-r plane:
test_coords_and_neighbors(T,transform(Point(points_2_test[i].x(),
-r,points_2_test[i].y())),
transform(Vector(0,-1,0)),tolerance,++k % 8 );
-r,points_2_test[i].y())),
transform(Vector(0,-1,0)),tolerance,++k % 8 );
//test point on y=r plane:
test_coords_and_neighbors(T,transform(Point(points_2_test[i].x(),
r,points_2_test[i].y())),
transform(Vector(0,1,0)),tolerance,++k % 8);
r,points_2_test[i].y())),
transform(Vector(0,1,0)),tolerance,++k % 8);
}
//test a sample point:
//with Delaunay triangulation filering:
test_coords_and_neighbors(T,transform(Point(points_2_data[n/2].x(),
points_2_data[n/2].y(), r)),
transform(Vector(0,0,1)),Coord_type(0),0);
points_2_data[n/2].y(), r)),
transform(Vector(0,0,1)),Coord_type(0),0);
//considering all points:
test_coords_and_neighbors(T,transform(Point(points_2_data[n/2].x(),
points_2_data[n/2].y(), r)),
transform(Vector(0,0,1)),Coord_type(0),4);
points_2_data[n/2].y(), r)),
transform(Vector(0,0,1)),Coord_type(0),4);
}

View File

@ -16,11 +16,12 @@
//
// Author(s) : Naceur MESKINI.
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Delaunay_triangulation_2.h>
#include <CGAL/_test_natural_neighbors_2.cpp>
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Delaunay_triangulation_2.h>
#include <iostream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;