songsenand
/
cgal
mirror of https://github.com/CGAL/cgal
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Merge pull request #4421 from MaelRL/Triangulation_on_sphere_2-GF 

new package: Triangulation on the sphere
This commit is contained in:
Sebastien Loriot 2021-04-21 15:24:45 +02:00 committed by GitHub
parent 5b253ad2fc 8226666398
commit b8bc56ca12
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
88 changed files with 39602 additions and 145 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

3
BGL/include/CGAL/boost/graph/parameters_interface.h
View File

@ -39,7 +39,6 @@ CGAL_add_named_parameter(face_to_face_map_t, face_to_face_map, face_to_face_map)
CGAL_add_named_parameter(implementation_tag_t, implementation_tag, implementation_tag)
CGAL_add_named_parameter(prevent_unselection_t, prevent_unselection, prevent_unselection)
CGAL_add_named_parameter(stream_precision_t, stream_precision, stream_precision)
CGAL_add_named_parameter(verbose_t, verbose, verbose)
// List of named parameters used for IO
@ -52,6 +51,8 @@ CGAL_add_named_parameter(vertex_color_map_t, vertex_color_map, vertex_color_map)
CGAL_add_named_parameter(vertex_texture_map_t, vertex_texture_map, vertex_texture_map)
CGAL_add_named_parameter(face_color_map_t, face_color_map, face_color_map)
CGAL_add_named_parameter(repair_polygon_soup_t, repair_polygon_soup, repair_polygon_soup)
CGAL_add_named_parameter(output_color_t, output_color, output_color)
CGAL_add_named_parameter(stream_precision_t, stream_precision, stream_precision)
// List of named parameters that we use in the package 'Mesh_3'
CGAL_add_named_parameter(vertex_feature_degree_t, vertex_feature_degree, vertex_feature_degree_map)

6
Cartesian_kernel/include/CGAL/Cartesian/function_objects.h
View File

@ -3692,9 +3692,9 @@ namespace CartesianKernelFunctors {
// p,q,r supposed to be non collinear
// tests whether s is on the same side of p,q as r
// returns :
// COLLINEAR if pqr collinear
// POSITIVE if qrp and qrs have the same orientation
// NEGATIVE if qrp and qrs have opposite orientations
// COLLINEAR if qps collinear
// POSITIVE if qpr and qps have the same orientation
// NEGATIVE if qpr and qps have opposite orientations
CGAL_kernel_exactness_precondition( ! cl(p, q, r) );
CGAL_kernel_exactness_precondition( cp(p, q, r, s) );
return coplanar_orientationC3(p.x(), p.y(), p.z(),

4
Circular_kernel_3/demo/Circular_kernel_3/Viewer.cpp
View File

@ -295,9 +295,6 @@ void Viewer::initialize_buffers()
buffers[6].release();
}
vao[2].release();
}
void Viewer::compute_elements()
@ -841,7 +838,6 @@ void Viewer::draw()
rendering_program.release();
vao[1].release();
}
void Viewer::init()

4
Circular_kernel_3/include/CGAL/Circular_kernel_3/Circular_arc_3.h
View File

@ -73,8 +73,8 @@ namespace CGAL {
: _full(false)
{
// l must pass through s and t, and s != t
CGAL_kernel_precondition(SK().has_on_3_object()(c,s));
CGAL_kernel_precondition(SK().has_on_3_object()(c,t));
// CGAL_kernel_precondition(SK().has_on_3_object()(c,s));
// CGAL_kernel_precondition(SK().has_on_3_object()(c,t));
CGAL_kernel_precondition(s != t);
base = Rep(c,s,t);
// we can optimize the computations of the sign (for the has_on functor),

1
Documentation/doc/Documentation/packages.txt
View File

@ -67,6 +67,7 @@
\package_listing{Triangulation_2}
\package_listing{TDS_2}
\package_listing{Triangulation_on_sphere_2}
\package_listing{Periodic_2_triangulation_2}
\package_listing{Hyperbolic_triangulation_2}
\package_listing{Periodic_4_hyperbolic_triangulation_2}

9
Documentation/doc/biblio/cgal_manual.bib
View File

@ -353,6 +353,15 @@ Boissonnat}
,url = "https://tel.archives-ouvertes.fr/tel-00552215/"
}
@inproceedings{cgal:ccplr-redtp-10,
title={Robust and efficient {Delaunay} triangulations of points on or close to a sphere},
author={Caroli, Manuel and de Castro, Pedro MM and Loriot, S{\'e}bastien and Rouiller, Olivier and Teillaud, Monique and Wormser, Camille},
booktitle={International Symposium on Experimental Algorithms},
pages={462--473},
year={2010},
organization={Springer}
}
@article{cgal:cc-rgbss-78
,author = {E.~Catmull and J.~Clark}
,title = {Recursively generated {B}-spline surfaces

6
Homogeneous_kernel/include/CGAL/Homogeneous/function_objects.h
View File

@ -3837,9 +3837,9 @@ namespace HomogeneousKernelFunctors {
// p,q,r supposed to be non collinear
// tests whether s is on the same side of p,q as r
// returns :
// COLLINEAR if pqr collinear
// POSITIVE if qrp and qrs have the same orientation
// NEGATIVE if qrp and qrs have opposite orientations
// COLLINEAR if qps collinear
// POSITIVE if qpr and qps have the same orientation
// NEGATIVE if qpr and qps have opposite orientations
CGAL_kernel_exactness_precondition( ! cl(p, q, r) );
CGAL_kernel_exactness_precondition( cp(p, q, r, s) );

12
Installation/CHANGES.md
View File

@ -6,8 +6,18 @@ Release History
Release date: June 2021
### Triangulations on the Sphere (new package)
- This package enables the construction and manipulation of Delaunay triangulations on the 2-sphere.
Triangulations are built incrementally and can be modified by insertion or removal of vertices.
Point location querying and primitives to build the dual Voronoi diagram are provided.
### General
- Cmake minimal version is now 3.14
- Cmake minimal version is now `3.14`.
### File Input / Output
- Point set, polygon soup, and polygon mesh file I/O functions have been harmonized and documented:
- Point set I/O functions can be found in the packages [Point_set_processing_3](https://doc.cgal.org/5.3/Manual/packages.html#PkgPolygonMeshProcessing), and [Point_set_3](https://doc.cgal.org/5.3/Manual/packages.html#PkgPointSet3).
- Polygon mesh I/O functions can be found in the package [BGL](https://doc.cgal.org/5.3/Manual/packages.html#PkgBGL).

54
Installation/include/CGAL/license/Triangulation_on_sphere_2.h Normal file
View File

@ -0,0 +1,54 @@
// Copyright (c) 2016 GeometryFactory SARL (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org)
//
// $URL$
// $Id$
// SPDX-License-Identifier: LGPL-3.0-or-later OR LicenseRef-Commercial
//
// Author(s) : Andreas Fabri
//
// Warning: this file is generated, see include/CGAL/licence/README.md
#ifndef CGAL_LICENSE_TRIANGULATION_ON_SPHERE_2_H
#define CGAL_LICENSE_TRIANGULATION_ON_SPHERE_2_H
#include <CGAL/config.h>
#include <CGAL/license.h>
#ifdef CGAL_TRIANGULATION_ON_SPHERE_2_COMMERCIAL_LICENSE
# if CGAL_TRIANGULATION_ON_SPHERE_2_COMMERCIAL_LICENSE < CGAL_RELEASE_DATE
# if defined(CGAL_LICENSE_WARNING)
CGAL_pragma_warning("Your commercial license for CGAL does not cover "
"this release of the 2D Triangulations on Sphere package.")
# endif
# ifdef CGAL_LICENSE_ERROR
# error "Your commercial license for CGAL does not cover this release \
of the 2D Triangulations on Sphere package. \
You get this error, as you defined CGAL_LICENSE_ERROR."
# endif // CGAL_LICENSE_ERROR
# endif // CGAL_TRIANGULATION_ON_SPHERE_2_COMMERCIAL_LICENSE < CGAL_RELEASE_DATE
#else // no CGAL_TRIANGULATION_ON_SPHERE_2_COMMERCIAL_LICENSE
# if defined(CGAL_LICENSE_WARNING)
CGAL_pragma_warning("\nThe macro CGAL_TRIANGULATION_ON_SPHERE_2_COMMERCIAL_LICENSE is not defined."
"\nYou use the CGAL 2D Triangulations on Sphere package under "
"the terms of the GPLv3+.")
# endif // CGAL_LICENSE_WARNING
# ifdef CGAL_LICENSE_ERROR
# error "The macro CGAL_TRIANGULATION_ON_SPHERE_2_COMMERCIAL_LICENSE is not defined.\
You use the CGAL 2D Triangulations on Sphere package under the terms of \
the GPLv3+. You get this error, as you defined CGAL_LICENSE_ERROR."
# endif // CGAL_LICENSE_ERROR
#endif // no CGAL_TRIANGULATION_ON_SPHERE_2_COMMERCIAL_LICENSE
#endif // CGAL_LICENSE_TRIANGULATION_ON_SPHERE_2_H

1
Installation/include/CGAL/license/gpl_package_list.txt
View File

@ -90,6 +90,7 @@ TDS_2 2D Triangulation Data Structure
TDS_3 3D Triangulation Data Structure
Three Three
Triangulation_2 2D Triangulation
Triangulation_on_sphere_2 2D Triangulation on Sphere
Triangulation_3 3D Triangulations
Triangulation dD Triangulations
Visibility_2 2D Visibility Computation

16
Intersections_3/include/CGAL/Intersections_3/internal/Triangle_3_Sphere_3_do_intersect.h
View File

@ -77,42 +77,42 @@ template <class K>
inline
typename K::Boolean
do_intersect(const typename K::Sphere_3 &sp,
const typename K::Ray_3 &lin,
const typename K::Ray_3 &ray,
const K & /* k */)
{
return squared_distance(sp.center(), lin) <= sp.squared_radius();
return squared_distance(sp.center(), ray) <= sp.squared_radius();
}
template <class K>
inline
typename K::Boolean
do_intersect(const typename K::Ray_3 &lin,
do_intersect(const typename K::Ray_3 &ray,
const typename K::Sphere_3 &sp,
const K & /* k */)
{
return squared_distance(sp.center(), lin) <= sp.squared_radius();
return squared_distance(sp.center(), ray) <= sp.squared_radius();
}
template <class K>
inline
typename K::Boolean
do_intersect(const typename K::Sphere_3 &sp,
const typename K::Segment_3 &lin,
const typename K::Segment_3 &seg,
const K & /* k */)
{
return squared_distance(sp.center(), lin) <= sp.squared_radius();
return squared_distance(sp.center(), seg) <= sp.squared_radius();
}
template <class K>
inline
typename K::Boolean
do_intersect(const typename K::Segment_3 &lin,
do_intersect(const typename K::Segment_3 &seg,
const typename K::Sphere_3 &sp,
const K & /* k */)
{
return squared_distance(sp.center(), lin) <= sp.squared_radius();
return squared_distance(sp.center(), seg) <= sp.squared_radius();
}
} // namespace internal

2
Kernel_23/doc/Kernel_23/Concepts/FunctionObjectConcepts.h
View File

@ -6905,7 +6905,7 @@ public:
/// @{
/*!
computes the unit normal of the vectors`q-p` and `r-p`.
computes the unit normal of the vectors `q-p` and `r-p`.
This requires that `Kernel::FT` supports the `sqrt` operation.
*/
Kernel::Vector_3 operator()(const Kernel::Point_3& p,

2
Periodic_3_triangulation_3/doc/Periodic_3_triangulation_3/CGAL/Periodic_3_triangulation_hierarchy_3.h
View File

@ -12,7 +12,7 @@ location queries.
implementation, only `Periodic_3_Delaunay_triangulation_3` is supported for.</I>
`PTr::Vertex` has to be a model of the concept
`Periodic_3TriangulationHierarchyVertexBase_3`.
`Periodic_3TriangulationDSVertexBase_3`.
`PTr::Geom_traits` has to be a model of the concept
`Periodic_3DelaunayTriangulationTraits_3`.

1
Periodic_3_triangulation_3/include/CGAL/Periodic_3_triangulation_3.h
View File

@ -22,6 +22,7 @@
#include <CGAL/basic.h>
#include <CGAL/internal/Periodic_3_triangulation_iterators_3.h>
#include <CGAL/result_of.h>
#include <CGAL/Periodic_3_triangulation_ds_cell_base_3.h>
#include <CGAL/Periodic_3_triangulation_ds_vertex_base_3.h>
#include <CGAL/Periodic_3_triangulation_traits_3.h>

39
Spatial_sorting/include/CGAL/Hilbert_sort_on_sphere_3.h
View File

@ -22,11 +22,13 @@
namespace CGAL {
template <class K, class Hilbert_policy >
class Hilbert_sort_on_sphere_3 {
typedef typename K::Point_3 Point_3;
static const double _sqrt_of_one_over_three;
template <class K,
class Hilbert_policy,
class P = typename K::Point_3>
class Hilbert_sort_on_sphere_3
{
typedef P Point_3;
typedef typename K::FT FT;
// Face 1, x > sqrt(1/3)
// Face 2, y > sqrt(1/3)
@ -42,28 +44,33 @@ class Hilbert_sort_on_sphere_3 {
typedef internal::Transform_coordinates_traits_3<K,-1,0,1,1> Face_5_traits_3; // -z +x
typedef internal::Transform_coordinates_traits_3<K,1,1,0,0> Face_6_traits_3; // +x +y
Hilbert_sort_2<Face_1_traits_3, Hilbert_policy > _hs_1_object;
Hilbert_sort_2<Face_2_traits_3, Hilbert_policy > _hs_2_object;
Hilbert_sort_2<Face_3_traits_3, Hilbert_policy > _hs_3_object;
Hilbert_sort_2<Face_4_traits_3, Hilbert_policy > _hs_4_object;
Hilbert_sort_2<Face_5_traits_3, Hilbert_policy > _hs_5_object;
Hilbert_sort_2<Face_6_traits_3, Hilbert_policy > _hs_6_object;
K _k;
Point_3 _p;
double _sq_r;
FT _sq_r;
const FT _sqrt_of_one_over_three;
public:
Hilbert_sort_on_sphere_3 (const K &k=K(),
double sq_r = 1.0,
const Point_3 &p = Point_3(0,0,0),
std::ptrdiff_t limit=1)
Hilbert_sort_on_sphere_3 (const K &k = K(),
FT sq_r = 1.0,
const Point_3 &p = Point_3(0,0,0),
std::ptrdiff_t limit=1)
: _hs_1_object(Face_1_traits_3(),limit),
_hs_2_object(Face_2_traits_3(),limit),
_hs_3_object(Face_3_traits_3(),limit),
_hs_4_object(Face_4_traits_3(),limit),
_hs_5_object(Face_5_traits_3(),limit),
_hs_6_object(Face_6_traits_3(),limit),
_k(k), _p(p), _sq_r(sq_r)
_k(k), _p(p), _sq_r(sq_r),
_sqrt_of_one_over_three(CGAL_NTS sqrt(FT(1)/FT(3)))
{
CGAL_precondition( sq_r > 0 );
}
@ -74,10 +81,10 @@ public:
typedef typename std::iterator_traits<RandomAccessIterator>::value_type Point;
std::vector< Point > vec[6];
const double mulcte = _sqrt_of_one_over_three * CGAL_NTS sqrt(_sq_r);
const double lxi = _p.x() - mulcte, lxs = _p.x() + mulcte;
const double lyi = _p.y() - mulcte, lys = _p.y() + mulcte;
const double lzs = _p.z() + mulcte;
const FT mulcte = _sqrt_of_one_over_three * CGAL_NTS sqrt(_sq_r);
const FT lxi = _p.x() - mulcte, lxs = _p.x() + mulcte;
const FT lyi = _p.y() - mulcte, lys = _p.y() + mulcte;
const FT lzs = _p.z() + mulcte;
for(RandomAccessIterator i = begin; i != end; ++i) {
const Point &p = *i;
@ -105,8 +112,6 @@ public:
*begin++ = vec[i][j];
}
};
template <class K, class Hilbert_policy >
const double Hilbert_sort_on_sphere_3<K,Hilbert_policy>::_sqrt_of_one_over_three = 0.57735026919;
} // namespace CGAL

37
Spatial_sorting/include/CGAL/Spatial_sort_traits_adapter_3.h
View File

@ -36,6 +36,39 @@ public:
typedef typename boost::property_traits<PointPropertyMap>::key_type Point_3;
typedef typename boost::call_traits<Point_3>::param_type Arg_type;
struct Compute_x_3
: public Base_traits::Compute_x_3
{
Compute_x_3(const PointPropertyMap& ppmap, const typename Base_traits::Compute_x_3& base):
Base_traits::Compute_x_3(base), ppmap_(ppmap){}
const PointPropertyMap& ppmap_;
bool operator()(Arg_type p) const {
return static_cast<const typename Base_traits::Compute_x_3*>(this)->operator()(get(ppmap_,p));
}
};
struct Compute_y_3
: public Base_traits::Compute_y_3
{
Compute_y_3(const PointPropertyMap& ppmap, const typename Base_traits::Compute_y_3& base):
Base_traits::Compute_y_3(base), ppmap_(ppmap){}
const PointPropertyMap& ppmap_;
bool operator()(Arg_type p) const {
return static_cast<const typename Base_traits::Compute_y_3*>(this)->operator()(get(ppmap_,p));
}
};
struct Compute_z_3
: public Base_traits::Compute_z_3
{
Compute_z_3(const PointPropertyMap& ppmap, const typename Base_traits::Compute_z_3& base):
Base_traits::Compute_z_3(base), ppmap_(ppmap){}
const PointPropertyMap& ppmap_;
bool operator()(Arg_type p) const {
return static_cast<const typename Base_traits::Compute_z_3*>(this)->operator()(get(ppmap_,p));
}
};
struct Less_x_3
: public Base_traits::Less_x_3
{
@ -69,6 +102,10 @@ public:
}
};
Compute_x_3 compute_x_3_object () const {return Compute_x_3(ppmap_, static_cast<const Gt*>(this)->compute_x_3_object() );}
Compute_y_3 compute_y_3_object () const {return Compute_y_3(ppmap_, static_cast<const Gt*>(this)->compute_y_3_object() );}
Compute_z_3 compute_z_3_object () const {return Compute_z_3(ppmap_, static_cast<const Gt*>(this)->compute_z_3_object() );}
Less_x_3 less_x_3_object () const {return Less_x_3(ppmap_, static_cast<const Gt*>(this)->less_x_3_object() );}
Less_y_3 less_y_3_object () const {return Less_y_3(ppmap_, static_cast<const Gt*>(this)->less_y_3_object() );}
Less_z_3 less_z_3_object () const {return Less_z_3(ppmap_, static_cast<const Gt*>(this)->less_z_3_object() );}

129
Spatial_sorting/include/CGAL/spatial_sort_on_sphere.h
View File

@ -16,29 +16,31 @@
#include <CGAL/spatial_sort.h>
#include <CGAL/hilbert_sort_on_sphere.h>
#include <CGAL/Hilbert_policy_tags.h>
namespace CGAL {
namespace internal {
template <class RandomAccessIterator, class Policy, class Kernel>
template <class RandomAccessIterator, class PolicyTag, class Kernel,
class FT = typename Kernel::FT,
class Point = typename Kernel::Point_3>
void spatial_sort_on_sphere (RandomAccessIterator begin, RandomAccessIterator end,
const Kernel &k,
Policy /*policy*/,
typename Kernel::Point_3 *,
double sq_r,
const typename Kernel::Point_3 &p,
const Hilbert_policy<PolicyTag> /*policy*/,
const FT sq_r,
const Point &p,
std::ptrdiff_t threshold_hilbert,
std::ptrdiff_t threshold_multiscale,
double ratio)
{
typedef Hilbert_sort_on_sphere_3<Kernel, Policy> Sort;
typedef Hilbert_sort_on_sphere_3<Kernel, Hilbert_policy<PolicyTag>, Point> Sort;
typedef std::iterator_traits<RandomAccessIterator> ITraits;
typedef typename ITraits::difference_type Diff_t;
boost::rand48 random;
boost::random_number_generator<boost::rand48, Diff_t> rng(random);
CGAL::cpp98::random_shuffle(begin,end, rng);
CGAL::cpp98::random_shuffle(begin, end, rng);
if (threshold_hilbert==0) threshold_hilbert=4;
if (threshold_multiscale==0) threshold_multiscale=16;
@ -50,87 +52,52 @@ void spatial_sort_on_sphere (RandomAccessIterator begin, RandomAccessIterator en
} // end of namespace internal
template <class RandomAccessIterator, class Policy, class Kernel>
template <class RandomAccessIterator, class PolicyTag,
class Kernel = typename CGAL::Kernel_traits<typename std::iterator_traits<RandomAccessIterator>::value_type>::Kernel,
class FT = typename Kernel::FT,
class Point = typename Kernel::Point_3>
void spatial_sort_on_sphere (RandomAccessIterator begin, RandomAccessIterator end,
const Hilbert_policy<PolicyTag> policy,
const FT sq_r = FT(1),
const Point &p = Point(0,0,0),
std::ptrdiff_t threshold_hilbert = 0,
std::ptrdiff_t threshold_multiscale = 0,
const double ratio = 0.)
{
internal::spatial_sort_on_sphere (begin, end, Kernel(), policy, sq_r, p,
threshold_hilbert, threshold_multiscale, ratio);
}
template <class RandomAccessIterator, class Kernel,
class FT = typename Kernel::FT,
class Point = typename Kernel::Point_3>
void spatial_sort_on_sphere (RandomAccessIterator begin, RandomAccessIterator end,
const Kernel &k,
Policy policy,
double sq_r=1.0,
const typename Kernel::Point_3 &p = typename Kernel::Point_3(0,0,0),
std::ptrdiff_t threshold_hilbert=0,
std::ptrdiff_t threshold_multiscale=0,
double ratio=0.0)
const FT sq_r = FT(1),
const Point &p = Point(0,0,0),
std::ptrdiff_t threshold_hilbert = 0,
std::ptrdiff_t threshold_multiscale = 0,
const double ratio = 0.)
{
typedef std::iterator_traits<RandomAccessIterator> ITraits;
typedef typename ITraits::value_type value_type;
internal::spatial_sort_on_sphere(begin, end, k, policy, static_cast<value_type *> (0),
sq_r,p, threshold_hilbert,threshold_multiscale,ratio);
internal::spatial_sort_on_sphere (begin, end, k,
Hilbert_sort_median_policy(), sq_r, p,
threshold_hilbert, threshold_multiscale, ratio);
}
template <class RandomAccessIterator>
template <class RandomAccessIterator,
class Kernel = typename CGAL::Kernel_traits<typename std::iterator_traits<RandomAccessIterator>::value_type>::Kernel,
class FT = typename Kernel::FT,
class Point = typename Kernel::Point_3>
void spatial_sort_on_sphere (RandomAccessIterator begin, RandomAccessIterator end,
Hilbert_sort_median_policy policy,
double sq_r=1.0,
const typename CGAL::Kernel_traits<typename std::iterator_traits<RandomAccessIterator>::value_type>::Kernel::Point_3 &p =
typename CGAL::Kernel_traits<typename std::iterator_traits<RandomAccessIterator>::value_type>::Kernel::Point_3(0,0,0),
std::ptrdiff_t threshold_hilbert=0,
std::ptrdiff_t threshold_multiscale=0,
double ratio=0.0)
const FT sq_r = FT(1),
const Point &p = Point(0,0,0),
std::ptrdiff_t threshold_hilbert = 0,
std::ptrdiff_t threshold_multiscale = 0,
const double ratio = 0.)
{
typedef std::iterator_traits<RandomAccessIterator> ITraits;
typedef typename ITraits::value_type value_type;
typedef CGAL::Kernel_traits<value_type> KTraits;
typedef typename KTraits::Kernel Kernel;
spatial_sort_on_sphere (begin, end, Kernel(), policy, sq_r, p,
threshold_hilbert,threshold_multiscale,ratio);
}
template <class RandomAccessIterator>
void spatial_sort_on_sphere (RandomAccessIterator begin, RandomAccessIterator end,
Hilbert_sort_middle_policy policy,
double sq_r=1.0,
const typename CGAL::Kernel_traits<typename std::iterator_traits<RandomAccessIterator>::value_type>::Kernel::Point_3 &p =
typename CGAL::Kernel_traits<typename std::iterator_traits<RandomAccessIterator>::value_type>::Kernel::Point_3(0,0,0),
std::ptrdiff_t threshold_hilbert=0,
std::ptrdiff_t threshold_multiscale=0,
double ratio=0.0)
{
typedef std::iterator_traits<RandomAccessIterator> ITraits;
typedef typename ITraits::value_type value_type;
typedef CGAL::Kernel_traits<value_type> KTraits;
typedef typename KTraits::Kernel Kernel;
spatial_sort_on_sphere (begin, end, Kernel(), policy, sq_r,p,
threshold_hilbert,threshold_multiscale,ratio);
}
template <class RandomAccessIterator, class Kernel>
void spatial_sort_on_sphere (RandomAccessIterator begin, RandomAccessIterator end,
const Kernel &k,
double sq_r = 1.0,
const typename Kernel::Point_3 &p = typename Kernel::Point_3(0,0,0),
std::ptrdiff_t threshold_hilbert=0,
std::ptrdiff_t threshold_multiscale=0,
double ratio=0.0)
{
spatial_sort_on_sphere (begin, end, k,
Hilbert_sort_median_policy(), sq_r, p,
threshold_hilbert,threshold_multiscale,ratio);
}
template <class RandomAccessIterator>
void spatial_sort_on_sphere (RandomAccessIterator begin, RandomAccessIterator end,
double sq_r = 1.0,
const typename CGAL::Kernel_traits<typename std::iterator_traits<RandomAccessIterator>::value_type>::Kernel::Point_3 &p =
typename CGAL::Kernel_traits<typename std::iterator_traits<RandomAccessIterator>::value_type>::Kernel::Point_3(0,0,0),
std::ptrdiff_t threshold_hilbert=0,
std::ptrdiff_t threshold_multiscale=0,
double ratio=0.0)
{
spatial_sort_on_sphere (begin, end,
Hilbert_sort_median_policy(), sq_r,p,
threshold_hilbert,threshold_multiscale,ratio);
internal::spatial_sort_on_sphere (begin, end, Kernel(),
Hilbert_sort_median_policy(), sq_r, p,
threshold_hilbert, threshold_multiscale, ratio);
}
} // end of namespace CGAL

25
TDS_2/doc/TDS_2/Concepts/TriangulationDSFaceBase_2.h
View File

@ -91,6 +91,11 @@ typedef TriangulationDataStructure_2::Vertex_handle Vertex_handle;
*/
typedef TriangulationDataStructure_2::Face_handle Face_handle;
/*!
*/
typedef TriangulationDataStructure_2::Face_data TDS_data;
/// @}
/// \name Creation
@ -187,5 +192,25 @@ void for_compact_container(void *p);
/// @}
/// \name Internal
/// \cgalAdvancedBegin
/// These functions are used internally by the triangulation data
/// structure. The user is not encouraged to use them directly as they
/// may change in the future.
/// \cgalAdvancedEnd
/// @{
/*!
*/
TDS_data& tds_data();
/*!
*/
const TDS_data& tds_data() const;
/// @}
}; /* end TriangulationDSFaceBase_2 */

72
TDS_2/doc/TDS_2/Concepts/TriangulationDataStructure_2.h
View File

@ -98,6 +98,11 @@ The face type, requirements for this type are described in concept `Triangulatio
*/
typedef unspecified_type Face;
/*!
%Face data type, requirements are described in `TriangulationDataStructure_2::Face_data`.
*/
typedef unspecified_type Face_data;
/*!
Handle to a vertex.
\cgalModels `Handle`
@ -689,12 +694,12 @@ Vertex_handle file_input( istream& is, bool skip_first=false);
/*!
reads a combinatorial triangulation data structure from `is` and assigns it to tthe triangulation data structure.
*/
istream& operator>> (istream& is, TriangulationDataStructure_3 & tds);
istream& operator>> (istream& is, TriangulationDataStructure_2 & tds);
/*!
writes `tds` into the stream `os`.
*/
ostream& operator<< (ostream& os, const TriangulationDataStructure_3 & tds);
ostream& operator<< (ostream& os, const TriangulationDataStructure_2 & tds);
/// @}
@ -992,3 +997,66 @@ std::ostream& operator<< (std::ostream& os, const TriangulationDataStructure_2::
/// @}
}; /* end TriangulationDataStructure_2::Face */
/*!
\ingroup PkgTDS2Concepts
\cgalConcept
Various algorithms using a triangulation data structure, such as Delaunay triangulations
or Alpha Shapes, must be able to associate a state to a face elemental.
For efficiency, this information must be stored directly within the face.
This class is only meant to store a state (Boolean). Consequently, the state must be the default
value (i.e. `false`) unless a setting function (`mark_in_conflict()`, etc.) has been called.
The three states are "in conflict", "on boundary", and "processed".
By default, a face is not in conflict, not on boundary, and not processed.
\sa `TriangulationDataStructure_2::Face`
*/
class TriangulationDataStructure_2::Face_data
{
public:
/// \name Setting
/// @{
/// Clears all flags: the face is neither in conflict, nor on the boundary, nor processed.
void clear();
/// Sets the "in conflict" state to `true`.
///
/// \post `is_in_conflict()` returns `true`
void mark_in_conflict();
/// Sets the "on boundary" state to `true`.
///
/// \post `is_on_boundary()` returns `true`
void mark_on_boundary();
/// Sets the "processed" state to `true`.
///
/// \post `processed()` returns `true`
void mark_processed();
/// @}
/// \name Access Functions
/// @{
/// Checks whether the face has default state (not in conflict, not on boundary, not processed).
bool is_clear();
/// Returns whether the face has been marked as "in conflict".
bool is_in_conflict();
/// Returns whether the face has been marked as "on boundary".
bool is_on_boundary();
/// Returns whether the face has been marked as "processed".
bool processed();
/// @}
}; /* end Face_data */

1
TDS_2/doc/TDS_2/PackageDescription.txt
View File

@ -52,6 +52,7 @@ These refining concepts and their models are described in Chapter
\cgalCRPSection{Concepts}
- `TriangulationDataStructure_2`
- `TriangulationDataStructure_2::Face`
- `TriangulationDataStructure_2::Face_data`
- `TriangulationDataStructure_2::Vertex`
- `TriangulationDSFaceBase_2`

16
TDS_2/include/CGAL/Triangulation_data_structure_2.h
View File

@ -67,6 +67,22 @@ public:
typedef Triangulation_data_structure_2<Vb, Fb2> Other;
};
class Face_data {
unsigned char conflict_state;
public:
Face_data() : conflict_state(0) {}
void clear() { conflict_state = 0; }
void mark_in_conflict() { conflict_state = 1; }
void mark_on_boundary() { conflict_state = 2; }
void mark_processed() { conflict_state = 1; }
bool is_clear() const { return conflict_state == 0; }
bool is_in_conflict() const { return conflict_state == 1; }
bool is_on_boundary() const { return conflict_state == 2; }
bool processed() const { return conflict_state == 1; }
};
typedef Vertex_base Vertex;
typedef Face_base Face;

5
TDS_2/include/CGAL/Triangulation_ds_face_base_2.h
View File

@ -30,6 +30,7 @@ public:
typedef TDS Triangulation_data_structure;
typedef typename TDS::Vertex_handle Vertex_handle;
typedef typename TDS::Face_handle Face_handle;
typedef typename TDS::Face_data TDS_data;
template <typename TDS2>
struct Rebind_TDS { typedef Triangulation_ds_face_base_2<TDS2> Other; };
@ -37,6 +38,7 @@ public:
private:
Vertex_handle V[3];
Face_handle N[3];
TDS_data _tds_data;
public:
Triangulation_ds_face_base_2();
@ -81,6 +83,9 @@ public:
void * for_compact_container() const {return N[0].for_compact_container(); }
void for_compact_container(void* p) { N[0].for_compact_container(p);}
// TDS internal data access functions.
TDS_data& tds_data() { return _tds_data; }
const TDS_data& tds_data() const { return _tds_data; }
static int ccw(int i) {return Triangulation_cw_ccw_2::ccw(i);}
static int cw(int i) {return Triangulation_cw_ccw_2::cw(i);}

1
TDS_3/doc/TDS_3/PackageDescription.txt
View File

@ -49,7 +49,6 @@ Section \ref TDS3secintro.)
\cgalClassifedRefPages
\cgalCRPSection{Concepts}
- `TriangulationDataStructure_3`
- `TriangulationDataStructure_3::Cell`
- `TriangulationDataStructure_3::Cell_data`

3
Triangulation_2/doc/Triangulation_2/CGAL/Delaunay_triangulation_2.h
View File

@ -330,8 +330,7 @@ Oriented_side side_of_oriented_circle(Face_handle f, const Point& p) const;
/*!
tests the validity of the triangulation as a
`Triangulation_2` and additionally tests the Delaunay property. This
method is mainly useful for debugging Delaunay triangulation
algorithms designed by the user.
method is mainly useful for debugging Delaunay triangulation algorithms.
*/
bool is_valid(bool verbose = false, int level = 0) const;

18
Triangulation_2/doc/Triangulation_2/CGAL/Triangulation_2.h
View File

@ -297,7 +297,7 @@ iterator over finite vertices.
typedef unspecified_type Finite_vertices_iterator;
/*!
iterator over the points corresponding the
iterator over the points corresponding to the
finite vertices of the triangulation.
*/
typedef unspecified_type Point_iterator;
@ -305,7 +305,7 @@ typedef unspecified_type Point_iterator;
/*!
range type for iterating over all faces (including infinite faces), with a nested
type `iterator` that has as value type `Face_handle`
type `iterator` that has as value type `Face_handle`.
*/
typedef Iterator_range<unspecified_type> All_face_handles;
@ -317,14 +317,14 @@ typedef Iterator_range<All_edges_iterator> All_edges;
/*!
range type for iterating over all vertices (including the infinite vertex), with a nested
type `iterator` that has as value type `Vertex_handle`
type `iterator` that has as value type `Vertex_handle`.
*/
typedef Iterator_range<unspecified_type> All_vertex_handles;
/*!
range type for iterating over finite faces, with a nested
type `iterator` that has as value type `Face_handle`
type `iterator` that has as value type `Face_handle`.
*/
typedef Iterator_range<unspecified_type> Finite_face_handles;
@ -336,7 +336,7 @@ typedef Iterator_range<Finite_edges_iterator> Finite_edges;
/*!
range type for iterating over finite vertices, with a nested
type `iterator` that has as value type `Vertex_handle`
type `iterator` that has as value type `Vertex_handle`.
*/
typedef Iterator_range<unspecified_type> Finite_vertex_handles;
@ -929,7 +929,7 @@ Point_iterator points_end() const;
/*!
returns a range of iterators over finite vertices.
\note While the value type of `Finite_vertices_iterator` is `Vertex`, the value type of
`Finite_vertex_handles::iterator` is `Vertex_handle`
`Finite_vertex_handles::iterator` is `Vertex_handle`.
*/
Finite_vertex_handles finite_vertex_handles() const;
@ -941,7 +941,7 @@ Finite_edges finite_edges() const;
/*!
returns a range of iterators over finite faces.
\note While the value type of `Finite_faces_iterator` is `Face`, the value type of
`Finite_face_handles::iterator` is `Face_handle`
`Finite_face_handles::iterator` is `Face_handle`.
*/
Finite_face_handles finite_face_handles() const;
@ -995,7 +995,7 @@ All_faces_iterator all_faces_end() const;
/*!
returns a range of iterators over all vertices.
\note While the value type of `All_vertices_iterator` is `Vertex`, the value type of
`All_vertex_handles::iterator` is `Vertex_handle`
`All_vertex_handles::iterator` is `Vertex_handle`.
*/
All_vertex_handles all_vertex_handles() const;
@ -1007,7 +1007,7 @@ All_edges all_edges() const;
/*!
returns a range of iterators over all faces.
\note While the value type of `All_faces_iterator` is `Face`, the value type of
`All_face_handles::iterator` is `Face_handle`
`All_face_handles::iterator` is `Face_handle`.
*/
All_face_handles all_face_handles() const;

16
Triangulation_2/doc/Triangulation_2/Triangulation_2.txt
View File

@ -10,7 +10,7 @@ namespace CGAL {
\image html tr1dt1.svg
\image latex tr1dt1.png
This chapter describes the two dimensional triangulations
This chapter describes the two-dimensional triangulations
of \cgal.
Section \ref Section_2D_Triangulations_Definitions recalls the
main definitions about triangulations.
@ -46,7 +46,7 @@ and vertices.
\section Section_2D_Triangulations_Definitions Definitions
A two dimensional triangulation can be roughly described as a set \f$ T\f$
A two-dimensional triangulation can be roughly described as a set \f$ T\f$
of triangular facets such that:
- two facets either are disjoint or share a lower dimensional
@ -85,8 +85,8 @@ A point \f$ p\f$ in the domain of \f$ T\f$ is said to singular
if its surrounding in \f$ U_T\f$
is neither a topological ball nor a topological disc.
Then, a two dimensional triangulation can be described as a
two dimensional simplicial complex that is pure,
Then, a two-dimensional triangulation can be described as a
two-dimensional simplicial complex that is pure,
connected and without singularity.
Each facet of a triangulation can be given an orientation
@ -100,7 +100,7 @@ that all pairs of incident facets have consistent orientations.
The data structure underlying \cgal triangulations
allows the user to represent the combinatorics of
any orientable two dimensional triangulations
any orientable two-dimensional triangulations
without boundaries.
On top of this data structure, the 2D triangulations classes
take care of the geometric embedding of the triangulation
@ -126,7 +126,7 @@ to design a face of any dimension,
and the two-dimensional faces of a triangulation
should be properly called <I>facets</I>.
However, following a common usage, we hereafter often call <I>faces</I>, the facets
of a two dimensional triangulation.
of a two-dimensional triangulation.
\section Section_2D_Triangulations_Representation Representation
@ -731,7 +731,7 @@ regular triangulation of \f$ { PW}\f$.
Alternatively, the regular triangulation
of the weighted points set \f$ { PW}\f$
can be obtained as the projection
on the two dimensional plane of the convex hull of the set of three
on the two-dimensional plane of the convex hull of the set of three
dimensional points
\f$ { P'}= \{ (p_i,p_i ^2 - w_i ) | i = 1, \ldots , n \}\f$.
@ -759,7 +759,7 @@ they have to be stored somewhere.
The regular triangulation store those points in special vertices, called
hidden vertices.
A hidden point can reappear as vertex of the triangulation
only when the two dimensional face that hides it
only when the two-dimensional face that hides it
is removed from the triangulation. To deal with this feature,
each face of a regular triangulation stores a list of hidden vertices.
The points in those vertices

5
Triangulation_2/include/CGAL/Triangulation_2.h
View File

@ -44,6 +44,7 @@
#include <boost/random/variate_generator.hpp>
#include <CGAL/boost/iterator/transform_iterator.hpp>
#include <boost/iterator/zip_iterator.hpp>
#include <boost/utility/result_of.hpp>
#ifndef CGAL_TRIANGULATION_2_DONT_INSERT_RANGE_OF_POINTS_WITH_INFO
#include <CGAL/internal/info_check.h>
@ -304,8 +305,10 @@ public:
typename boost::result_of<const Construct_point_2(const Point&)>::type
construct_point(const Point& p) const { return geom_traits().construct_point_2_object()(p); }
const Point& point(Face_handle c, int i) const;
const Point& point(Face_handle f, int i) const;
const Point& point(Vertex_handle v) const;
Segment segment(Face_handle f, int i) const;
Segment segment(const Edge& e) const;
Segment segment(const Edge_circulator& ec) const;

20
Triangulation_on_sphere_2/benchmark/CMakeLists.txt Normal file
View File

@ -0,0 +1,20 @@
# Created by the script cgal_create_cmake_script
# This is the CMake script for compiling a CGAL application.
project( Triangulation_on_sphere_2_Benchmarks )
cmake_minimum_required(VERSION 3.1...3.15)
find_package(CGAL REQUIRED COMPONENTS Core )
if ( CGAL_FOUND )
create_single_source_cgal_program( "bench_dtos2.cpp" )
create_single_source_cgal_program( "generate_points.cpp" )
else()
message(STATUS "This program requires the CGAL library, and will not be compiled.")
endif()

107
Triangulation_on_sphere_2/benchmark/bench_dtos2.cpp Normal file
View File

@ -0,0 +1,107 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Delaunay_triangulation_on_sphere_traits_2.h>
#include <CGAL/Delaunay_triangulation_on_sphere_2.h>
#include <CGAL/Projection_on_sphere_traits_3.h>
#include <CGAL/algorithm.h>
#include <CGAL/convex_hull_3.h>
#include <CGAL/Delaunay_triangulation_3.h>
#include <CGAL/point_generators_2.h>
#include <CGAL/point_generators_3.h>
#include <CGAL/Surface_mesh.h>
#include <CGAL/squared_distance_3.h>
#include <CGAL/Timer.h>
#include <boost/iterator/transform_iterator.hpp>
#include <cmath>
#include <fstream>
#include <vector>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Surface_mesh<K> Surface_mesh;
typedef K::Segment_3 Segment_3;
typedef CGAL::Delaunay_triangulation_3<K> Delaunay;
typedef CGAL::Delaunay_triangulation_on_sphere_traits_2<K> Gt;
typedef CGAL::Projection_on_sphere_traits_3<K> Gt2;
typedef CGAL::Delaunay_triangulation_on_sphere_2<Gt> DTOS;
typedef CGAL::Delaunay_triangulation_on_sphere_2<Gt2> DTOS2;
typedef K::Point_3 Point;
typedef CGAL::Delaunay_triangulation_3<K, CGAL::Fast_location> Delaunay_fast;
typedef CGAL::Creator_uniform_3<double, Point> Creator;
int main(int, char**)
{
CGAL::Timer time;
const std::size_t nu_of_pts = 1e7;
const double radius = 5184.152;
CGAL::Random random;
std::cout << "Seed is " << random.get_seed() << std::endl;
CGAL::Random_points_on_sphere_3<Point, Creator> on_sphere(radius, random);
std::vector<Point> points;
points.reserve(nu_of_pts);
for(std::size_t count=0; count<nu_of_pts; ++count)
points.push_back(*on_sphere++);
std::cout << points.size() << " points" << std::endl;
// Delaunay_traits
DTOS dtos;
dtos.set_radius(radius);
std::cout << " ***STARTING***" << std::endl;
time.start();
dtos.insert(points.begin(), points.end());
time.stop();
assert(dtos.number_of_vertices() == nu_of_pts);
std::cout << "Triangulation sphere: "
<< dtos.number_of_vertices() << " vertices in " << time.time() << " sec" << std::endl;
//Triangulation with points on the sphere (projection_traits)
Gt2 traits(K::Point_3(0, 0, 0), radius);
DTOS2 dtos2(traits);
Gt2::Construct_point_on_sphere_2 cst = traits.construct_point_on_sphere_2_object();
time.reset();
time.start();
dtos2.insert(boost::make_transform_iterator(points.begin(), cst),
boost::make_transform_iterator(points.end(), cst));
time.stop();
std::cout << "Triangulation w/ sphere projection traits: "
<< dtos2.number_of_vertices() << " vertices in " << time.time() << " sec" << std::endl;
// Surface_mesh sm;
// time.reset();
// time.start();
// CGAL::convex_hull_3(points.begin(), points.end(), sm);
// time.stop();
// std::cout << "Convex hull 3D: " << time.time() << " " << std::endl;
time.reset();
time.start();
Delaunay T;
T.insert(Point(0, 0, 0));
T.insert(points.begin(), points.end());
time.stop();
std::cout << "Delaunay 3D with origin: "
<< T.number_of_vertices() << " vertices in " << time.time() << " sec" << std::endl;
time.reset();
time.start();
Delaunay_fast T_fast_on2;
T_fast_on2.insert(Point(0, 0, 0));
T_fast_on2.insert(points.begin(), points.end());
time.stop();
std::cout << "Delaunay 3D with origin, fast location: " << time.time() << std::endl;
return EXIT_SUCCESS;
}

38
Triangulation_on_sphere_2/benchmark/generate_points.cpp Normal file
View File

@ -0,0 +1,38 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Random.h>
#include <CGAL/point_generators_3.h>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
int main(int argc, char* argv[])
{
int no_of_pts;
double radius;
if(argc > 1)
{
no_of_pts = pow(10, atoi(argv[1]));
radius = atoi(argv[2]);
}
else
{
no_of_pts = 1000000;
radius = 1;
}
CGAL::Random random(7);
typedef CGAL::Creator_uniform_3<double, K::Point_3> Creator;
std::cout << no_of_pts << std::endl;
CGAL::Random_points_on_sphere_3<K::Point_3, Creator> on_sphere(radius, random);
for(int count=0; count<no_of_pts; ++count)
{
K::Point_3 p = *on_sphere;
++on_sphere;
std::cout << p << std::endl;
}
return EXIT_SUCCESS;
}

49
Triangulation_on_sphere_2/demo/Triangulation_on_sphere_2/CMakeLists.txt Normal file
View File

@ -0,0 +1,49 @@
# Created by the script cgal_create_cmake_script
# This is the CMake script for compiling a CGAL application.
cmake_minimum_required(VERSION 3.1...3.15)
project (Triangulation_on_sphere_2_Demo)
set(CMAKE_ALLOW_LOOSE_LOOP_CONSTRUCTS true)
# Find includes in corresponding build directories
set(CMAKE_INCLUDE_CURRENT_DIR ON)
# Instruct CMake to run moc automatically when needed.
set(CMAKE_AUTOMOC ON)
if(NOT POLICY CMP0070 AND POLICY CMP0053)
# Only set CMP0053 to OLD with CMake<3.10, otherwise there is a warning.
cmake_policy(SET CMP0053 OLD)
endif()
if(POLICY CMP0071)
cmake_policy(SET CMP0071 NEW)
endif()
# Find CGAL and CGAL Qt5
find_package(CGAL REQUIRED OPTIONAL_COMPONENTS Qt5)
# Find Qt5 itself
find_package(Qt5 QUIET COMPONENTS Script OpenGL Gui Svg)
find_package(Eigen3 3.1.0) #(requires 3.1.0 or greater)
include(CGAL_Eigen3_support)
if(CGAL_Qt5_FOUND AND Qt5_FOUND AND TARGET CGAL::Eigen3_support)
# Include this package's headers first
include_directories(BEFORE ./ ./include)
# ui file, created wih Qt Designer
qt5_wrap_ui( uis Mainwindow.ui )
#qt5_generate_moc( main.cpp Mainwindow.moc)
add_executable ( Triangulation_on_sphere_2_Demo main.cpp Viewer.cpp ${uis})
add_to_cached_list( CGAL_EXECUTABLE_TARGETS Triangulation_on_sphere_2_Demo )
target_link_libraries( Triangulation_on_sphere_2_Demo PRIVATE CGAL::CGAL CGAL::CGAL_Qt5 CGAL::Eigen3_support)
else()
message(STATUS "NOTICE: This demo requires CGAL, the QGLViewer, OpenGL and Qt4, and will not be compiled.")
endif()

68
Triangulation_on_sphere_2/demo/Triangulation_on_sphere_2/Mainwindow.ui Normal file
View File

@ -0,0 +1,68 @@
<?xml version="1.0" encoding="UTF-8"?>
<ui version="4.0">
<class>MainWindow</class>
<widget class="QMainWindow" name="MainWindow">
<property name="geometry">
<rect>
<x>0</x>
<y>0</y>
<width>800</width>
<height>600</height>
</rect>
</property>
<property name="windowTitle">
<string>2D Triangulation on 3D Sphere</string>
</property>
<widget class="QWidget" name="centralwidget">
<layout class="QGridLayout" name="gridLayout">
<item row="0" column="0">
<widget class="Viewer" name="viewer" native="true"/>
</item>
</layout>
</widget>
<widget class="QMenuBar" name="menubar">
<property name="geometry">
<rect>
<x>0</x>
<y>0</y>
<width>800</width>
<height>20</height>
</rect>
</property>
<widget class="QMenu" name="menu_File">
<property name="title">
<string>&amp;File</string>
</property>
<addaction name="action_Open"/>
<addaction name="action_Quit"/>
</widget>
<addaction name="menu_File"/>
</widget>
<widget class="QStatusBar" name="statusbar"/>
<action name="action_Quit">
<property name="text">
<string>&amp;Quit</string>
</property>
<property name="shortcut">
<string>Ctrl+Q</string>
</property>
</action>
<action name="action_Open">
<property name="text">
<string>&amp;Open</string>
</property>
<property name="shortcut">
<string>Ctrl+O</string>
</property>
</action>
</widget>
<customwidgets>
<customwidget>
<class>Viewer</class>
<extends>QWidget</extends>
<header>Viewer.h</header>
</customwidget>
</customwidgets>
<resources/>
<connections/>
</ui>

448
Triangulation_on_sphere_2/demo/Triangulation_on_sphere_2/Viewer.cpp Normal file
View File

@ -0,0 +1,448 @@
#include "Viewer.h"
#include <CGAL/point_generators_3.h>
#include <CGAL/squared_distance_3.h>
#include <CGAL/Exact_spherical_kernel_3.h>
#include <CGAL/Qt/CreateOpenGLContext.h>
#include "create_sphere.h"
#include <vector>
Viewer::Viewer(QWidget* parent )
: CGAL::QGLViewer(parent)
{
draw_balls = true;
draw_inputs = false;
trivial_center.push_back(0.);
trivial_center.push_back(0.);
trivial_center.push_back(0.);
radius_ = 1.0;
center_ = Point_3(0,0,0);
}
void Viewer::compile_shaders()
{
initializeOpenGLFunctions();
for(int i = 0; i<SIZE_OF_VBO; ++i)
if(! buffers[i].create())
{
std::cerr<<"VBO Creation FAILED"<<std::endl;
return;
}
for(int i = 0; i<SIZE_OF_VAO; ++i)
if(!vao[i].create())
{
std::cerr<<"VAO Creation FAILED"<<std::endl;
return;
}
//The sphere
//Vertex source code
const char vertex_source[] =
{
"#version 120 \n"
"attribute highp vec4 vertex;\n"
"attribute highp vec3 normal;\n"
"attribute highp vec4 center;\n"
"uniform highp mat4 mvp_matrix;\n"
"uniform highp mat4 mv_matrix; \n"
"varying highp vec4 fP; \n"
"varying highp vec3 fN; \n"
"void main(void)\n"
"{\n"
" fP = mv_matrix * (vertex+vec4(center.xyz, 0.0)); \n"
" fN = mat3(mv_matrix)* normal; \n"
" gl_Position = mvp_matrix * (vertex+vec4(center.xyz, 0.0));\n"
"}"
};
//Vertex source code
const char fragment_source[] =
{
"#version 120 \n"
"varying highp vec4 fP; \n"
"varying highp vec3 fN; \n"
"uniform highp vec4 color; \n"
"uniform highp vec4 light_pos; \n"
"uniform highp vec4 light_diff; \n"
"uniform highp vec4 light_spec; \n"
"uniform highp vec4 light_amb; \n"
"uniform float spec_power ; \n"
"void main(void) { \n"
" highp vec3 L = light_pos.xyz - fP.xyz; \n"
" highp vec3 V = -fP.xyz; \n"
" highp vec3 N = normalize(fN); \n"
" L = normalize(L); \n"
" V = normalize(V); \n"
" highp vec3 R = reflect(-L, N); \n"
" highp vec4 diffuse = max(dot(N,L), 0.0) * light_diff * color; \n"
" highp vec4 specular = pow(max(dot(R,V), 0.0), spec_power) * light_spec; \n"
"gl_FragColor = light_amb*color + diffuse ; \n"
"} \n"
"\n"
};
QOpenGLShader *vertex_shader = new QOpenGLShader(QOpenGLShader::Vertex);
if(!vertex_shader->compileSourceCode(vertex_source))
{
std::cerr<<"Compiling vertex source FAILED"<<std::endl;
}
QOpenGLShader *fragment_shader= new QOpenGLShader(QOpenGLShader::Fragment);
if(!fragment_shader->compileSourceCode(fragment_source))
{
std::cerr<<"Compiling fragmentsource FAILED"<<std::endl;
}
if(!rendering_program.addShader(vertex_shader))
{
std::cerr<<"adding vertex shader FAILED"<<std::endl;
}
if(!rendering_program.addShader(fragment_shader))
{
std::cerr<<"adding fragment shader FAILED"<<std::endl;
}
if(!rendering_program.link())
{
std::cerr<<"linking Program FAILED"<<std::endl;
}
//The edges
//Vertex source code
const char vertex_source_edges[] =
{
"#version 120 \n"
"attribute highp vec4 vertex;\n"
"uniform highp mat4 mvp_matrix;\n"
"void main(void)\n"
"{\n"
" gl_Position = mvp_matrix * vertex;\n"
"}"
};
//Vertex source code
const char fragment_source_edges[] =
{
"#version 120 \n"
"void main(void) { \n"
"gl_FragColor = vec4(0.0, 0.0, 0.0, 1.0); \n"
"} \n"
"\n"
};
QOpenGLShader *vertex_shader_edges = new QOpenGLShader(QOpenGLShader::Vertex);
if(!vertex_shader_edges->compileSourceCode(vertex_source_edges))
{
std::cerr<<"Compiling vertex source FAILED"<<std::endl;
}
QOpenGLShader *fragment_shader_edges = new QOpenGLShader(QOpenGLShader::Fragment);
if(!fragment_shader_edges->compileSourceCode(fragment_source_edges))
{
std::cerr<<"Compiling fragmentsource FAILED"<<std::endl;
}
if(!rendering_program_edges.addShader(vertex_shader_edges))
{
std::cerr<<"adding vertex shader FAILED"<<std::endl;
}
if(!rendering_program_edges.addShader(fragment_shader_edges))
{
std::cerr<<"adding fragment shader FAILED"<<std::endl;
}
if(!rendering_program_edges.link())
{
std::cerr<<"linking Program FAILED"<<std::endl;
}
}
void Viewer::initialize_buffers()
{
rendering_program.bind();
//The big white sphere
vao[SPHERE_VAO].bind();
//points of the sphere
buffers[SPHERE_POINTS].bind();
buffers[SPHERE_POINTS].allocate(pos_sphere.data(),
static_cast<int>(pos_sphere.size()*sizeof(float)));
vertexLocation[0] = rendering_program.attributeLocation("vertex");
rendering_program.enableAttributeArray(vertexLocation[0]);
rendering_program.setAttributeBuffer(vertexLocation[0],GL_FLOAT,0,3);
buffers[SPHERE_POINTS].release();
//normals of the sphere
buffers[SPHERE_NORMALS].bind();
buffers[SPHERE_NORMALS].allocate(normals.data(),
static_cast<int>(normals.size()*sizeof(float)));
normalsLocation[0] = rendering_program.attributeLocation("normal");
rendering_program.enableAttributeArray(normalsLocation[0]);
rendering_program.setAttributeBuffer(normalsLocation[0],GL_FLOAT,0,3);
buffers[SPHERE_NORMALS].release();
//center of the sphere
buffers[SPHERE_CENTER].bind();
buffers[SPHERE_CENTER].allocate(trivial_center.data(),
static_cast<int>(trivial_center.size()*sizeof(float)));
trivialCenterLocation = rendering_program.attributeLocation("center");
rendering_program.enableAttributeArray(trivialCenterLocation);
rendering_program.setAttributeBuffer(trivialCenterLocation,GL_FLOAT,0,3);
buffers[SPHERE_CENTER].release();
glVertexAttribDivisor(trivialCenterLocation, 1);
glVertexAttribDivisor(normalsLocation[0], 0);
vao[SPHERE_VAO].release();
//The little green spheres
vao[POINTS_VAO].bind();
//points of the spheres
buffers[POINTS_POINTS].bind();
buffers[POINTS_POINTS].allocate(pos_sphere_inter.data(),
static_cast<int>(pos_sphere_inter.size()*sizeof(float)));
vertexLocation[2] = rendering_program.attributeLocation("vertex");
rendering_program.enableAttributeArray(vertexLocation[2]);
rendering_program.setAttributeBuffer(vertexLocation[2],GL_FLOAT,0,3);
buffers[POINTS_POINTS].release();
//normals of the sphere
buffers[POINTS_NORMALS].bind();
buffers[POINTS_NORMALS].allocate(normals_inter.data(),
static_cast<int>(normals_inter.size()*sizeof(float)));
normalsLocation[2] = rendering_program.attributeLocation("normal");
rendering_program.enableAttributeArray(normalsLocation[2]);
rendering_program.setAttributeBuffer(normalsLocation[2],GL_FLOAT,0,3);
buffers[POINTS_NORMALS].release();
//center of the sphere
buffers[POINTS_CENTERS].bind();
buffers[POINTS_CENTERS].allocate(pos_points.data(),
static_cast<int>(pos_points.size()*sizeof(float)));
centerLocation = rendering_program.attributeLocation("center");
rendering_program.enableAttributeArray(centerLocation);
rendering_program.setAttributeBuffer(centerLocation,GL_FLOAT,0,3);
buffers[POINTS_CENTERS].release();
glVertexAttribDivisor(centerLocation, 1);
glVertexAttribDivisor(normalsLocation[2], 0);
vao[POINTS_VAO].release();
rendering_program.release();
//The circles
rendering_program_edges.bind();
vao[EDGES_VAO].bind();
buffers[EDGES_POINTS].bind();
buffers[EDGES_POINTS].allocate(pos_lines.data(),
static_cast<int>(pos_lines.size()*sizeof(float)));
vertexLocation[1] = rendering_program_edges.attributeLocation("vertex");
rendering_program_edges.enableAttributeArray(vertexLocation[1]);
rendering_program_edges.setAttributeBuffer(vertexLocation[1],GL_FLOAT,0,3);
buffers[EDGES_POINTS].release();
vao[EDGES_VAO].release();
rendering_program_edges.release();
}
void Viewer::compute_elements()
{
pos_sphere_inter.clear();
normals_inter.clear();
pos_sphere.clear();
normals.clear();
trivial_center.clear();
trivial_center.push_back(static_cast<float>(center_.x()));
trivial_center.push_back(static_cast<float>(center_.y()));
trivial_center.push_back(static_cast<float>(center_.z()));
create_flat_sphere(static_cast<float>(radius_), pos_sphere, normals, 9);
create_flat_sphere(0.001f*static_cast<float>(radius_), pos_sphere_inter, normals_inter, 36);
pos_points.resize(0);
pos_lines.resize(0);
std::cout << inputs.size() << " inputs to points" << std::endl;
// draw points as small spheres
for(const auto& p : inputs)
{
pos_points.push_back(p.x());
pos_points.push_back(p.y());
pos_points.push_back(p.z());
}
std::cout << subsampled_arcs.size() << " edge arcs" << std::endl;
// draw the primal edges
std::size_t tot_size = 0;
for(const auto& arc : subsampled_arcs)
tot_size += arc.size();
pos_lines.reserve(2 * 3 * tot_size);
for(const auto& arc : subsampled_arcs)
{
for(std::size_t i=0; i<arc.size()-1; ++i)
{
pos_lines.push_back(arc[i].x());
pos_lines.push_back(arc[i].y());
pos_lines.push_back(arc[i].z());
pos_lines.push_back(arc[i+1].x());
pos_lines.push_back(arc[i+1].y());
pos_lines.push_back(arc[i+1].z());
}
}
setSceneRadius(radius_*1.3);
setSceneCenter({center_.x(), center_.y(), center_.z()});
showEntireScene();
}
void Viewer::attrib_buffers(CGAL::QGLViewer* viewer)
{
QMatrix4x4 mvpMatrix;
QMatrix4x4 mvMatrix;
double mat[16];
viewer->camera()->getModelViewProjectionMatrix(mat);
for(int i=0; i < 16; i++)
mvpMatrix.data()[i] = (float)mat[i];
viewer->camera()->getModelViewMatrix(mat);
for(int i=0; i < 16; i++)
mvMatrix.data()[i] = (float)mat[i];
// define material
QVector4D ambient(0.1f, 0.1f, 0.1f, 1.0f);
QVector4D diffuse( 0.9f,
0.9f,
0.9f,
0.0f );
QVector4D specular( 0.0f,
0.0f,
0.0f,
0.0f );
QVector4D position( -1.2f, 1.2f, .9797958971f, 1.0f );
GLfloat shininess = 1.0f;
mvpLocation = rendering_program.uniformLocation("mvp_matrix");
mvLocation = rendering_program.uniformLocation("mv_matrix");
colorLocation = rendering_program.uniformLocation("color");
lightLocation[0] = rendering_program.uniformLocation("light_pos");
lightLocation[1] = rendering_program.uniformLocation("light_diff");
lightLocation[2] = rendering_program.uniformLocation("light_spec");
lightLocation[3] = rendering_program.uniformLocation("light_amb");
lightLocation[4] = rendering_program.uniformLocation("spec_power");
rendering_program.setUniformValue(lightLocation[0], position);
rendering_program.setUniformValue(lightLocation[1], diffuse);
rendering_program.setUniformValue(lightLocation[2], specular);
rendering_program.setUniformValue(lightLocation[3], ambient);
rendering_program.setUniformValue(lightLocation[4], shininess);
rendering_program.setUniformValue(mvpLocation, mvpMatrix);
rendering_program.setUniformValue(mvLocation, mvMatrix);
rendering_program_edges.setUniformValue("mvp_matrix", mvpMatrix);
}
void Viewer::draw()
{
glEnable(GL_DEPTH_TEST);
QColor color;
rendering_program.bind();
if(draw_balls)
{
vao[SPHERE_VAO].bind();
attrib_buffers(this);
color.setRgbF(1.0f, 1.0f, 1.0f);
rendering_program.setUniformValue(colorLocation, color);
glDrawArraysInstanced(GL_TRIANGLES, 0, static_cast<GLsizei>(pos_sphere.size()/3), 1);
vao[SPHERE_VAO].release();
}
if(draw_inputs)
{
vao[POINTS_VAO].bind();
attrib_buffers(this);
color.setRgbF(0.0f, 1.0f, 0.0f);
rendering_program.setUniformValue(colorLocation, color);
glDrawArraysInstanced(GL_TRIANGLES, 0,
static_cast<GLsizei>(pos_sphere_inter.size()/3),
static_cast<GLsizei>(pos_points.size()/3));
vao[POINTS_VAO].release();
}
rendering_program.release();
rendering_program_edges.bind();
vao[EDGES_VAO].bind();
attrib_buffers(this);
glDrawArrays(GL_LINES, 0, static_cast<GLsizei>(pos_lines.size()/3));
vao[EDGES_VAO].release();
rendering_program_edges.release();
}
void Viewer::init()
{
glDrawArraysInstanced = (PFNGLDRAWARRAYSINSTANCEDARBPROC)this->context()->getProcAddress("glDrawArraysInstancedARB");
glVertexAttribDivisor = (PFNGLVERTEXATTRIBDIVISORARBPROC)this->context()->getProcAddress("glVertexAttribDivisorARB");
if(!glDrawArraysInstanced || !glDrawArraysInstanced)
{
qDebug()<<"Missing OpenGL extensions. Demo not available.";
return;
}
compile_shaders();
compute_elements();
initialize_buffers();
glEnable(GL_BLEND);
glHint(GL_LINE_SMOOTH_HINT, GL_NICEST);
glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);
glEnable(GL_POINT_SMOOTH);
glEnable(GL_LINE_SMOOTH);
}
QString Viewer::helpString() const
{
QString text("<h2>S i m p l e V i e w e r</h2>");
text += "Use the mouse to move the camera around the object. ";
text += "You can respectively revolve around, zoom and translate with the three mouse buttons. ";
text += "Left and middle buttons pressed together rotate around the camera view direction axis<br><br>";
text += "Pressing <b>Alt</b> and one of the function keys (<b>F1</b>..<b>F12</b>) defines a camera keyFrame. ";
text += "Simply press the function key again to restore it. Several keyFrames define a ";
text += "camera path. Paths are saved when you quit the application and restored at next start.<br><br>";
text += "Press <b>F</b> to display the frame rate, <b>A</b> for the world axis, ";
text += "<b>Alt+Return</b> for full screen mode and <b>Control+S</b> to save a snapshot. ";
text += "See the <b>Keyboard</b> tab in this window for a complete shortcut list.<br><br>";
text += "Double clicks automates single click actions: A left button double click aligns the closer axis with the camera (if close enough). ";
text += "A middle button double click fits the zoom of the camera and the right button re-centers the scene.<br><br>";
text += "A left button double click while holding right button pressed defines the camera <i>Revolve Around Point</i>. ";
text += "See the <b>Mouse</b> tab and the documentation web pages for details.<br><br>";
text += "Press <b>Escape</b> to exit the viewer.";
return text;
}
void Viewer::keyPressEvent(QKeyEvent *e)
{
if(e->key() == Qt::Key_P)
{
draw_inputs = !draw_inputs;
update();
return;
}
if(e->key() == Qt::Key_S)
{
draw_balls = !draw_balls;
update();
return;
}
QGLViewer::keyPressEvent(e);
}

173
Triangulation_on_sphere_2/demo/Triangulation_on_sphere_2/Viewer.h Normal file
View File

@ -0,0 +1,173 @@
#ifndef TRIANGULATION_ON_SPHERE_2_VIEWER_H
#define TRIANGULATION_ON_SPHERE_2_VIEWER_H
#include <QMap>
#include <CGAL/Qt/qglviewer.h>
#include <QOpenGLFunctions_2_1>
#include <QOpenGLVertexArrayObject>
#include <QOpenGLBuffer>
#include <QOpenGLShaderProgram>
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Delaunay_triangulation_on_sphere_traits_2.h>
#include <CGAL/Projection_on_sphere_traits_3.h>
#include <CGAL/Delaunay_triangulation_on_sphere_2.h>
#include <list>
#include <vector>
class Viewer
: public CGAL::QGLViewer
{
typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel;
typedef Kernel::FT FT;
typedef Kernel::Point_3 Point_3;
typedef Kernel::Segment_3 Segment_3;
typedef CGAL::Projection_on_sphere_traits_3<Kernel> Projection_traits;
typedef CGAL::Delaunay_triangulation_on_sphere_2<Projection_traits> Projected_DToS2;
typedef std::list<std::vector<Point_3> > Subsampled_arcs;
template <class TOS>
void build_the_boundary(const TOS& tos)
{
typedef typename TOS::Geom_traits::SK Spherical_kernel;
if(tos.dimension() < 2)
return;
typename TOS::Solid_edges_iterator it = tos.solid_edges_begin();
for (; it != tos.solid_edges_end(); ++it)
{
typename TOS::Edge e = *it;
const bool diametral_edge = CGAL::collinear(tos.construct_point(tos.point(e.first, (e.second+1)%3)),
tos.construct_point(tos.point(e.first, (e.second+2)%3)),
tos.geom_traits().center());
subsampled_arcs.emplace_front();
// primal
if(!diametral_edge)
{
typename TOS::Arc_on_sphere_2 as = tos.segment_on_sphere(e);
std::vector<typename Spherical_kernel::Point_3> discretization_points;
CGAL::Triangulations_on_sphere_2::internal::subsample_arc_on_sphere_2<Spherical_kernel>(
as, std::back_inserter(discretization_points), min_edge_size);
subsampled_arcs.begin()->reserve(discretization_points.size());
for(const typename Spherical_kernel::Point_3& spt : discretization_points)
subsampled_arcs.begin()->push_back(Point_3(spt.x(), spt.y(), spt.z()));
}
else
{
const Segment_3 s = tos.segment(e);
subsampled_arcs.begin()->push_back(s.source());
subsampled_arcs.begin()->push_back(s.target());
}
}
}
public:
Viewer(QWidget* parent = 0);
template <typename TOS, typename Iterator>
void open(Iterator begin, Iterator end,
const TOS& tos)
{
draw_balls = true;
draw_inputs = false;
subsampled_arcs.clear();
inputs.clear();
pos_points.clear();
pos_lines.clear();
pos_sphere_inter.clear();
normals_inter.clear();
normals_lines.clear();
radius_ = tos.geom_traits().radius();
center_ = tos.geom_traits().center();
min_edge_size = 0.01 * radius_;
std::copy(begin, end, std::back_inserter(inputs));
build_the_boundary(tos);
compute_elements();
initialize_buffers();
qDebug() << "Finished loading";
}
GLuint dl_nb;
protected :
virtual void draw();
virtual void init();
virtual QString helpString() const;
virtual void keyPressEvent(QKeyEvent *e);
private:
//Shaders elements
int vertexLocation[3];
int normalsLocation[3];
int centerLocation;
int trivialCenterLocation;
int mvpLocation;
int mvLocation;
int colorLocation;
int lightLocation[5];
std::vector<float> pos_points;
std::vector<float> pos_lines;
std::vector<float> pos_sphere;
std::vector<float> pos_sphere_inter;
std::vector<float> normals;
std::vector<float> normals_inter;
std::vector<float> trivial_center;
std::vector<float> normals_lines;
double min_edge_size;
Subsampled_arcs subsampled_arcs;
std::vector<Point_3> inputs;
bool draw_inputs;
bool draw_balls;
enum VBO {
SPHERE_POINTS =0,
SPHERE_NORMALS,
SPHERE_CENTER,
POINTS_POINTS,
POINTS_NORMALS,
POINTS_CENTERS,
EDGES_POINTS,
SIZE_OF_VBO
};
enum VAO {
SPHERE_VAO = 0,
POINTS_VAO,
EDGES_VAO,
SIZE_OF_VAO
};
QOpenGLBuffer buffers[SIZE_OF_VBO];
QOpenGLVertexArrayObject vao[SIZE_OF_VAO];
QOpenGLShaderProgram rendering_program;
QOpenGLShaderProgram rendering_program_edges;
typedef void (APIENTRYP PFNGLDRAWARRAYSINSTANCEDARBPROC) (GLenum mode, GLint first, GLsizei count, GLsizei primcount);
typedef void (APIENTRYP PFNGLVERTEXATTRIBDIVISORARBPROC) (GLuint index, GLuint divisor);
PFNGLDRAWARRAYSINSTANCEDARBPROC glDrawArraysInstanced;
PFNGLVERTEXATTRIBDIVISORARBPROC glVertexAttribDivisor;
double radius_;
Point_3 center_;
void initialize_buffers();
void compute_elements();
void attrib_buffers(CGAL::QGLViewer*);
void compile_shaders();
};
#endif // TRIANGULATION_ON_SPHERE_2_VIEWER_H

355
Triangulation_on_sphere_2/demo/Triangulation_on_sphere_2/create_sphere.h Normal file
View File

@ -0,0 +1,355 @@
#ifndef POLYHEDRON_DEMO_CREATE_SPHERE_H
#define POLYHEDRON_DEMO_CREATE_SPHERE_H
#include <CGAL/number_type_config.h>
#include <cmath>
#include <vector>
template <class FLOAT>
void create_flat_sphere(FLOAT R,
std::vector<FLOAT>& positions_spheres,
std::vector<FLOAT>& normals_spheres,
int prec = 36)
{
// The smaller they are, the more precise the sphere will be
// Must be divisors of 180 and 360
const float rings = float(prec) / 2;
const float sectors = float(prec);
const float to_rad = static_cast<float>(CGAL_PI / 180.);
float T, P;
float x[4],y[4],z[4];
//Top of the sphere
for(float t=0; t<360.f; t+=sectors)
{
positions_spheres.push_back(0);
positions_spheres.push_back(0);
positions_spheres.push_back(R);
normals_spheres.push_back(0);
normals_spheres.push_back(0);
normals_spheres.push_back(1);
P = rings*to_rad;
T = t*to_rad;
x[1] = sin(P) * cos(T);
y[1] = sin(P) * sin(T);
z[1] = cos(P);
positions_spheres.push_back(R * x[1]);
positions_spheres.push_back(R * y[1]);
positions_spheres.push_back(R * z[1]);
normals_spheres.push_back(x[1]);
normals_spheres.push_back(y[1]);
normals_spheres.push_back(z[1]);
//
P = rings*to_rad;
T = (t+sectors)*to_rad;
x[2] = sin(P) * cos(T);
y[2] = sin(P) * sin(T);
z[2] = cos(P);
positions_spheres.push_back(R * x[2]);
positions_spheres.push_back(R * y[2]);
positions_spheres.push_back(R * z[2]);
normals_spheres.push_back(x[2]);
normals_spheres.push_back(y[2]);
normals_spheres.push_back(z[2]);
}
//Body of the sphere
for(float p=rings; p<180.f-rings; p+=rings)
{
for(float t=0; t<360.f; t+=sectors)
{
//A
P = p*to_rad;
T = t*to_rad;
x[0] = sin(P) * cos(T);
y[0] = sin(P) * sin(T);
z[0] = cos(P);
positions_spheres.push_back(R * x[0]);
positions_spheres.push_back(R * y[0]);
positions_spheres.push_back(R * z[0]);
normals_spheres.push_back(x[0]);
normals_spheres.push_back(y[0]);
normals_spheres.push_back(z[0]);
//B
P = (p+rings)*to_rad;
T = t*to_rad;
x[1] = sin(P) * cos(T);
y[1] = sin(P) * sin(T);
z[1] = cos(P);
positions_spheres.push_back(R * x[1]);
positions_spheres.push_back(R * y[1]);
positions_spheres.push_back(R * z[1]);
normals_spheres.push_back(x[1]);
normals_spheres.push_back(y[1]);
normals_spheres.push_back(z[1]);
//C
P = p*to_rad;
T = (t+sectors)*to_rad;
x[2] = sin(P) * cos(T);
y[2] = sin(P) * sin(T);
z[2] = cos(P);
positions_spheres.push_back(R * x[2]);
positions_spheres.push_back(R * y[2]);
positions_spheres.push_back(R * z[2]);
normals_spheres.push_back(x[2]);
normals_spheres.push_back(y[2]);
normals_spheres.push_back(z[2]);
//D
P = (p+rings)*to_rad;
T = (t+sectors)*to_rad;
x[3] = sin(P) * cos(T);
y[3] = sin(P) * sin(T);
z[3] = cos(P);
positions_spheres.push_back(R * x[3]);
positions_spheres.push_back(R * y[3]);
positions_spheres.push_back(R * z[3]);
normals_spheres.push_back(x[3]);
normals_spheres.push_back(y[3]);
normals_spheres.push_back(z[3]);
positions_spheres.push_back(R * x[1]);
positions_spheres.push_back(R * y[1]);
positions_spheres.push_back(R * z[1]);
normals_spheres.push_back(x[1]);
normals_spheres.push_back(y[1]);
normals_spheres.push_back(z[1]);
positions_spheres.push_back(R * x[2]);
positions_spheres.push_back(R * y[2]);
positions_spheres.push_back(R * z[2]);
normals_spheres.push_back(x[2]);
normals_spheres.push_back(y[2]);
normals_spheres.push_back(z[2]);
}
}
//Bottom of the sphere
for(float t=0; t<360.f; t+=sectors)
{
positions_spheres.push_back(0);
positions_spheres.push_back(0);
positions_spheres.push_back(-R);
normals_spheres.push_back(0);
normals_spheres.push_back(0);
normals_spheres.push_back(-1);
P = (180-rings)*to_rad;
T = t*to_rad;
x[1] = sin(P) * cos(T);
y[1] = sin(P) * sin(T);
z[1] = cos(P);
positions_spheres.push_back(R * x[1]);
positions_spheres.push_back(R * y[1]);
positions_spheres.push_back(R * z[1]);
normals_spheres.push_back(x[1]);
normals_spheres.push_back(y[1]);
normals_spheres.push_back(z[1]);
P = (180-rings)*to_rad;
T = (t+sectors)*to_rad;
x[2] = sin(P) * cos(T);
y[2] = sin(P) * sin(T);
z[2] = cos(P);
positions_spheres.push_back(R * x[2]);
positions_spheres.push_back(R * y[2]);
positions_spheres.push_back(R * z[2]);
normals_spheres.push_back(x[2]);
normals_spheres.push_back(y[2]);
normals_spheres.push_back(z[2]);
}
}
template <class FLOAT>
void create_flat_and_wire_sphere(FLOAT R,
std::vector<FLOAT>& positions_spheres,
std::vector<FLOAT>& normals_spheres,
std::vector<FLOAT>& positions_wire_spheres,
int prec = 36)
{
// The smaller they are, the more precise the sphere will be
// Must be a divisor of 360 and 180.
const float rings = float(prec)/2;
const float sectors = float(prec);
const float to_rad = static_cast<float>(CGAL_PI / 180.0);
create_flat_sphere(R, positions_spheres, normals_spheres);
float T, P;
float x[4],y[4],z[4];
//Top of the sphere
for(float t=0; t<360; t+=sectors)
{
positions_wire_spheres.push_back(0);
positions_wire_spheres.push_back(0);
positions_wire_spheres.push_back(R);
P = rings*to_rad;
T = t*to_rad;
x[1] = sin(P) * cos(T);
y[1] = sin(P) * sin(T);
z[1] = cos(P);
positions_wire_spheres.push_back(R * x[1]);
positions_wire_spheres.push_back(R * y[1]);
positions_wire_spheres.push_back(R * z[1]);
positions_wire_spheres.push_back(R * x[1]);
positions_wire_spheres.push_back(R * y[1]);
positions_wire_spheres.push_back(R * z[1]);
//
P = rings*to_rad;
T = (t+sectors)*to_rad;
x[2] = sin(P) * cos(T);
y[2] = sin(P) * sin(T);
z[2] = cos(P);
positions_wire_spheres.push_back(R * x[2]);
positions_wire_spheres.push_back(R * y[2]);
positions_wire_spheres.push_back(R * z[2]);
positions_wire_spheres.push_back(R * x[2]);
positions_wire_spheres.push_back(R * y[2]);
positions_wire_spheres.push_back(R * z[2]);
positions_wire_spheres.push_back(0);
positions_wire_spheres.push_back(0);
positions_wire_spheres.push_back(R);
}
//Body of the sphere
for(float p=rings; p<180.f-rings; p+=rings)
{
for(float t=0; t<360.f; t+=sectors)
{
//A
P = p*to_rad;
T = t*to_rad;
x[0] = sin(P) * cos(T);
y[0] = sin(P) * sin(T);
z[0] = cos(P);
//B
P = (p+rings)*to_rad;
T = t*to_rad;
x[1] = sin(P) * cos(T);
y[1] = sin(P) * sin(T);
z[1] = cos(P);
//C
P = p*to_rad;
T = (t+sectors)*to_rad;
x[2] = sin(P) * cos(T);
y[2] = sin(P) * sin(T);
z[2] = cos(P);
//D
P = (p+rings)*to_rad;
T = (t+sectors)*to_rad;
x[3] = sin(P) * cos(T);
y[3] = sin(P) * sin(T);
z[3] = cos(P);
positions_wire_spheres.push_back(R * x[0]);
positions_wire_spheres.push_back(R * y[0]);
positions_wire_spheres.push_back(R * z[0]);
positions_wire_spheres.push_back(R * x[1]);
positions_wire_spheres.push_back(R * y[1]);
positions_wire_spheres.push_back(R * z[1]);
positions_wire_spheres.push_back(R * x[1]);
positions_wire_spheres.push_back(R * y[1]);
positions_wire_spheres.push_back(R * z[1]);
positions_wire_spheres.push_back(R * x[3]);
positions_wire_spheres.push_back(R * y[3]);
positions_wire_spheres.push_back(R * z[3]);
positions_wire_spheres.push_back(R * x[3]);
positions_wire_spheres.push_back(R * y[3]);
positions_wire_spheres.push_back(R * z[3]);
positions_wire_spheres.push_back(R * x[2]);
positions_wire_spheres.push_back(R * y[2]);
positions_wire_spheres.push_back(R * z[2]);
positions_wire_spheres.push_back(R * x[2]);
positions_wire_spheres.push_back(R * y[2]);
positions_wire_spheres.push_back(R * z[2]);
positions_wire_spheres.push_back(R * x[0]);
positions_wire_spheres.push_back(R * y[0]);
positions_wire_spheres.push_back(R * z[0]);
}
}
//Bottom of the sphere
for(float t=0; t<360.f; t+=sectors)
{
P = (180-rings)*to_rad;
T = t*to_rad;
x[1] = sin(P) * cos(T);
y[1] = sin(P) * sin(T);
z[1] = cos(P);
P = (180-rings)*to_rad;
T = (t+sectors)*to_rad;
x[2] = sin(P) * cos(T);
y[2] = sin(P) * sin(T);
z[2] = cos(P);
positions_wire_spheres.push_back(0);
positions_wire_spheres.push_back(0);
positions_wire_spheres.push_back(-R);
positions_wire_spheres.push_back(R * x[1]);
positions_wire_spheres.push_back(R * y[1]);
positions_wire_spheres.push_back(R * z[1]);
positions_wire_spheres.push_back(R * x[1]);
positions_wire_spheres.push_back(R * y[1]);
positions_wire_spheres.push_back(R * z[1]);
positions_wire_spheres.push_back(R * x[2]);
positions_wire_spheres.push_back(R * y[2]);
positions_wire_spheres.push_back(R * z[2]);
positions_wire_spheres.push_back(R * x[2]);
positions_wire_spheres.push_back(R * y[2]);
positions_wire_spheres.push_back(R * z[2]);
positions_wire_spheres.push_back(0);
positions_wire_spheres.push_back(0);
positions_wire_spheres.push_back(-R);
}
}
#endif //POLYHEDRON_DEMO_DRAW_SPHERE_H

127
Triangulation_on_sphere_2/demo/Triangulation_on_sphere_2/main.cpp Normal file
View File

@ -0,0 +1,127 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Delaunay_triangulation_on_sphere_traits_2.h>
#include <CGAL/Projection_on_sphere_traits_3.h>
#include <CGAL/Delaunay_triangulation_on_sphere_2.h>
#include <CGAL/utility.h>
#include <CGAL/Qt/DemosMainWindow.h>
#include <QApplication>
#include <QMainWindow>
#include <QFileDialog>
#include <QInputDialog>
#include <CGAL/Three/Three.h>
#include <boost/iterator/transform_iterator.hpp>
#include <list>
#include <fstream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel;
typedef Kernel::FT FT;
typedef Kernel::Point_3 Point_3;
typedef Kernel::Segment_3 Segment_3;
typedef CGAL::Projection_on_sphere_traits_3<Kernel> Projection_traits;
typedef CGAL::Delaunay_triangulation_on_sphere_2<Projection_traits> Projected_DToS2;
typedef std::list<std::vector<Point_3> > Subsampled_arcs;
#include "Viewer.h"
#include "ui_Mainwindow.h"
template <class Output_iterator>
void read_points(const char* file_path, Output_iterator out)
{
std::ifstream input(file_path);
if(!input)
{
std::cerr << "Error while reading " << file_path << std::endl;
std::exit(EXIT_FAILURE);
}
double x,y,z;
while(input >> x >> y >> z)
*out++ = Point_3(x, y, z);
}
class MainWindow
: public CGAL::Qt::DemosMainWindow,
public Ui::MainWindow
{
Q_OBJECT
public:
MainWindow()
{
setupUi(this);
}
public slots:
void open(QString filename)
{
std::vector<Point_3> lst_pt;
read_points(filename.toUtf8().data(), std::back_inserter(lst_pt));
const Point_3 center(0,0,0);
const QString infos = QInputDialog::getText(nullptr, "Sphere", "Enter the sphere's information : (Radius center.x center.y center.z):",
QLineEdit::Normal,"100.0 0.0 0.0 0.0");
QStringList list = infos.split(QRegExp("\\s+"), CGAL_QT_SKIP_EMPTY_PARTS);
if (list.isEmpty()) return;
if (list.size()!=4){
QMessageBox *msgBox = new QMessageBox;
msgBox->setWindowTitle("Error");
msgBox->setText("ERROR : Input should consists of 4 doubles: The radius first, then the coordinates of the center.");
msgBox->exec();
return;
}
double coords[4];
for(int j=0; j<4; ++j)
{
bool ok;
coords[j] = list.at(j).toDouble(&ok);
if(!ok)
{
QMessageBox *msgBox = new QMessageBox;
msgBox->setWindowTitle("Error");
msgBox->setText("ERROR : Input is invalid.");
msgBox->exec();
return;
}
}
Projection_traits traits(Point_3(coords[1], coords[2], coords[3]), coords[0]);
Projected_DToS2 dtos(lst_pt.begin(), lst_pt.end(), traits);
std::cout << dtos.number_of_vertices() << " vertices" << std::endl;
// Instantiate the viewer
viewer->open(lst_pt.begin(), lst_pt.end(), dtos);
}
void on_action_Quit_triggered()
{
close();
}
void on_action_Open_triggered()
{
QString filename = QFileDialog::getOpenFileName(this);
if(!filename.isNull())
open(filename);
}
};
int main(int argc, char** argv)
{
// Read command lines arguments
QApplication application(argc,argv);
MainWindow mainWindow;
mainWindow.show();
// Run main loop
return application.exec();
}
#include "main.moc"

260
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/CGAL/Delaunay_triangulation_on_sphere_2.h Normal file
View File

@ -0,0 +1,260 @@
namespace CGAL {
/*!
\ingroup PkgTriangulationOnSphere2TriangulationClasses
The class `Delaunay_triangulation_on_sphere_2` is designed to represent
the Delaunay triangulation of a set of points on the 2-sphere.
A Delaunay triangulation of a set of points is a triangulation of the set of points that fulfills
the following <I>empty circle property</I> (also called <I>Delaunay property</I>): the circumscribing
sphere of any simplex of the triangulation contains no point of the set in its interior.
For a point set with no case of co-circularity of more than three points,
the Delaunay triangulation is unique, and defined as the dual of the Voronoi diagram of the points.
The setting of 3D points on the 2-sphere is particular in that the empty circle property
can be reduced to a single 3D orientation test \cgalCite{cgal:ccplr-redtp-10}.
\tparam Traits is the geometric traits; it must be a model of `DelaunayTriangulationOnSphereTraits_2`.
\tparam TDS is the triangulation data structure, which must be a model of `TriangulationDataStructure_2`
whose vertex base must be a model of `TriangulationOnSphereVertexBase_2` and whose face base
must be a model of `TriangulationOnSphereFaceBase_2`. \cgal provides a default instantiation
for this parameter, which is the class
`CGAL::Triangulation_data_structure_2 < CGAL::Triangulation_on_sphere_vertex_base_2<Traits>,
CGAL::Triangulation_on_sphere_face_base_2<Traits> >`.
\sa `CGAL::Triangulation_on_sphere_2<Traits, TDS>`
*/
template< typename Traits, typename TDS >
class Delaunay_triangulation_on_sphere_2
: public Triangulation_on_sphere_2<Traits, TDS>
{
public:
/// \name Types
/// @{
/*!
The traits class.
*/
typedef Traits Geom_traits;
/*!
The point type representing a point on the sphere.
*/
typedef Traits::Point_on_sphere_2 Point;
/// @}
/// \name Creation
/// @{
/*!
introduces an empty triangulation and sets the center and radius of the sphere to `c` and `r` respectively.
*/
Delaunay_triangulation_on_sphere_2(const Point_3& c, const FT r);
/*!
introduces an empty triangulation.
\note The values for the center and radius must be either already set in the traits
(if `gt` is provided) or must be set after the construction, using the function `set_center_and_radius()`.
*/
Delaunay_triangulation_on_sphere_2(const Geom_traits& gt = Geom_traits());
/*!
introduces an empty triangulation, sets the center and radius of the sphere to `c` and `r` respectively,
and inserts the point range `[first; beyond[`.
\tparam PointOnSphereIterator must be a model of `InputIterator` with value type `Point_on_sphere_2` or `Point_3`.
*/
template <typename PointOnSphereIterator>
Delaunay_triangulation_on_sphere_2(PointOnSphereIterator first, PointOnSphereIterator beyond,
const Point_3& center, const FT radius);
/*!
introduces an empty triangulation whose center and radius are set according to values within the traits
and inserts the point range `[first;beyond[`.
\warning It is the user's responsability to ensure that the center and radius are set as intended in `gt`.
\tparam PointOnSphereIterator must be a model of `InputIterator` with value type `Point_on_sphere_2` or `Point_3`.
*/
template <typename PointOnSphereIterator>
Delaunay_triangulation_on_sphere_2(PointOnSphereIterator first, PointOnSphereIterator beyond,
const Geom_traits& gt = Geom_traits());
/*!
Copy constructor. All the vertices and faces are duplicated.
*/
Delaunay_triangulation_on_sphere_2(const Delaunay_triangulation_on_sphere_2<Traits,TDS> &tr);
/// @}
public:
/// \name Predicates
///
/// @{
/*!
returns the side of `p` with respect to the circle circumscribing the triangle associated with `f`.
*/
Oriented_side side_of_oriented_circle(Face_handle f, const Point& p) const;
/// @}
public:
/// \name Insertion and Removal
///
/// Methods for the insertion and removal of points on the sphere.
/// In the degenerate setting of co-circular points, the Delaunay triangulation is known
/// not to be uniquely defined. In this case, \cgal chooses a particular Delaunay triangulation
/// using a symbolic perturbation scheme \cgalCite{cgal:dt-pvr3d-03}.
///
/// @{
/*!
inserts the point `p`.
If the point `p` coincides with an already existing vertex, this vertex is returned
and the triangulation remains unchanged.
The optional parameter `f` is used to give a hint about the location of `p`.
*/
Vertex_handle insert(const Point& p, Face_handle f = Face_handle());
/*!
inserts the point `p` at the location given by `(lt, loc, li)`.
\sa `Triangulation_on_sphere_2::locate()`
*/
Vertex_handle insert(const Point& p, Locate_type& lt, Face_handle loc, int li );
/*!
Equivalent to `insert(p)`.
*/
Vertex_handle push_back(const Point& p);
/*!
inserts the points in the range `[first, beyond)` and returns the number of inserted points.
\tparam PointOnSphereIterator must be a model of `InputIterator` with value type `Point`.
*/
template <class PointOnSphereIterator>
size_type insert(PointOnSphereIterator first, PointOnSphereIterator beyond);
/*!
removes the vertex `v` from the triangulation.
*/
void remove(Vertex_handle v);
/// @}
/// \name Queries
/// @{
/*!
outputs the faces and boundary edges of the conflict zone of point `p` into output iterators.
This function outputs in the container pointed to by `fit` the faces which are in conflict with point `p`,
i. e., the faces whose circumcircle contains `p`. It outputs in the container pointed to by `eit`
the boundary of the zone in conflict with `p`. The boundary edges of the conflict zone
are output in counter-clockwise order and each edge is described through its incident face
which is not in conflict with `p`. The function returns in a `std::pair` the resulting output iterators.
\tparam OutItFaces is an output iterator with `Face_handle` as value type.
\tparam OutItBoundaryEdges is an output iterator with `Edge` as value type.
\warning This function makes uses of the member `tds_data` (see the concept `TriangulationDSFaceBase_2`)
of the face to mark faces in conflict. It is the responsability of the user to make sure the flags are cleared.
\pre `dimension() == 2`.
*/
template <typename OutputItFaces, typename OutputItBoundaryEdges>
std::pair<OutputItFaces, OutputItBoundaryEdges>
get_conflicts_and_boundary(const Point& p,
OutputItFaces fit,
OutputItBoundaryEdges eit,
Face_handle start) const;
/// @}
/// \name Voronoi Diagram
///
/// The following member functions provide the elements of the dual Voronoi diagram.
/// Two different embeddings are possible: a "straight" embedding, using line segments living
/// in Euclidean 3D sphere, and a "curved" embedding, using arc segments on the sphere.
///
/// Note that the following operations are constructions, which should be kept in mind in the choice
/// of the underlying kernel.
///
/// @{
// Straight
/*!
returns the center of the circle circumscribed to face `f`
*/
Point_3 dual(const Face_handle f) const;
/*!
returns the line segment with endpoints the circumcenters of the faces incident to the edge `e`.
*/
Segment_3 dual(const Edge& e) const;
/*!
returns the line segment with endpoints the circumcenters of the faces incident to the edge `*ec`.
*/
Segment_3 dual(const Edge_circulator ec) const;
/*!
returns the line segment with endpoints the circumcenters of the faces incident to the edge `*ei`.
*/
Segment_3 dual(const All_edges_iterator ei) const;
// Curved
/*!
returns the intersection of the dual of the face `f` and the sphere.
\pre `dimension() == 2` and `f` is a solid face.
*/
Point dual_on_sphere(const Face_handle f) const;
/*!
returns the arc of great circle with endpoints the circumcenters of the faces incident to the edge `e`.
\pre `dimension() == 2` and `e` is not a ghost edge.
*/
Arc_on_sphere_2 dual_on_sphere(const Edge& e) const;
/*!
returns the arc of great circle with endpoints the circumcenters of the faces incident to the edge `*ec`.
\pre `dimension() == 2` and `*ec` is not a ghost edge.
*/
Arc_on_sphere_2 dual_on_sphere(const Edge_circulator ec) const;
/*!
returns the arc of great circle with endpoints the circumcenters of the faces incident to the edge `*ei`.
\pre `dimension() == 2` and `*ei` is not a ghost edge.
*/
Arc_on_sphere_2 dual_on_sphere(const All_edges_iterator ei) const;
/// @}
/// \name Miscellaneous
///
/// @{
/*!
tests the validity of the triangulation as a `Triangulation_on_sphere_2`
and additionally tests the Delaunay property. This method is mainly useful
for debugging Delaunay triangulation algorithms.
*/
bool is_valid(bool verbose = false, int level = 0) const;
/// @}
};
} // namespace CGAL

119
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/CGAL/Delaunay_triangulation_on_sphere_traits_2.h Normal file
View File

@ -0,0 +1,119 @@
namespace CGAL {
/*!
\ingroup PkgTriangulationOnSphere2TriangulationClasses
The class `Delaunay_triangulation_on_sphere_traits_2` is a model
of the concept `DelaunayTriangulationOnSphereTraits_2`.
The `Point_on_sphere_2` type is implemented as a kernel `Point_3` type.
If the kernel template parameter `LK` does not enable exact representation of points on sphere
(i.e., at least a mean to represent algebraic coordinates), then it cannot be guaranteed
that all points on the sphere are in a convex position and thus some points might be hidden upon insertion.
It is possible to ensure that no point can be hidden by enforcing a tiny gap between points \cgalCite{cgal:ccplr-redtp-10}.
In Lemma 4.1 of this publication, it is in particular proven that if points lie within a distance
\f$ \delta \f$ of the sphere, no point is hidden as long as points are separated by at least \f$ 2 \sqrt{R\delta} \f$
with \f$ R \f$ the radius of the sphere.
Thus, if `LK` offers exact representation, then \f$ \delta \f$ is set to \f$ 0 \f$.
Otherwise, \f$ \delta \f$ is set to the maximal distance between two consecutive `LK::FT`
for the coordinates of a point that is (theoretically) on the sphere.
This bound \f$ \delta \f$ is then used in the functions `is_on_sphere()` and `are_points_too_close()`
using the relation above to ensure that a point being inserted is either marked as "too close"
(see \link CGAL::Triangulation_on_sphere_2::Locate_type `CGAL::Triangulation_on_sphere_2<Traits,TDS>::Locate_type` \endlink)
and thus not inserted, or guaranteed to not be hidden upon insertion.
\tparam LK a linear kernel type; it must be a model of `Kernel`.
\tparam SK a spherical kernel type; it must be a model of `SphericalKernel` refining `LK`.
\cgalModels `DelaunayTriangulationOnSphereTraits_2`
\sa `CGAL::Projection_on_sphere_traits_3`
*/
template <typename LK,
typename SK = CGAL::Spherical_kernel_3<
LK, CGAL::Algebraic_kernel_for_spheres_2_3<typename LK::FT> > >
class Delaunay_triangulation_on_sphere_traits_2
{
public:
/// The field number type
typedef typename LK::FT FT;
///
typedef typename LK::Point_3 Point_on_sphere_2;
/// An arc of a great circle, used to represent a curved segment on the sphere (Voronoi or Delaunay edge).
typedef typename SK::Circular_arc_3 Arc_on_sphere_2;
///
typedef typename LK::Point_3 Point_3;
///
typedef typename LK::Segment_3 Segment_3;
///
typedef typename LK::Triangle_3 Triangle_3;
/// \name Predicates
///
/// @{
/// Internally uses `LK::Coplanar_orientation_3`
typedef unspecified_type Collinear_are_strictly_ordered_on_great_circle_2;
///
typedef typename LK::Compare_xyz_3 Compare_on_sphere_2;
/// If the kernel cannot represent algebraic coordinates exactly, there is a tolerance
/// around the sphere, and thus different points of \f$ \mathbb{R}^3\f$ can actually
/// correspond to the same point on the sphere.
/// This functor checks if two points project onto the same point on the sphere.
typedef unspecified_type Equal_on_sphere_2;
///
typedef typename LK::Orientation_3 Side_of_oriented_circle_on_sphere_2;
/// Internally uses `LK::Orientation_3`
typedef unspecified_type Orientation_on_sphere_2;
/// @}
/// \name Constructions
///
/// @{
/// Internally uses `SK::Construct_circular_arc_3`
typedef typename unspecified_type Construct_arc_on_sphere_2;
///
typedef typename LK::Construct_circumcenter_3 Construct_circumcenter_3;
///
typedef typename LK::Construct_point_3 Construct_point_3;
///
typedef typename LK::Construct_segment_3 Construct_segment_3;
///
typedef typename LK::Construct_triangle_3 Construct_triangle_3;
/// @}
public:
/// \name Precision predicates
///
/// @{
/// returns whether `p` is exactly on the sphere if `LK` can represent algebraic coordinates,
/// or whether `p` is within an automatically computed small distance otherwise.
bool is_on_sphere(const Point_on_sphere_2& p) const;
/// returns `false` if `LK` can represent algebraic coordinates, or whether the distance
/// between `p` and `q` is lower than \f$ 2 \sqrt{R\delta} \f$ otherwise.
bool are_points_too_close(const Point_on_sphere_2& p, const Point_on_sphere_2& q) const;
/// @}
};
} // namespace CGAL

69
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/CGAL/Projection_on_sphere_traits_3.h Normal file
View File

@ -0,0 +1,69 @@
namespace CGAL {
/*!
\ingroup PkgTriangulationOnSphere2TriangulationClasses
The class `Projection_on_sphere_traits_3` is a model of the concept `DelaunayTriangulationOnSphereTraits_2`.
It implements the `Point_on_sphere_2` type as a custom point type which represents the projection
of a point living in the 3D Euclidean space onto the sphere along the segment between said point
and the center of the sphere.
\tparam LK a linear kernel type; must be a model of `Kernel`.
\tparam SK a spherical kernel type; must be a model of `SphericalKernel`.
\cgalModels `DelaunayTriangulationOnSphereTraits_2`
\sa `CGAL::Delaunay_triangulation_on_sphere_traits_2`
*/
template <typename LK,
typename SK = CGAL::Spherical_kernel_3<
LK, CGAL::Algebraic_kernel_for_spheres_2_3<typename LK::FT> > >
class Projection_on_sphere_traits_3
{
public:
/// The field number type
typedef typename LK::FT FT;
/// The point on the sphere type
typedef unspecified_type Point_on_sphere_2;
///
typedef typename SK::Circular_arc_3 Arc_on_sphere_2;
///
typedef typename LK::Point_3 Point_3;
///
typedef typename LK::Segment_3 Segment_3;
///
typedef typename LK::Triangle_3 Triangle_3;
/// \name Predicates
///
/// @{
/// Points are equal if they have the same projection onto the sphere.
typedef unspecified_type Equal_on_sphere_2;
/// @}
public:
/// \name Precision predicates
///
/// @{
/// Due to the chosen point representation, any point is theoretically on the sphere, and this function
/// always returns `true`.
bool is_on_sphere(const Point_on_sphere_2& p) const;
/// returns `false` if `LK` can represent algebraic coordinates, or whether the distance
/// between `p` and `q` is lower than \f$ 2 \sqrt{R\delta} \f$ otherwise (see the traits class
/// `CGAL::Delaunay_triangulation_on_sphere_traits_2`).
bool are_points_too_close(const Point_on_sphere_2& p, const Point_on_sphere_2& q) const;
/// @}
};
} // namespace CGAL

634
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/CGAL/Triangulation_on_sphere_2.h Normal file
View File

@ -0,0 +1,634 @@
namespace CGAL {
/*!
\ingroup PkgTriangulationOnSphere2TriangulationClasses
The class `Triangulation_on_sphere_2` is the basic class designed to represent a triangulation
of a point set on a sphere: its vertices coincide with the points of the set.
\warning This triangulation supports neither the insertion nor the removal of vertices,
see `CGAL::Delaunay_triangulation_on_sphere_2` for such purposes.
This triangulation class is very similar to `CGAL::Triangulation_2` as both classes represent
triangulations of a 2-manifold domain without boundary. A significant difference is that
in the case of Euclidean 2D triangulation, it is necessary to introduce so-called <i>infinite
faces</i> to complete the convex hull into an actual 2-manifold without boundary that the triangulation
data structure can represent. This is not necessary for triangulations on the sphere,
which are already perfectly adapted to the triangulation data structure.
There is an exception to the previous statement: in the degenerate configuration
where all points lie on the same hemisphere, the triangulation has a border.
Internally, the triangulation data structure must however remain a 2-manifold at all time,
and to ensure this property, fictitious faces called <i>ghost faces</i> are added. We call faces that
are not ghost faces <em>solid faces</em>, and edges of such faces <em>solid edges</em>.
\tparam Traits is the geometric traits; it must be a model of the concept `TriangulationOnSphereTraits_2`.
\tparam TDS is the triangulation data structure; it must be a model of the concept `TriangulationDataStructure_2`,
whose vertex base must be a model of `TriangulationOnSphereVertexBase_2` and whose face base
must be a model of `TriangulationOnSphereFaceBase_2`. By default, the triangulation data structure
is instantiated by `Triangulation_data_structure_2<Triangulation_on_sphere_vertex_base_2<Gt>, Triangulation_on_sphere_face_base_2<Gt> >`.
\sa `CGAL::Delaunay_triangulation_on_sphere_2<Traits, TDS>`
*/
template <typename Traits, typename TDS>
class Triangulation_on_sphere_2
: public Triangulation_cw_ccw_2
{
public:
/// \name Types
/// @{
/*!
The traits class.
*/
typedef Traits Geom_traits;
/*!
The triangulation data structure type.
*/
typedef TDS Triangulation_data_structure;
/*!
Size type (an unsigned integral type).
*/
typedef Triangulation_data_structure::size_type size_type;
/*!
The number type.
*/
typedef Traits::FT FT;
/*!
The point type representing a point on the sphere.
*/
typedef Traits::Point_on_sphere_2 Point;
/*!
The 3D point type.
*/
typedef Traits::Point_3 Point_3;
/*!
The 3D segment type.
*/
typedef Traits::Segment_3 Segment_3;
/*!
The 3D triangle type.
*/
typedef Traits::Triangle_3 Triangle_3;
/*!
An arc of a great circle, used to represent a curved segment (Voronoi or Delaunay edge).
*/
typedef Traits::Arc_on_sphere_2 Arc_on_sphere_2;
public:
/*!
The vertex type.
*/
typedef TDS::Vertex Vertex;
/*!
The edge type.
*/
typedef TDS::Edge Edge;
/*!
The face type.
*/
typedef TDS::Face Face;
public:
/// \name Handles, Iterators, and Circulators
///
/// The vertices and faces of the triangulation are accessed through
/// handles, iterators and circulators. The handles are models
/// of the concept `Handle` which basically offers the two dereference
/// operators `*` and `->`. The handles are also model of the concepts
/// `LessThanComparable` and `Hashable`, that is they can be used as keys
/// in containers such as `std::map` and `std::unordered_map`.
/// The iterators and circulators are all bidirectional and non-mutable.
/// The circulators and iterators are convertible to handles with the same value type,
/// so that whenever a handle appear in the parameter list of a function,
/// an appropriate iterator or circulator can be passed as well.
///
/// The edges of the triangulation can also be visited through
/// iterators and circulators, the edge circulators and iterators are
/// also bidirectional and non mutable.
///
/*!
Handle to a vertex.
*/
typedef TDS::Vertex_handle Vertex_handle;
/*!
Handle to a face.
*/
typedef TDS::Face_handle Face_handle;
/*!
Iterator over all vertices.
*/
typedef TDS::Vertex_iterator Vertices_iterator;
/*!
Range type for iterating over all vertices, with a nested
type `iterator` that has as value type `Vertex_handle`.
*/
typedef Iterator_range<unspecified_type> Vertex_handles;
/*!
Iterator over all edges.
*/
typedef TDS::Edge_iterator All_edges_iterator;
/*!
Range type for iterating over all edges (including non-solid ones).
*/
typedef Iterator_range<All_edges_iterator> All_edges;
/*!
Iterator over all faces.
*/
typedef TDS::Face_iterator All_faces_iterator;
/*!
Range type for iterating over all faces (including ghost faces), with a nested
type `iterator` that has as value type `Face_handle`.
*/
typedef Iterator_range<All_faces_iterator> All_face_handles;
/*!
Iterator over all solid edges.
*/
typedef unspecified_type Solid_edges_iterator;
/*!
Range type for iterating over all solid edges.
*/
typedef Iterator_range<Solid_edges_iterator> Solid_edges;
/*!
Iterator over all solid faces.
*/
typedef unspecified_type Solid_faces_iterator;
/*!
Range type for iterating over solid faces, with a nested
type `iterator` that has as value type `Face_handle`.
*/
typedef Iterator_range<unspecified_type> Solid_face_handles;
/*!
Circulator over all vertices incident to a given vertex.
*/
typedef unspecified_type Vertex_circulator;
/*!
Circulator over all edges incident to a given vertex.
*/
typedef unspecified_type Edge_circulator;
/*!
Circulator over all faces incident to a given vertex.
*/
typedef unspecified_type Face_circulator;
/*!
Iterator over the points corresponding the vertices of the triangulation.
*/
typedef unspecified_type Point_iterator;
/*!
Range type for iterating over the points of the finite vertices.
*/
typedef Iterator_range<Point_iterator> Points;
/// @}
public:
/// \name Creation
/// @{
/*!
constructs an empty triangulation.
\note The values for the center and radius must be either already set in the traits
(if `gt` is passed) or must be set after the construction, using the function `set_center_and_radius()`.
*/
Triangulation_on_sphere_2(const Traits& gt = Traits());
/*!
constructs an empty triangulation and sets the center and radius to `c` and `r` respectively.
*/
Triangulation_on_sphere_2(const Point_3& c, const FT r);
/*!
Copy constructor. All the vertices and faces are duplicated.
After the copy, `*this` and `tr` refer to different triangulations:
if `tr` is modified, `*this` is not.
*/
Triangulation_on_sphere_2(const Triangulation_on_sphere_2& tr);
/*!
Assignment operator. This performs a deep copy of the triangulation, duplicating both vertices and faces.
*/
Triangulation_on_sphere_2& operator=(Triangulation_on_sphere_2<Traits,TDS> tr);
/*!
The triangulations `tr` and `*this` are swapped.
*/
void swap(Triangulation_on_sphere_2& tr);
/*!
deletes all faces and vertices, resulting in an empty triangulation.
*/
void clear();
/// @}
public:
/// \name Ghost Predicates
/// @{
/*!
returns `true` if `f` is a ghost face, and `false` otherwise.
\pre `dimension() == 2`
*/
bool is_ghost(const Face_handle f) const;
/*!
returns `true` if `e` is a ghost edge, that is if both its incident faces are ghost faces,
and `false` otherwise.
\pre `dimension() == 2`
*/
bool is_ghost(const Edge& e) const;
/// @}
public:
/// \name Modifying the domain
///
/// The following functions can be used to modify the geometry of the spherical domain.
///
/// \warning The triangulation is cleared in the process.
///
/// @{
/*! */
void set_center(const Point_3& c);
/*! */
void set_radius(const FT radius);
/*! */
void set_center_and_radius(const Point_3& c, const FT radius);
/// @}
public:
/// \name Access Functions
/// @{
/*!
returns:
- `-2` if the triangulation is empty
- `-1` if the triangulation contains a single vertex
- `0` if the triangulation contains exactly two vertices
- `1` if the triangulation contains three (or more) coplanar vertices
- `2` if the triangulation contains at least four non-coplanar vertices
Note that a triangulation of dimension `1` is just a polygon drawn on a circle. The polygon is
not triangulated itself. Thus the triangulation of dimension one consists of one polygon and
has no faces.
*/
int dimension() const;
/*!
returns a const reference to the triangulation traits object.
*/
const Geom_traits& geom_traits() const;
/*!
returns a const reference to the triangulation data structure.
*/
const Triangulation_data_structure& tds() const;
/*!
returns the number of vertices.
*/
size_type number_of_vertices() const;
/*!
returns the number of faces. Note that this includes ghost faces.
*/
size_type number_of_faces() const;
/*!
returns the number of ghost faces.
*/
size_type number_of_ghost_faces() const;
/*!
returns the number of solid faces.
*/
size_type number_of_solid_faces() const;
/// @}
/// \name Geometric Access Functions
/// @{
/*!
returns the geometric position of the vertex `v`.
*/
const Point& point(const Vertex_handle v);
/*!
returns the geometric position of the `i`-th vertex of the face `f`.
*/
const Point& point(const Face_handle f, const int i);
/*!
returns the 3D line segment formed by the vertices of the edge `e`.
\pre `t.dimension()` \f$ \geq1\f$.
*/
Segment_3 segment(const Edge& e) const;
/*!
returns the 3D line segment formed by the vertices of the edge `(f, i)`.
\pre `t.dimension()` \f$ \geq1\f$.
*/
Segment_3 segment(const Face_handle f, int i) const;
/*!
returns the great circle arc formed by the vertices of the edge `e`.
\pre `t.dimension()` \f$ \geq1\f$.
*/
Arc_on_sphere_2 segment_on_sphere(const Edge& e) const;
/*!
returns the great circle arc formed by the vertices of the edge `(f, i)`.
\pre `t.dimension()` \f$ \geq1\f$.
*/
Arc_on_sphere_2 segment_on_sphere(const Face_handle f, int i) const;
/*!
returns the 3D triangle formed by the three vertices of the face `f`.
\pre `t.dimension()` \f$ \geq2\f$.
*/
Triangle_3 triangle(const Face_handle f) const;
/// @}
/// \name All Face, Edge and Vertex Iterators
///
/// The following iterators allow respectively to visit all faces, edges and
/// vertices of the triangulation. These iterators are non mutable,
/// bidirectional and their value types are respectively `Face`,
/// `Edge` and `Vertex`. They are all invalidated by any change in the triangulation.
///
/// @{
/*!
Starts at an arbitrary vertex.
*/
Vertices_iterator vertices_begin() const;
/*!
Past-the-end iterator
*/
Vertices_iterator vertices_end() const;
/*!
returns a range of iterators over all vertices.
\note While the value type of `Vertices_iterator` is `Vertex`, the value type of
`Vertex_handles::iterator` is `Vertex_handle`.
*/
Vertex_handles vertex_handles() const;
/*!
Starts at an arbitrary edge.
*/
All_edges_iterator all_edges_begin() const;
/*!
Past-the-end iterator
*/
All_edges_iterator all_edges_end() const;
/*!
returns a range of iterators over all edges.
*/
All_edges all_edges() const;
/*!
Starts at an arbitrary face.
*/
All_faces_iterator all_faces_begin() const;
/*!
Past-the-end iterator
*/
All_faces_iterator all_faces_end() const;
/*!
returns a range of iterators over all faces.
\note While the value type of `All_faces_iterator` is `Face`, the value type of
`All_face_handles::iterator` is `Face_handle`.
*/
All_face_handles all_face_handles() const;
/*!
Starts at an arbitrary face.
*/
Solid_faces_iterator solid_faces_begin() const;
/*!
Past-the-end iterator
*/
Solid_faces_iterator solid_faces_end() const;
/*!
returns a range of iterators over all solid faces.
\note While the value type of `Solid_faces_iterator` is `Face`, the value type of
`Solid_face_handles::iterator` is `Face_handle`.
*/
Solid_face_handles solid_face_handles() const;
/*!
Starts at an arbitrary face.
*/
Solid_edges_iterator solid_edges_begin() const;
/*!
Past-the-end iterator
*/
Solid_edges_iterator solid_edges_end() const;
/*!
returns a range of iterators over all solid edges.
*/
Solid_edges solid_edges() const;
/*!
returns a range of iterators over all the points of the triangulations.
*/
Points points() const;
/// @}
/// \name Face, Edge and Vertex Circulators
///
/// The triangulation also provides circulators that allows to visit respectively
/// all faces or edges incident to a given vertex or all vertices adjacent to a given vertex.
/// These circulators are non-mutable and bidirectional. The `operator++` moves the
/// circulator counterclockwise around the vertex while the `operator-` moves clockwise.
/// A face circulator is invalidated by any modification of the face pointed to.
/// An edge or a vertex circulator are invalidated by any modification of one of the two
/// faces incident to the edge pointed to.
///
/// @{
/*!
Starts at an arbitrary vertex adjacent to the vertex `v`.
*/
Vertex_circulator adjacent_vertices(Vertex_handle v) const;
/*!
Starts at the first vertex of `f` adjacent to `v` in counterclockwise order around `v`.
\pre The face `f` is incident to the vertex `v`.
*/
Vertex_circulator adjacent_vertices(Vertex_handle v, Face_handle f);
/*!
Starts at an arbitrary edge incident to the vertex `v`.
*/
Edge_circulator incident_edges(Vertex_handle v) const;
/*!
Starts at the first edge of `f` incident to `v`, in counterclockwise order around `v`.
\pre The face `f` is incident to the vertex `v`.
*/
Edge_circulator incident_edges(Vertex_handle v, Face_handle f) const;
/*!
Starts at an arbitrary face incident to the vertex `v`.
Note that this may contain ghost faces.
*/
Face_circulator incident_faces(Vertex_handle v) const;
/*!
Starts at face `f`.
\pre The face `f` is incident to the vertex `v`.
*/
Face_circulator incident_faces(Vertex_handle v, Face_handle f) const;
/// @}
/// \name Combinatorial Predicates
///
/// The class `Triangulation_on_sphere_2` provides methods to test the
/// presence of a particular feature (edge or face) in the triangulation.
///
/// @{
/*!
returns `true` if there exists an edge (ghost or solid) having `va` and `vb` as vertices.
*/
bool is_edge(Vertex_handle va, Vertex_handle vb);
/*!
returns `true` if there exists an edge (ghost or solid) having `va` and `vb` as vertices.
If `true` is returned, the edge with vertices `va` and `vb` is the edge `e=(fr,i)`
where `fr` is a handle to the face incident to `e` and on the right side of `e` oriented from `va` to `vb`.
*/
bool is_edge(Vertex_handle va, Vertex_handle vb, Face_handle& fr, int& i);
/*!
returns `true` if there exists a face (ghost or solid) having `v1`, `v2` and `v3` as vertices.
*/
bool is_face(Vertex_handle v1, Vertex_handle v2, Vertex_handle v3);
/*!
returns `true` if there exists a face (ghost or solid) having `v1`, `v2` and `v3` as vertices.
If `true` is returned, `fr` is a handle to the face with `v1`, `v2` and `v3` as vertices.
*/
bool is_face(Vertex_handle v1, Vertex_handle v2, Vertex_handle v3, Face_handle& fr);
/// @}
/// \name Queries
///
/// The class `Triangulation_on_sphere_2` provides methods to locate a given
/// point with respect to a triangulation. It also provides methods to locate
/// a point with respect to a given face of the triangulation.
///
/// @{
/*!
specifies which case occurs when locating a query point in the triangulation.
*/
enum Locate_type { VERTEX=0, /*!< when the point coincides with a vertex of the triangulation */
EDGE, /*!< when the point is in the relative interior of an edge */
FACE, /*!< when the point is in the interior of a face */
OUTSIDE_CONVEX_HULL, /*!< when the point is on the same 3D plane as the existing vertices, but not on an existing edge */
OUTSIDE_AFFINE_HULL, /*!< when the insertion of the point would increase the dimension of the triangulation. */
NOT_ON_SPHERE, /*!< when the point is not on the sphere */
TOO_CLOSE /*!< when the point is too close to a vertex of the triangulation (see `TriangulationOnSphereTraits_2` for more details) */
};
/*!
locates the point `query` in the triangulation, and returns information on this location.
If the point is (according to the traits) not on the sphere or is too close to an existing vertex,
or if the dimension of the triangulation is not 2, or if the point is outside the affine hull,
then the returned `Face_handle` is `nullptr`.
Otherwise, a face intersected by the ray spawned from the center of the sphere
and passing through the point is returned.
The optional `Face_handle` argument, if provided, is used as a hint
of where the locate process should start its search.
*/
Face_handle locate(const Point& query, Face_handle f = Face_handle()) const;
/*!
Same as above. Additionally, the parameters `lt` and `li` describe where the query point is located.
The variable `lt` is set to the locate type of the query.
If `lt==VERTEX`, the variable `li` is set to the index of the vertex,
and if `lt==EDGE` `li` is set to the index of the vertex opposite to the edge.
Note that `li` has no meaning when the query type is `FACE`, `OUTSIDE_CONVEX_HULL`,
or `OUTSIDE_AFFINE_HULL`.
*/
Face_handle locate(const Point& query, Locate_type& lt, int& li, Face_handle h = Face_handle()) const;
/// @}
/// \name Miscellaneous
///
/// @{
/*!
tests the validity of the triangulation as a `Triangulation_on_sphere_2`.
This function tests the validity of the underlying data structure (using the function
`TriangulationDataStructure_2::is_valid()`), and the validity of the triangulation itself
(consistency between the dimension and the number of simplices, face orientation checks, etc.).
This method is mainly useful for debugging Delaunay triangulation algorithms.
*/
bool is_valid(bool verbose = false, int level = 0) const;
/// @}
};
} // namespace CGAL

36
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/CGAL/Triangulation_on_sphere_face_base_2.h Normal file
View File

@ -0,0 +1,36 @@
namespace CGAL {
/*!
\ingroup PkgTriangulationOnSphere2VertexFaceClasses
\cgalModels `TriangulationOnSphereFaceBase_2`
The class `Triangulation_on_sphere_face_base_2` is a model for the concept
`TriangulationOnSphereFaceBase_2`. It is the default face base class
for triangulations on the sphere.
`Triangulation_on_sphere_face_base_2` can be simply plugged into the triangulation data
structure of a triangulation, or used as a base class to derive
other base face classes tuned for specific applications.
\tparam Traits must be a geometric traits class.
The geometric traits is actually not used by the class.
\tparam Fb has to be a model of the concept `TriangulationDSFaceBase_2`
and will serve as a base class for `Triangulation_on_sphere_face_base_2`.
By default this parameter is instantiated by `Triangulation_ds_face_base_2<>`.
\sa `CGAL::Triangulation_ds_face_base_2<Tds>`
\sa `CGAL::Triangulation_face_base_with_info_2<Traits,Tds>`
\sa `CGAL::Triangulation_on_sphere_vertex_base_2<Traits,Vb>`
\sa `CGAL::Triangulation_on_sphere_2<Traits,Tds>`
*/
template< typename Traits, typename Fb >
class Triangulation_on_sphere_face_base_2
: public Fb
{
public:
};
} // namespace CGAL

47
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/CGAL/Triangulation_on_sphere_vertex_base_2.h Normal file
View File

@ -0,0 +1,47 @@
namespace CGAL {
/*!
\ingroup PkgTriangulationOnSphere2VertexFaceClasses
\cgalModels `TriangulationOnSphereVertexBase_2`
The class `Triangulation_on_sphere_vertex_base_2` is the default model
for the concept `TriangulationOnSphereVertexBase_2` and
the default vertex base class for triangulations on the sphere.
`Triangulation_on_sphere_vertex_base_2` can be simply plugged into the triangulation data
structure of a triangulation, or used as a base class to derive
other base vertex classes tuned for specific applications.
\tparam Traits must be a geometric traits class which provides the type `Point_on_sphere_2`.
It is strongly recommended to instantiate this traits class with the model used
for the triangulation traits class. This ensures that the point type defined by
`Triangulation_on_sphere_vertex_base_2` is the same as the point type defined by
the triangulation.
\tparam Vb must be a model of the concept `TriangulationDSVertexBase_2`.
By default this parameter is instantiated by `Triangulation_ds_vertex_base_2<>`.
\sa `CGAL::Triangulation_ds_vertex_base_2<Tds>`
\sa `CGAL::Triangulation_on_sphere_face_base_2<Traits,Fb>`
\sa `CGAL::Triangulation_vertex_base_with_info_2<Info,Traits,Vb>`
\sa `CGAL::Triangulation_on_sphere_2<Traits,Tds>`
*/
template< typename Traits, typename Vb >
class Triangulation_on_sphere_vertex_base_2
: public Vb
{
public:
/// \name Types
/// @{
/*!
*/
typedef Traits::Point_on_sphere_2 Point;
/// @}
};
} // namespace CGAL

65
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/Concepts/DelaunayTriangulationOnSphereTraits_2.h Normal file
View File

@ -0,0 +1,65 @@
/*!
\ingroup PkgTriangulationOnSphere2Concepts
\cgalConcept
\cgalRefines TriangulationOnSphereTraits_2
The concept `DelaunayTriangulationOnSphereTraits_2` describes the set of requirements
to be fulfilled by any class used to instantiate the first template
parameter of the class `CGAL::Delaunay_triangulation_on_sphere_2<Traits, Tds>`.
To the requirements listed within the concept `TriangulationOnSphereTraits_2`,
this concept adds types and functors requirements related to build the dual on the sphere.
\cgalHasModel `CGAL::Delaunay_triangulation_on_sphere_traits_2`
\cgalHasModel `CGAL::Projection_on_sphere_traits_3`
*/
class DelaunayTriangulationOnSphereTraits_2
{
public:
/// \name Constructions
///
/// @{
/// Construction object type. Must provide the operator:
///
/// `Point_on_sphere_2 operator()(Point_on_sphere_2 p, Point_on_sphere_2 q, Point_on_sphere_2 r)`
///
/// which returns the intersection of the dual of the face defined by the three points `p`, `q`, and `r`,
/// and the sphere, on the positive side of the plane defined by `p`, `q`, and `r`. The dual of the face
/// is the line orthogonal to the face, passing through the center of the smallest circumscribing sphere
/// of the face.
///
/// \pre the center of the sphere is on the negative side of the plane defined by `p`, `q`, and `r`.
///
/// \note This type is only required for the computation of dual objects (Voronoi vertices and edges)
/// and a dummy type can be used otherwise.
typedef unspecified_type Construct_circumcenter_on_sphere_2;
/// Construction object type. Must provide the operator:
///
/// `Point_3 operator()(Point_3 p, Point_3 q, Point_3 r)`
///
/// which returns the center of the circle circumscribed to the face with vertices `p`, `q`, and `r`.
///
/// \note This type is only required for the computation of dual objects (Voronoi vertices and edges)
/// and a dummy type can be used otherwise.
typedef unspecified_type Construct_circumcenter_3;
/// @}
public:
/// \name Operations
///
/// The following functions give access to the predicate and constructor objects.
///
/// @{
///
Construct_circumcenter_on_sphere_2 construct_circumcenter_on_sphere_2_object();
///
Construct_circumcenter_3 construct_circumcenter_3_object();
/// @}
};

35
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/Concepts/TriangulationOnSphereFaceBase_2.h Normal file
View File

@ -0,0 +1,35 @@
/*!
\ingroup PkgTriangulationOnSphere2Concepts
\cgalConcept
\cgalRefines `TriangulationDSFaceBase_2`
The concept `TriangulationOnSphereFaceBase_2` describes the requirements for a face class
of a triangulation data structure that is itself plugged into a triangulation on the sphere.
The data structure concept `TriangulationDataStructure_2` was primarily designed
to serve as a data structure for the 2D triangulation classes of \cgal, which are triangulations
embedded in the 2D Euclidean plane.
However, its genericity makes it usable for any orientable triangulated surface without boundary
regardless of the dimensionality of the space the triangulation is embedded in, and thus
it is a valid data structure for the triangulations on the sphere of this package.
A departing feature is that if the vertices of the triangulation are all lying on one hemisphere,
then the triangulation is not an orientable triangulated surface without boundary.
In this case, fictitious faces are added to the triangulation, called <i>ghost faces</i>,
such that the triangulation is a topological sphere.
\cgalHasModel `CGAL::Triangulation_on_sphere_face_base_2<Traits, Fb>`
\sa `TriangulationDataStructure_2`
\sa `TriangulationOnSphereVertexBase_2`
*/
class TriangulationOnSphereFaceBase_2
{
public:
/// provides read access to a Boolean used to indicate if the face is a ghost face.
bool is_ghost();
/// provides write access to a Boolean used to indicate if the face is a ghost face.
void set_ghost(const bool b);
};

267
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/Concepts/TriangulationOnSphereTraits_2.h Normal file
View File

@ -0,0 +1,267 @@
/*!
\ingroup PkgTriangulationOnSphere2Concepts
\cgalConcept
\cgalRefines SpatialSortingTraits_3
The concept `TriangulationOnSphereTraits_2` describes the set of requirements
to be fulfilled by any class used to instantiate the first template
parameter of the class `CGAL::Triangulation_on_sphere_2<Traits, Tds>`.
This concept provides the types of the geometric primitives used in the
triangulation and the function object types for the required predicates on those primitives.
In particular, the traits class is expected to contain information about the sphere (center and radius)
on which the points of the triangulation lie, as well as to provide means to change this information.
\cgalHasModel `CGAL::Delaunay_triangulation_on_sphere_traits_2`
\cgalHasModel `CGAL::Projection_on_sphere_traits_3`
\sa `DelaunayTriangulationOnSphereTraits_2`
*/
class TriangulationOnSphereTraits_2
{
public:
/// \name Types
/// @{
/*!
The number type; must be a model of `FieldNumberType`.
*/
typedef unspecified_type FT;
/*!
The point type of the triangulation, representing a point on the sphere.
*/
typedef unspecified_type Point_on_sphere_2;
/*!
An arc of a great circle, used to represent a curved segment (a Voronoi or a Delaunay edge).
*/
typedef unspecified_type Arc_on_sphere_2;
/*!
The 3D point type.
*/
typedef unspecified_type Point_3;
/*!
The 3D segment type.
*/
typedef unspecified_type Segment_3;
/*!
The 3D triangle type.
*/
typedef unspecified_type Triangle_3;
/// @}
public:
/// \name Predicates
/// @{
/// Predicate object type. Must provide the operator:
///
/// `bool operator()(Point_on_sphere_2 p, Point_on_sphere_2 q, Point_on_sphere_2 r)`
///
/// which returns `true` if `r` strictly lies on the shortest path between `p` and `q`
/// on the great circle passing through `p` and `q`.
///
/// If `p` and `q` are diametrically opposed, `true` is returned.
typedef unspecified_type Collinear_are_strictly_ordered_on_great_circle_2;
/// Predicate object type. Must provide the operator:
///
/// `bool operator()(Point_on_sphere_2 p, Point_on_sphere_2 q)`
///
/// which returns `true` or `false` based on a consistent, user-defined strict total order of points.
typedef unspecified_type Compare_on_sphere_2;
/// Predicate object type. Must provide the operator:
///
/// `bool operator()(Point_on_sphere_2 p, Point_on_sphere_2 q)`
///
/// which returns `true` if, and only if, `p` and `q` are equal.
typedef unspecified_type Equal_on_sphere_2;
/// Predicate object type. Must provide the operator:
///
/// `Orientation operator()(Point_on_sphere_2 p, Point_on_sphere_2 q, Point_on_sphere_2 r)`
///
/// which returns `CGAL::POSITIVE`, if `r` lies on the left hemisphere while walking
/// the shortest path between `p` and `q` on the great circle defined by `p` and `q`;
/// returns `CGAL::NEGATIVE` if `r` lies on the right hemisphere; and returns `CGAL::COPLANAR`
/// if `r` lies on the great circle.
///
/// If `p` and `q` are diametrically opposed, `CGAL::COPLANAR` is returned.
typedef unspecified_type Orientation_on_sphere_2;
/// Predicate object type. Must provide the operator:
///
/// `bool operator()(Point_on_sphere_2 p, Point_on_sphere_2 q, Point_on_sphere_2 r, Point_on_sphere_2 s)`
///
/// which returns `CGAL::POSITIVE`, if `s` lies on the positive side of the oriented plane `h`
/// defined by `p`, `q`, and `r`; returns `CGAL::NEGATIVE` if `s` lies on the negative side of `h`;
/// and returns `CGAL::COPLANAR` if `s` lies on `h`.
typedef unspecified_type Side_of_oriented_circle_on_sphere_2;
/// @}
/// \name Constructions
///
/// @{
/// Construction object. Must provide the operator:
///
/// `Arc_on_sphere_2 operator()(Point_on_sphere_2 p, Point_on_sphere_2 q)`
///
/// which introduces an arc of great circle, with source `p` and target `q`.
///
/// The circular arc constructed from a source, and a target, is defined as the set of points
/// of the great circle that lie between the source `p` and the target `q`,
/// when traversing the great circle counterclockwise seen from the side of the plane
/// of the great circle pointed by its <I>positive</I> normal vectors.
///
/// In this definition, we say that a normal vector \f$ (a,b,c)\f$ is <I>positive</I> if
/// \f$ (a,b,c)>(0,0,0)\f$ (i.e.\ \f$ (a>0) || (a==0) \&\& (b>0) || (a==0)\&\&(b==0)\&\&(c>0)\f$).
///
/// \pre `p` and `q` are not diametrically opposed.
typedef unspecified_type Construct_arc_on_sphere_2;
/// Construction object. Must provide the operator:
///
/// `Point_on_sphere_2 operator()(Point_3 p)`,
///
/// which expresses the point `p`, a 3D point living on the sphere,
/// using the (user-defined) point-on-sphere representation.
typedef unspecified_type Construct_point_on_sphere_2;
/// Construction object. Must provide the operator:
///
/// `Point_3 operator()(Point_on_sphere_2 p)`,
///
/// which expresses the point `p` on the sphere as a point in the embedding 3D space.
typedef unspecified_type Construct_point_3;
/// Construction object. Must provide the operator:
///
/// `Segment_3 operator()(Point_3 p, Point_3 q)`,
///
/// which constructs a segment from two 3D points.
typedef unspecified_type Construct_segment_3;
/// Construction object. Must provide the operator:
///
/// `Triangle_3 operator()(Point_3 p, Point_3 q, Point_3 r)`,
///
/// which constructs a triangle from three 3D points.
typedef unspecified_type Construct_triangle_3;
/// @}
public:
/// \name Operations
///
/// The following functions provide access to the predicate and constructor objects.
///
/// @{
// Predicates
///
Compare_on_sphere_2 compare_on_sphere_2_object();
///
Equal_on_sphere_2 equal_on_sphere_2_object();
///
Orientation_on_sphere_2 orientation_on_sphere_2_object();
///
Side_of_oriented_circle_on_sphere_2 side_of_oriented_circle_on_sphere_2_object();
///
Collinear_are_strictly_ordered_on_great_circle_2 collinear_are_strictly_ordered_on_great_circle_2_object();
// Constructions
///
Construct_arc_on_sphere_2 construct_arc_on_sphere_2_object();
///
Construct_point_on_sphere_2 construct_point_on_sphere_2_object();
///
Construct_point_3 construct_point_3_object();
///
Construct_segment_3 construct_segment_3_object();
///
Construct_triangle_3 construct_triangle_3_object();
/// @}
public:
/// \name Configuration of the spherical domain
///
/// @{
/// sets the center of the sphere.
///
/// \note This function is meant to be used only by the triangulation class as modifying the domain
/// requires clearing the triangulation. Users can change the domain using
/// `CGAL::Triangulation_on_sphere_2::set_center_and_radius()`.
void set_center(Point_3 center);
/// returns the center of the sphere.
Point_3 center() const;
/// sets the radius of the sphere.
///
/// \note This function is meant to be used only by the triangulation class as modifying the domain
/// requires clearing the triangulation. Users can change the domain using
/// `CGAL::Triangulation_on_sphere_2::set_center_and_radius()`.
void set_radius(FT radius);
/// returns the radius of the sphere.
FT radius() const;
/// @}
public:
/// \name Precision predicates
///
/// When a kernel does not offer arbitrary precision - or at least an exact representation of
/// algebraic coordinates, which is usually the case for the sake of computational speed,
/// then most points that are theoretically on the sphere are not in practice exactly on the sphere.
/// When this is the case, it cannot be guaranteed that all points are in a convex position,
/// which can then result in points being hidden upon insertion, see \cgalCite{cgal:ccplr-redtp-10}.
/// To ensure that no point is hidden, a trick is to enforce a gap between points (\cgalCite{cgal:ccplr-redtp-10},
/// Lemma 4.1). This gap is based on a fixed maximal allowed distance \f$ \delta \f$ between a point
/// and the theoretical sphere. The gap must then be at least \f$ 2 \sqrt{R\delta} \f$, where \f$ R \f$
/// is the radius of the sphere.
///
/// The following two predicates serve to check if the point is on the sphere, that is
/// if for a model-defined \f$ \delta \f$ the distance to the sphere is smaller
/// than \f$ \delta \f$ and if the distance between two points is greater than \f$ 2 \sqrt{R\delta} \f$.
/// It is also of course possible to construct traits classes with arbitrary precision or
/// with a representation that ensures that points are exactly on the sphere and in this case,
/// \f$ \delta = 0 \f$ can be chosen and the predicates below are trivial.
///
/// \sa `CGAL::Delaunay_triangulation_on_sphere_traits_2`
///
/// @{
/// returns whether the point `p` is on the sphere or not.
///
/// \note A point that is not on the sphere will not be inserted in the triangulation.
bool is_on_sphere(Point_on_sphere_2 p);
/// returns whether the points `p` and `q` are too close (see
/// <a href="https://doc.cgal.org/latest/Triangulation_on_sphere_2/classTriangulationOnSphereTraits__2.html#amgrpbc40b59d2ea2be0fc60516e5f7bb0f92">
/// Precision predicates</a>).
bool are_points_too_close(Point_on_sphere_2 p, Point_on_sphere_2 q);
/// @}
};

83
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/Concepts/TriangulationOnSphereVertexBase_2.h Normal file
View File

@ -0,0 +1,83 @@
/*!
\ingroup PkgTriangulationOnSphere2Concepts
\cgalConcept
\cgalRefines `TriangulationDSVertexBase_2`
The concept `TriangulationOnSphereVertexBase_2` describes the requirements for the vertex base class
of a triangulation data structure to be plugged in a triangulation on the sphere.
It refines the concept `TriangulationDSVertexBase_2`, adding geometric information:
the vertex base of a triangulation stores a point.
\cgalHasModel `CGAL::Triangulation_on_sphere_vertex_base_2<Traits,Vb>`
\sa `TriangulationDataStructure_2`
\sa `TriangulationOnSphereFaceBase_2`
*/
class TriangulationOnSphereVertexBase_2
{
public:
/// \name Types
/// @{
/*!
Must be the equivalent to the point type `Point_on_sphere_2` of the concept `TriangulationOnSphereTraits_2`.
*/
typedef unspecified_type Point;
/*!
Must be the equivalent to the point type `Face_handle` of the triangulation.
*/
typedef unspecified_type Face_handle;
/// @}
/// \name Creation
/// @{
/*!
constructs a vertex whose geometric embedding is the point `p`.
*/
TriangulationOnSphereVertexBase_2(Point p);
/*!
constructs a vertex whose geometric embedding is the point `p` and pointing to face `f`.
`Face_handle` is the type defined the triangulation data structure used in the Delaunay triangulation
on the sphere.
*/
TriangulationOnSphereVertexBase_2(Point p, Face_handle f);
/// @}
/// \name Access Functions
/// @{
/*!
returns the point.
*/
Point point() const;
/// @}
/// \name Setting
/// @{
/*!
sets the point.
*/
void set_point(Point p);
/*!
inputs the non-combinatorial information given by the vertex:
the point and other possible information.
*/
std::istream& operator>>(std::istream& is, TriangulationOnSphereVertexBase_2& v);
/*!
outputs the non combinatorial operation given by the vertex: the
point and other possible information.
*/
std::ostream& operator<<(std::ostream& os, const TriangulationOnSphereVertexBase_2& v);
/// @}
};

7
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/Doxyfile.in Normal file
View File

@ -0,0 +1,7 @@
@INCLUDE = ${CGAL_DOC_PACKAGE_DEFAULTS}
PROJECT_NAME = "CGAL ${CGAL_DOC_VERSION} - 2D Triangulations on the Sphere"
HTML_EXTRA_FILES = ${CGAL_PACKAGE_DOC_DIR}/fig/geometric_embedding.png \
${CGAL_PACKAGE_DOC_DIR}/fig/ghost_faces.png \
${CGAL_PACKAGE_DOC_DIR}/fig/header.png

81
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/PackageDescription.txt Normal file
View File

@ -0,0 +1,81 @@
/// \defgroup PkgTriangulationOnSphere2Ref 2D Triangulation on Sphere Reference
/// \defgroup PkgTriangulationOnSphere2Concepts Concepts
/// \ingroup PkgTriangulationOnSphere2Ref
/// \defgroup PkgTriangulationOnSphere2TriangulationClasses Triangulation Classes
/// \ingroup PkgTriangulationOnSphere2Ref
/// \defgroup PkgTriangulationOnSphere2VertexFaceClasses Vertex and Face Classes
/// \ingroup PkgTriangulationOnSphere2Ref
/*!
\addtogroup PkgTriangulationOnSphere2Ref
\cgalPkgDescriptionBegin{2D Triangulations on Sphere,PkgTriangulationOnSphere2}
\cgalPkgPicture{triangulation_on_sphere.png}
\cgalPkgSummaryBegin
\cgalPkgAuthors{Mael Rouxel-Labbé, Monique Teillaud, and Claudia Werner}
\cgalPkgDesc{This package enables the construction and manipulation of Delaunay triangulations on the 2-sphere.
Triangulations are built incrementally and can be modified by insertion or removal of vertices.
Point location querying and primitives to build the dual Voronoi diagram are provided.}
\cgalPkgManuals{Chapter_2D_Triangulations_on_sphere,PkgTriangulationOnSphere2Ref}
\cgalPkgSummaryEnd
\cgalPkgShortInfoBegin
\cgalPkgSince{5.3}
\cgalPkgDependsOn{\ref PkgTDS2}
\cgalPkgBib{cgal:rtw-tos2}
\cgalPkgLicense{\ref licensesGPL "GPL"}
\cgalPkgDemo{Delaunay Triangulation, delaunay_Triangulation_on_sphere_2.zip}
\cgalPkgShortInfoEnd
\cgalPkgDescriptionEnd
A triangulation is a 2-dimensional simplicial complex that is pure,
connected, and without singularities. It can be viewed as a collection of triangular faces,
such that two faces either have an empty intersection or share an edge or a vertex.
This package handles triangulations of points on the 2-sphere, and is thus tightly linked
to the package \ref PkgTriangulation2, which handles triangulations of points in Euclidean 2D space,
since both domains are 2-manifolds. Concepts and classes are therefore similar
to those of the package \ref PkgTriangulation2.
The basic elements of the representation are vertices and faces.
Each triangular face gives access to its three incident vertices
and to its three adjacent faces.
Each vertex gives access to one of its incident faces
and, through that face, to the circular list of its incident faces.
The edges are not explicitly represented, they are only represented
through the adjacency relations of two faces.
The triangulation classes of \cgal depend on two template parameters.
The first template parameter stands for a geometric traits class which is assumed to provide
the geometric objects (points, segments, ...) forming the triangulation
and the geometric predicates on those objects.
The second template parameter stands for a model of triangulation data structure
acting as a container for faces and vertices while taking care
of the combinatorial aspects of the triangulation.
The concepts and models relative to the triangulation data structure
are described in Chapter \ref PkgTDS2Ref "2D Triangulation Data Structure".
\cgalClassifedRefPages
\cgalCRPSection{Concepts}
- `TriangulationOnSphereTraits_2`
- `DelaunayTriangulationOnSphereTraits_2`
- `TriangulationOnSphereFaceBase_2`
- `TriangulationOnSphereVertexBase_2`
\cgalCRPSection{Traits Classes}
- `CGAL::Delaunay_triangulation_on_sphere_traits_2<LK, SK>`
- `CGAL::Projection_on_sphere_traits_3<LK, SK>`
\cgalCRPSection{Triangulation Classes}
- `CGAL::Triangulation_on_sphere_face_base_2<Traits, Fb>`
- `CGAL::Triangulation_on_sphere_vertex_base_2<Traits, Vb>`
- `CGAL::Triangulation_on_sphere_2<Traits, TDS>`
- `CGAL::Delaunay_triangulation_on_sphere_2<Traits, TDS>`
\cgalCRPSection{Enum}
- \link CGAL::Triangulation_on_sphere_2::Locate_type `CGAL::Triangulation_on_sphere_2<Traits,TDS>::Locate_type` \endlink
*/

241
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/Triangulation_on_sphere_2.txt Normal file
View File

@ -0,0 +1,241 @@
namespace CGAL {
/*!
\mainpage User Manual
\anchor Chapter_2D_Triangulations_on_sphere
\cgalAutoToc
\author Mael Rouxel-Labbé, Monique Teillaud, and Claudia Werner
\cgalFigureAnchor{fig_TOS2_header}
<center>
<img src="header.png" style="max-width:70%;"/>
</center>
This chapter describes two-dimensional triangulations on the sphere of \cgal. It is organized as follows:
<ul>
<li>Section \ref Section_2D_ToS_Definitions introduces the main definitions about triangulations on the sphere.</li>
<li>Section \ref Section_2D_ToS_Implementation details the way two-dimensional triangulations on the sphere
are implemented and represented in \cgal.</li>
<li>Section \ref Section_2D_ToS_Examples lists a few basic examples to illustrate the capabilities
of the package.</li>
</ul>
\section Section_2D_ToS_Definitions Definitions
Denote by \f$ \mathbb{S(c, r)}\f$ the two-dimensional sphere embedded in the Euclidean space \f$ \mathbb{R}^3\f$,
with center `c` and radius `r`, that is
\f$ \mathbb{S(c, r)} = \left\{ x \in \mathbf{R}^3 : \left\| x - c \right\| = r \right\} \f$.
When the parameters `c` and `r` are not important, they will be omitted and \f$ \mathbb{S}\f$
will be used directly.
Given a set \f$ \mathcal{P}\f$ of points on \f$ \mathbb{S(c, r)}\f$, a <em>two-dimension triangulation</em>
of \f$ \mathcal{P}\f$ can be described as a two-dimensional simplicial complex that is pure,
connected, and without singularity whose vertices are exactly the points in \f$ \mathcal{P}\f$
(see the complete definition in the package
\link Section_2D_Triangulations_Definitions 2D Triangulations\endlink).
In \f$ \mathbb{R}^2\f$, a <em>Delaunay</em> triangulation is a two-dimension triangulation
that satisfies the <em>empty circle property</em> (also called <em>Delaunay property</em>):
the circumscribing circle of any facet of the triangulation contains no point in its interior.
This definition naturally extends to the two-dimensional sphere, as illustrated in the figure below.
\cgalFigureBegin{fig_del_def, delaunay_on_sphere_def.svg}
Delaunay property: the circumscribing circle (in green) on the sphere
of the Delaunay face \f$ p_1 p_2 p_3\f$ is empty.
\cgalFigureEnd
A particularity of Delaunay triangulations on \f$ \mathbb{S}\f$ is that the Delaunay property
test can be reduced to a simple orientation test in the three-dimensional space: indeed, the circumscribing
circle of a Delaunay face is the intersection of the plane passing through the three vertices
of the face with \f$ \mathbb{S}\f$, and determining whether a fourth point `q` is within the
circumscribing circle or not is equivalent to checking whether `q` is above or below the said supporting plane.
\subsection Section_2D_ToS_Imprecision Representation of a Point on the Sphere
The theoretical definition of Delaunay triangulations in the previous section assumes that points
in \f$ \mathcal{P}\f$ lie exactly on the sphere. In a real world however, one usually does not
manipulate points that lie exactly on \f$ \mathbb{S}\f$ because the chosen number type
is not able to represent square roots exactly, such as \f$ \mathtt{float}\f$, or \f$ \mathtt{double}\f$.
Some specific number types are able to represent exactly all points on the sphere, for example
the class `CORE::Expr`, but this comes with the drawback of a substantially higher computational cost.
This lack of exact representation is an obvious problem as the interpretation of the Delaunay property
as an orientation test does not stand if points do not lie in a convex position.
This gap between the theoretical and the practical settings was addressed by Caroli et al. \cgalCite{cgal:ccplr-redtp-10} :
the solution is to use a <em>regular</em> triangulation, which is a generalization of the Delaunay triangulation
to sets of <em>weighted points</em>
(see \link Section_2D_Triangulations_Regular 2D Regular Triangulations\endlink for more information).
A weighted point \f$(p,w)\f$ of \f$ \mathbb{R}^2\f$ can naturally be seen as as a circle with center \f$ p\f$
and radius \f$ r\f$ such that \f$ r^2 = w\f$ and similarly to Delaunay triangulations and the definition
of regular triangulations in \f$ \mathbb{R}^2\f$ can be naturally extended to circles on \f$ \mathbb{S}\f$.
Given a set of points \f$ \mathcal{P}_3\f$ living in \f$ \mathbb{R}^3\f$, Caroli et al. compute
a regular triangulation of the set of weighted points \f$ \mathcal{P}_w\f$ where the weighted
points of \f$ \mathcal{P}_w\f$ are the projection of the points of \f$ \mathcal{P}_3\f$ onto \f$ \mathbb{S}\f$,
and whose weights are based on the projection distance. Furthermore, they show that if the Euclidean distance
between the points of \f$ \mathcal{P}_w\f$ is sufficiently large, and if the points are sufficiently
close to the sphere, then the Delaunay triangulation of the Delaunay triangulation of \f$ \mathcal{P}_3\f$
is exactly the regular triangulation of \f$ \mathcal{P}_w\f$.
As a consequence is that one can manipulate 3D points that are not exactly on the sphere,
and it possible to build the Delaunay triangulation without taking the weights into account, as long as the
points are sufficiently far from each other.
Caroli et al. et provide bounds on this separation criterion for the case of the \f$ \mathtt{double}\f$
number type: two points must be at least \f$ 2^{-25}r\f$ apart. This is a condition that
is generally satisfied for most inputs: if \f$ r\f$ were for example the radius of the Earth,
roughly 6300 kms, then the point separation requirement would be of the order of 1 meter.
\section Section_2D_ToS_Implementation Implementation
The main class of this package is the class `CGAL::Delaunay_triangulation_on_sphere_2`;
it represents a Delaunay triangulation on the sphere, and provides insertion and removal of vertices,
as well as tools to draw the triangulation and its dual, the Voronoi diagram.
The base of `CGAL::Delaunay_triangulation_on_sphere_2` is the class `CGAL::Triangulation_on_sphere_2`.
It represents a triangulation on the sphere, but does <em>not</em> support insertion or removal
of vertices. Both classes are built on top of a data structure called the triangulation data structure.
The triangulation data structure can be thought of as a container for the faces and vertices
of the triangulation. This data structure also takes care of all the combinatorial aspects of the triangulation.
These triangulation classes are intentionally very similar to `CGAL::Delaunay_triangulation_2` and
`CGAL::Triangulation_2` as both classes represent triangulations of a 2-manifold domain without boundary.
As such, many details about the implementation are not repeated here, and complementary information
can be found in the Sections \ref Section_2D_Triangulations_Representation,
\ref Section_2D_Triangulations_Software_Design, and \ref Section_2D_Triangulations_Basic
of the package \ref PkgTriangulation2.
However, a significant departure from Euclidean 2D triangulations is the following: since
the triangulation data structure represents a 2-manifold without boundary, it is necessary
for 2D triangulations to introduce so-called <i>infinite faces</i> to complete the "real"
triangulation (see \ref Triangulation_2TheSetofFaces).
On the sphere, this trick is not necessary as the triangulation itself is already a 2-manifold
without boundary.
\subsection Section_2D_ToS_Ghost Ghost Faces
There is an exception to the previous statement: in degenerate configurations where all points
lie on the same hemisphere, the Delaunay triangulation on the sphere theoretically has a border.
Internally however, the triangulation data structure must remain a 2-manifold at all time.
To ensure this property, fictitious faces referred to as <i>ghost faces</i> are added.
These faces are characterized by the fact that the center of the sphere does not (strictly) lie
on the positive side of the supporting plane of the face. Conversely, faces that
are not ghost faces are called <em>solid faces</em>, and edges of such faces are <em>solid edges</em>.
\cgalFigureAnchor{fig_TOS2_ghost}
<center>
<img src="ghost_faces.png" style="max-width:70%;"/>
</center>
\cgalFigureCaptionBegin{fig_TOS2_ghost}
Solid (blue) and ghost (orange) faces of a Delaunay triangulation on the sphere. View from above (leftmost),
and below (middle and rightmost).
\cgalFigureCaptionEnd
\subsection Section_2D_ToS_Traits Traits Classes and Choice of Kernel
Two traits classes are offered with this package as models of the concept `DelaunayTriangulationOnSphereTraits_2`:
<ul>
<li>`CGAL::Delaunay_triangulation_on_sphere_traits_2`: This is the simplest possible traits,
which represents points on sphere directly as 3D points. If its kernel template parameter
cannot represent all points on the sphere exactly, it employs the solution
described in Section \ref Section_2D_ToS_Imprecision to determine whether a point
should be considered on the sphere, or too close to an existing vertex.
</li>
<li>`CGAL::Projection_on_sphere_traits_3`: This traits class utilizes a custom internal point type
for points on the sphere: given a point `p` in 3D space, this traits class manipulates directly
its projection on the sphere (that is, the intersection of the sphere and the segment with endpoints
`p` and the center of the sphere). Consequently, all points to be inserted are on the sphere.
This traits enable manipulating points that are not on the sphere, but whose triangulation
on the sphere is still interesting, such as geographical coordinates with altitude.
</li>
</ul>
Both these classes are templated by a kernel. The choice of this kernel is important to ensure
a correct result: for the construction of triangulations to be safe, the kernel should provide
exact predicates (for example, `CGAL::Exact_predicates_inexact_constructions_kernel`). In addition,
some auxiliary functions such as those used in the construction of the Voronoi diagram require creating
new geometric points; as such, a kernel offering an exact representation of points on the sphere
and exact constructions (for example, `CGAL::Exact_predicates_exact_constructions_kernel_with_sqrt`)
should be used if one wants to avoid otherwise inevitable numerical approximations.
\subsection Section_2D_ToS_Lowdim Dimension of a Triangulation on the Sphere
Mostly by convention with other triangulations, the dimension of a triangulation on the sphere
is defined as follows:
<ul>
<li>`-2`, if the triangulation is empty;</li>
<li>`-1`, if the triangulation contains a single vertex;</li>
<li>`0`, if the triangulation contains exactly two vertices;</li>
<li>`1`, if the triangulation contains at least three coplanar vertices (which do not necessarily
lie on a great circle);</li>
<li>`2`, if the triangulation contains at least four non-coplanar vertices.</li>
</ul>
Note that a triangulation of dimension `1` is just a polygon drawn on a circle. The polygon is
not triangulated itself. Thus the triangulation of dimension `1` consists of a planar polygon and
has no faces.
\subsection Section_2D_ToS_Embedding Geometric Embeddings
The descriptions of Delaunay (and regular) triangulations over \f$ \mathbb{S}\f$ have so far been mostly
combinatorial. The question of the geometrical embedding of the simplices of the triangulation
is also interesting. The two natural embedding of edges and faces of a triangulation of a set
of points on \f$ \mathbb{S}\f$ are to use either <em>straight</em> simplex, that is using
three-dimensional segments and triangles for the edges and faces of the triangulation, or
to use a <em>curved</em> embedding, where the edges are arc segments of great circles over \f$ \mathbb{S}\f$.
In the latter choice, the geometrical embedding of the face is defined implicitely by its three edges.
Both choices are available to users, for example using either `Triangulation_on_sphere_2::segment()`
or `Triangulation_on_sphere_2::segment_on_sphere()`. Similar choices are available in the construction
of the dual, the Voronoi diagram.
\cgalFigureAnchor{fig_TOS2_embedding}
<center>
<img src="geometric_embedding.png" style="max-width:70%;"/>
</center>
\cgalFigureCaptionBegin{fig_TOS2_embedding}
Curved (left) and straight (right) geometric embeddings of the edges of the Delaunay triangulation
of the French (including overseas territories) post offices.
\cgalFigureCaptionEnd
\section Section_2D_ToS_Examples Examples
\subsection Section_2D_ToS_Ex_Basic Basic Example
The following example uses the simplest traits class provided in this package. It demonstrates
how to iteratively insert a few points, and how the dimension and the data structure
of the triangulation evolve with these insertions.
\cgalExample{Triangulation_on_sphere_2/triang_on_sphere.cpp}
\subsection Section_2D_ToS_Ex_Exact Using an Exact Kernel
The following example illustrates the limitation of using a kernel with inexact representation
of points on the sphere and the rejection of a point for being too close to an already existing
vertex. A kernel providing exact representation is also shown to be able to insert
these two extremely close points.
\cgalExample{Triangulation_on_sphere_2/triang_on_sphere_exact.cpp}
\subsection Section_2D_ToS_Ex_Project Using the Projection Traits Class
In this example, the use of the projection traits class is illustrated.
\cgalExample{Triangulation_on_sphere_2/triang_on_sphere_proj.cpp}
\subsection Section_2D_ToS_Ex_Range Insertion of a Range and Spatial Sorting
The following example demonstrates how to insert a range of points at once. This enables
an internal algorithm to sort the set of points to ensure locality when inserting points,
which greatly speeds up the insertion.
\cgalExample{Triangulation_on_sphere_2/triang_on_sphere_range.cpp}
\section Section_2D_ToS_Design Design and Implementation History
Prototype code implementing the publication of Caroli et al. \cgalCite{cgal:ccplr-redtp-10}
was developed over the two internships of Oliver Rouiller and Claudia Werner at Inria,
under the supervision of Monique Teillaud and with the help of Sébastien Loriot.
Based on this prototype, Mael Rouxel-Labbé developed the initial version of this package.
*/
} /* namespace CGAL */

11
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/dependencies Normal file
View File

@ -0,0 +1,11 @@
Manual
Kernel_23
Circular_kernel_3
STL_Extension
Algebraic_foundations
Circulator
Spatial_sorting
Stream_support
TDS_2
Triangulation_2
Number_types

6
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/examples.txt Normal file
View File

@ -0,0 +1,6 @@
/*!
\example Triangulation_on_sphere_2/triang_on_sphere.cpp
\example Triangulation_on_sphere_2/triang_on_sphere_exact.cpp
\example Triangulation_on_sphere_2/triang_on_sphere_proj.cpp
\example Triangulation_on_sphere_2/triang_on_sphere_range.cpp
*/

BIN
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/fig/delaunay_on_sphere_def.png Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 46 KiB

86
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/fig/delaunay_on_sphere_def.svg Normal file
View File

@ -0,0 +1,86 @@
<?xml version="1.0" encoding="UTF-8"?>
<svg xmlns="http://www.w3.org/2000/svg" xmlns:xlink="http://www.w3.org/1999/xlink" width="294pt" height="284pt" viewBox="0 0 294 284" version="1.1">
<defs>
<g>
<symbol overflow="visible" id="glyph0-0">
<path style="stroke:none;" d=""/>
</symbol>
<symbol overflow="visible" id="glyph0-1">
<path style="stroke:none;" d="M 9.109375 -9.96875 C 9.109375 -10.09375 9.015625 -10.09375 8.984375 -10.09375 C 8.921875 -10.09375 8.90625 -10.078125 8.734375 -9.875 C 8.65625 -9.765625 8.0625 -9.015625 8.046875 -9.015625 C 7.578125 -9.9375 6.625 -10.09375 6.03125 -10.09375 C 4.203125 -10.09375 2.546875 -8.4375 2.546875 -6.8125 C 2.546875 -5.734375 3.203125 -5.109375 3.90625 -4.859375 C 4.0625 -4.8125 4.90625 -4.578125 5.328125 -4.46875 C 6.0625 -4.28125 6.25 -4.21875 6.5625 -3.90625 C 6.609375 -3.828125 6.90625 -3.5 6.90625 -2.828125 C 6.90625 -1.484375 5.671875 -0.109375 4.234375 -0.109375 C 3.0625 -0.109375 1.75 -0.609375 1.75 -2.21875 C 1.75 -2.5 1.8125 -2.84375 1.84375 -2.984375 C 1.84375 -3.03125 1.859375 -3.09375 1.859375 -3.125 C 1.859375 -3.1875 1.84375 -3.25 1.71875 -3.25 C 1.59375 -3.25 1.578125 -3.234375 1.515625 -2.984375 L 0.78125 -0.046875 C 0.78125 -0.03125 0.734375 0.15625 0.734375 0.171875 C 0.734375 0.296875 0.84375 0.296875 0.875 0.296875 C 0.9375 0.296875 0.953125 0.28125 1.125 0.078125 L 1.78125 -0.78125 C 2.125 -0.265625 2.875 0.296875 4.203125 0.296875 C 6.046875 0.296875 7.75 -1.484375 7.75 -3.28125 C 7.75 -3.890625 7.609375 -4.421875 7.0625 -4.953125 C 6.75 -5.25 6.5 -5.328125 5.171875 -5.671875 C 4.21875 -5.921875 4.09375 -5.96875 3.828125 -6.203125 C 3.578125 -6.4375 3.40625 -6.78125 3.40625 -7.265625 C 3.40625 -8.484375 4.625 -9.703125 5.984375 -9.703125 C 7.390625 -9.703125 8.046875 -8.84375 8.046875 -7.484375 C 8.046875 -7.109375 7.96875 -6.734375 7.96875 -6.671875 C 7.96875 -6.546875 8.09375 -6.546875 8.140625 -6.546875 C 8.265625 -6.546875 8.28125 -6.578125 8.328125 -6.8125 Z M 9.109375 -9.96875 "/>
</symbol>
<symbol overflow="visible" id="glyph1-0">
<path style="stroke:none;" d=""/>
</symbol>
<symbol overflow="visible" id="glyph1-1">
<path style="stroke:none;" d="M 5.609375 -5.390625 C 5.328125 -5.390625 5.203125 -5.390625 5 -5.21875 C 4.921875 -5.15625 4.765625 -4.9375 4.765625 -4.703125 C 4.765625 -4.421875 4.984375 -4.25 5.25 -4.25 C 5.59375 -4.25 5.984375 -4.53125 5.984375 -5.109375 C 5.984375 -5.796875 5.328125 -6.328125 4.328125 -6.328125 C 2.453125 -6.328125 0.578125 -4.28125 0.578125 -2.234375 C 0.578125 -0.984375 1.34375 0.140625 2.8125 0.140625 C 4.765625 0.140625 6 -1.375 6 -1.5625 C 6 -1.65625 5.90625 -1.71875 5.859375 -1.71875 C 5.8125 -1.71875 5.796875 -1.703125 5.671875 -1.578125 C 4.75 -0.359375 3.390625 -0.140625 2.84375 -0.140625 C 1.84375 -0.140625 1.53125 -1 1.53125 -1.71875 C 1.53125 -2.21875 1.78125 -3.609375 2.296875 -4.59375 C 2.671875 -5.265625 3.4375 -6.03125 4.34375 -6.03125 C 4.53125 -6.03125 5.328125 -6.015625 5.609375 -5.390625 Z M 5.609375 -5.390625 "/>
</symbol>
<symbol overflow="visible" id="glyph2-0">
<path style="stroke:none;" d=""/>
</symbol>
<symbol overflow="visible" id="glyph2-1">
<path style="stroke:none;" d="M 0.609375 1.828125 C 0.515625 2.25 0.453125 2.359375 -0.125 2.359375 C -0.296875 2.359375 -0.4375 2.359375 -0.4375 2.640625 C -0.4375 2.671875 -0.4375 2.78125 -0.265625 2.78125 C -0.09375 2.78125 0.109375 2.75 0.296875 2.75 L 0.921875 2.75 C 1.21875 2.75 1.953125 2.78125 2.25 2.78125 C 2.34375 2.78125 2.515625 2.78125 2.515625 2.53125 C 2.515625 2.359375 2.40625 2.359375 2.171875 2.359375 C 1.5 2.359375 1.46875 2.265625 1.46875 2.15625 C 1.46875 1.984375 2.109375 -0.484375 2.1875 -0.8125 C 2.359375 -0.421875 2.734375 0.140625 3.484375 0.140625 C 5.109375 0.140625 6.859375 -1.96875 6.859375 -4.078125 C 6.859375 -5.390625 6.109375 -6.328125 5.03125 -6.328125 C 4.109375 -6.328125 3.34375 -5.4375 3.1875 -5.234375 C 3.0625 -5.953125 2.515625 -6.328125 1.9375 -6.328125 C 1.515625 -6.328125 1.1875 -6.125 0.921875 -5.578125 C 0.65625 -5.0625 0.453125 -4.1875 0.453125 -4.125 C 0.453125 -4.078125 0.515625 -4 0.609375 -4 C 0.734375 -4 0.75 -4.015625 0.828125 -4.34375 C 1.046875 -5.1875 1.3125 -6.03125 1.890625 -6.03125 C 2.21875 -6.03125 2.34375 -5.8125 2.34375 -5.375 C 2.34375 -5.03125 2.296875 -4.890625 2.234375 -4.640625 Z M 3.09375 -4.46875 C 3.203125 -4.875 3.59375 -5.296875 3.828125 -5.5 C 3.984375 -5.640625 4.453125 -6.03125 5 -6.03125 C 5.640625 -6.03125 5.921875 -5.40625 5.921875 -4.65625 C 5.921875 -3.96875 5.515625 -2.359375 5.15625 -1.609375 C 4.8125 -0.828125 4.140625 -0.140625 3.484375 -0.140625 C 2.515625 -0.140625 2.359375 -1.375 2.359375 -1.4375 C 2.359375 -1.484375 2.375 -1.59375 2.390625 -1.65625 Z M 3.09375 -4.46875 "/>
</symbol>
<symbol overflow="visible" id="glyph3-0">
<path style="stroke:none;" d=""/>
</symbol>
<symbol overflow="visible" id="glyph3-1">
<path style="stroke:none;" d="M 2.9375 -1.640625 L 2.9375 -0.78125 C 2.9375 -0.421875 2.90625 -0.3125 2.171875 -0.3125 L 1.96875 -0.3125 L 1.96875 0 C 2.375 -0.03125 2.890625 -0.03125 3.3125 -0.03125 C 3.734375 -0.03125 4.25 -0.03125 4.671875 0 L 4.671875 -0.3125 L 4.453125 -0.3125 C 3.71875 -0.3125 3.703125 -0.421875 3.703125 -0.78125 L 3.703125 -1.640625 L 4.6875 -1.640625 L 4.6875 -1.953125 L 3.703125 -1.953125 L 3.703125 -6.484375 C 3.703125 -6.6875 3.703125 -6.75 3.53125 -6.75 C 3.453125 -6.75 3.421875 -6.75 3.34375 -6.625 L 0.28125 -1.953125 L 0.28125 -1.640625 Z M 2.984375 -1.953125 L 0.5625 -1.953125 L 2.984375 -5.671875 Z M 2.984375 -1.953125 "/>
</symbol>
<symbol overflow="visible" id="glyph3-2">
<path style="stroke:none;" d="M 2.890625 -3.515625 C 3.703125 -3.78125 4.28125 -4.46875 4.28125 -5.265625 C 4.28125 -6.078125 3.40625 -6.640625 2.453125 -6.640625 C 1.453125 -6.640625 0.6875 -6.046875 0.6875 -5.28125 C 0.6875 -4.953125 0.90625 -4.765625 1.203125 -4.765625 C 1.5 -4.765625 1.703125 -4.984375 1.703125 -5.28125 C 1.703125 -5.765625 1.234375 -5.765625 1.09375 -5.765625 C 1.390625 -6.265625 2.046875 -6.390625 2.40625 -6.390625 C 2.828125 -6.390625 3.375 -6.171875 3.375 -5.28125 C 3.375 -5.15625 3.34375 -4.578125 3.09375 -4.140625 C 2.796875 -3.65625 2.453125 -3.625 2.203125 -3.625 C 2.125 -3.609375 1.890625 -3.59375 1.8125 -3.59375 C 1.734375 -3.578125 1.671875 -3.5625 1.671875 -3.46875 C 1.671875 -3.359375 1.734375 -3.359375 1.90625 -3.359375 L 2.34375 -3.359375 C 3.15625 -3.359375 3.53125 -2.6875 3.53125 -1.703125 C 3.53125 -0.34375 2.84375 -0.0625 2.40625 -0.0625 C 1.96875 -0.0625 1.21875 -0.234375 0.875 -0.8125 C 1.21875 -0.765625 1.53125 -0.984375 1.53125 -1.359375 C 1.53125 -1.71875 1.265625 -1.921875 0.984375 -1.921875 C 0.734375 -1.921875 0.421875 -1.78125 0.421875 -1.34375 C 0.421875 -0.4375 1.34375 0.21875 2.4375 0.21875 C 3.65625 0.21875 4.5625 -0.6875 4.5625 -1.703125 C 4.5625 -2.515625 3.921875 -3.296875 2.890625 -3.515625 Z M 2.890625 -3.515625 "/>
</symbol>
<symbol overflow="visible" id="glyph3-3">
<path style="stroke:none;" d="M 1.265625 -0.765625 L 2.328125 -1.796875 C 3.875 -3.171875 4.46875 -3.703125 4.46875 -4.703125 C 4.46875 -5.84375 3.578125 -6.640625 2.359375 -6.640625 C 1.234375 -6.640625 0.5 -5.71875 0.5 -4.828125 C 0.5 -4.28125 1 -4.28125 1.03125 -4.28125 C 1.203125 -4.28125 1.546875 -4.390625 1.546875 -4.8125 C 1.546875 -5.0625 1.359375 -5.328125 1.015625 -5.328125 C 0.9375 -5.328125 0.921875 -5.328125 0.890625 -5.3125 C 1.109375 -5.96875 1.65625 -6.328125 2.234375 -6.328125 C 3.140625 -6.328125 3.5625 -5.515625 3.5625 -4.703125 C 3.5625 -3.90625 3.078125 -3.125 2.515625 -2.5 L 0.609375 -0.375 C 0.5 -0.265625 0.5 -0.234375 0.5 0 L 4.203125 0 L 4.46875 -1.734375 L 4.234375 -1.734375 C 4.171875 -1.4375 4.109375 -1 4 -0.84375 C 3.9375 -0.765625 3.28125 -0.765625 3.0625 -0.765625 Z M 1.265625 -0.765625 "/>
</symbol>
<symbol overflow="visible" id="glyph3-4">
<path style="stroke:none;" d="M 2.9375 -6.375 C 2.9375 -6.625 2.9375 -6.640625 2.703125 -6.640625 C 2.078125 -6 1.203125 -6 0.890625 -6 L 0.890625 -5.6875 C 1.09375 -5.6875 1.671875 -5.6875 2.1875 -5.953125 L 2.1875 -0.78125 C 2.1875 -0.421875 2.15625 -0.3125 1.265625 -0.3125 L 0.953125 -0.3125 L 0.953125 0 C 1.296875 -0.03125 2.15625 -0.03125 2.5625 -0.03125 C 2.953125 -0.03125 3.828125 -0.03125 4.171875 0 L 4.171875 -0.3125 L 3.859375 -0.3125 C 2.953125 -0.3125 2.9375 -0.421875 2.9375 -0.78125 Z M 2.9375 -6.375 "/>
</symbol>
</g>
</defs>
<g id="surface10774">
<path style="fill:none;stroke-width:0.8;stroke-linecap:butt;stroke-linejoin:round;stroke:rgb(0%,0%,0%);stroke-opacity:1;stroke-miterlimit:10;" d="M 435.953125 421.160156 C 435.953125 498.929688 372.910156 561.972656 295.144531 561.972656 C 217.375 561.972656 154.332031 498.929688 154.332031 421.160156 C 154.332031 343.394531 217.375 280.351562 295.144531 280.351562 C 372.910156 280.351562 435.953125 343.394531 435.953125 421.160156 Z M 435.953125 421.160156 " transform="matrix(1,0,0,-1,-153,563)"/>
<path style=" stroke:none;fill-rule:evenodd;fill:rgb(0%,0%,0%);fill-opacity:1;" d="M 143.796875 141.847656 C 143.796875 139.449219 140.195312 139.449219 140.195312 141.847656 C 140.195312 144.25 143.796875 144.25 143.796875 141.847656 Z M 143.796875 141.847656 "/>
<g style="fill:rgb(0%,0%,0%);fill-opacity:1;">
<use xlink:href="#glyph0-1" x="283.177808" y="199.223094"/>
</g>
<g style="fill:rgb(0%,0%,0%);fill-opacity:1;">
<use xlink:href="#glyph1-1" x="148.795439" y="151.211264"/>
</g>
<g style="fill:rgb(0%,0%,0%);fill-opacity:1;">
<use xlink:href="#glyph2-1" x="72.319298" y="182.176426"/>
</g>
<g style="fill:rgb(0%,0%,0%);fill-opacity:1;">
<use xlink:href="#glyph3-1" x="79.369298" y="184.328426"/>
</g>
<g style="fill:rgb(0%,0%,0%);fill-opacity:1;">
<use xlink:href="#glyph2-1" x="61.102381" y="97.365032"/>
</g>
<g style="fill:rgb(0%,0%,0%);fill-opacity:1;">
<use xlink:href="#glyph3-2" x="68.152381" y="99.517032"/>
</g>
<g style="fill:rgb(0%,0%,0%);fill-opacity:1;">
<use xlink:href="#glyph2-1" x="232.164693" y="67.862648"/>
</g>
<g style="fill:rgb(0%,0%,0%);fill-opacity:1;">
<use xlink:href="#glyph3-3" x="239.214693" y="70.014648"/>
</g>
<g style="fill:rgb(0%,0%,0%);fill-opacity:1;">
<use xlink:href="#glyph2-1" x="226.685146" y="44.49094"/>
</g>
<g style="fill:rgb(0%,0%,0%);fill-opacity:1;">
<use xlink:href="#glyph3-4" x="233.735146" y="46.64294"/>
</g>
<path style="fill:none;stroke-width:0.8;stroke-linecap:butt;stroke-linejoin:round;stroke:rgb(100%,0%,0%);stroke-opacity:1;stroke-miterlimit:10;" d="M 230.960938 392.6875 C 275.921875 433.171875 325.492188 468.226562 378.644531 497.132812 " transform="matrix(1,0,0,-1,-153,563)"/>
<path style="fill:none;stroke-width:0.8;stroke-linecap:butt;stroke-linejoin:round;stroke:rgb(100%,0%,0%);stroke-opacity:1;stroke-miterlimit:10;" d="M 230.960938 392.6875 C 227.496094 419.308594 227.441406 446.261719 230.796875 472.894531 " transform="matrix(1,0,0,-1,-153,563)"/>
<path style="fill:none;stroke-width:0.8;stroke-linecap:butt;stroke-linejoin:round;stroke:rgb(100%,0%,0%);stroke-opacity:1;stroke-miterlimit:10;" d="M 230.796875 472.894531 C 271 507.796875 323.949219 524.339844 376.871094 518.535156 " transform="matrix(1,0,0,-1,-153,563)"/>
<path style="fill:none;stroke-width:0.8;stroke-linecap:butt;stroke-linejoin:round;stroke:rgb(100%,0%,0%);stroke-opacity:1;stroke-miterlimit:10;" d="M 230.796875 472.894531 C 276.140625 496.8125 328.039062 505.320312 378.644531 497.132812 " transform="matrix(1,0,0,-1,-153,563)"/>
<path style="fill:none;stroke-width:0.8;stroke-linecap:butt;stroke-linejoin:round;stroke:rgb(0%,0%,0%);stroke-opacity:1;stroke-dasharray:4;stroke-miterlimit:10;" d="M 435.957031 420.902344 C 435.066406 470.46875 372.523438 510.960938 295.578125 511.789062 C 218.632812 512.613281 155.835938 473.46875 154.632812 423.925781 " transform="matrix(1,0,0,-1,-153,563)"/>
<path style="fill:none;stroke-width:0.8;stroke-linecap:butt;stroke-linejoin:round;stroke:rgb(0%,100%,0%);stroke-opacity:1;stroke-miterlimit:10;" d="M 379.625 506.074219 C 378.8125 535.113281 340.398438 558.65625 293.824219 558.65625 C 247.246094 558.65625 210.148438 535.113281 210.960938 506.074219 C 211.769531 477.035156 250.183594 453.492188 296.761719 453.492188 C 343.335938 453.492188 380.433594 477.035156 379.625 506.074219 Z M 379.625 506.074219 " transform="matrix(1,0,0,-1,-153,563)"/>
<path style="fill:none;stroke-width:0.8;stroke-linecap:butt;stroke-linejoin:round;stroke:rgb(100%,0%,0%);stroke-opacity:1;stroke-miterlimit:10;" d="M 376.871094 518.535156 C 378.171875 511.476562 378.765625 504.308594 378.644531 497.132812 " transform="matrix(1,0,0,-1,-153,563)"/>
<path style=" stroke:none;fill-rule:evenodd;fill:rgb(0%,0%,0%);fill-opacity:1;" d="M 225.671875 44.464844 C 225.671875 42.0625 222.074219 42.0625 222.074219 44.464844 C 222.074219 46.863281 225.671875 46.863281 225.671875 44.464844 Z M 225.671875 44.464844 "/>
<path style=" stroke:none;fill-rule:evenodd;fill:rgb(0%,0%,0%);fill-opacity:1;" d="M 227.445312 65.867188 C 227.445312 63.46875 223.84375 63.46875 223.84375 65.867188 C 223.84375 68.269531 227.445312 68.269531 227.445312 65.867188 Z M 227.445312 65.867188 "/>
<path style=" stroke:none;fill-rule:evenodd;fill:rgb(0%,0%,0%);fill-opacity:1;" d="M 79.597656 90.105469 C 79.597656 87.707031 75.996094 87.707031 75.996094 90.105469 C 75.996094 92.503906 79.597656 92.503906 79.597656 90.105469 Z M 79.597656 90.105469 "/>
<path style="fill:none;stroke-width:0.8;stroke-linecap:butt;stroke-linejoin:round;stroke:rgb(0%,0%,0%);stroke-opacity:1;stroke-miterlimit:10;" d="M 154.632812 423.925781 C 154.3125 399.796875 168.953125 376.445312 195.304688 359.042969 C 221.65625 341.640625 257.542969 331.628906 295.007812 331.226562 C 332.472656 330.824219 368.417969 340.066406 394.878906 356.898438 C 421.339844 373.734375 436.125 396.769531 435.957031 420.902344 " transform="matrix(1,0,0,-1,-153,563)"/>
<path style=" stroke:none;fill-rule:evenodd;fill:rgb(0%,0%,0%);fill-opacity:1;" d="M 79.761719 170.3125 C 79.761719 167.914062 76.160156 167.914062 76.160156 170.3125 C 76.160156 172.710938 79.761719 172.710938 79.761719 170.3125 Z M 79.761719 170.3125 "/>
</g>
</svg>
Side by Side

After

Width:  |  Height:  |  Size: 14 KiB

BIN
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/fig/geometric_embedding.png Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 488 KiB

BIN
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/fig/ghost_faces.png Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 236 KiB

BIN
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/fig/header.png Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 570 KiB

BIN
Triangulation_on_sphere_2/doc/Triangulation_on_sphere_2/fig/triangulation_on_sphere.png Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 16 KiB

19
Triangulation_on_sphere_2/examples/Triangulation_on_sphere_2/CMakeLists.txt Normal file
View File

@ -0,0 +1,19 @@
cmake_minimum_required(VERSION 3.1...3.15)
project( Triangulation_on_sphere_2_Examples )
find_package(CGAL REQUIRED COMPONENTS Core)
if ( CGAL_FOUND )
create_single_source_cgal_program( "triang_on_sphere.cpp" )
create_single_source_cgal_program( "triang_on_sphere_range.cpp" )
create_single_source_cgal_program( "triang_on_sphere_exact.cpp" )
create_single_source_cgal_program( "triang_on_sphere_proj.cpp" )
create_single_source_cgal_program( "triang_on_sphere_geo.cpp" )
else()
message(STATUS "This program requires the CGAL library, and will not be compiled.")
endif()

9031
Triangulation_on_sphere_2/examples/Triangulation_on_sphere_2/data/poste_france.xyz Normal file
View File

File diff suppressed because it is too large Load Diff

20950
Triangulation_on_sphere_2/examples/Triangulation_on_sphere_2/data/radar.xyz Normal file
View File

File diff suppressed because it is too large Load Diff

46
Triangulation_on_sphere_2/examples/Triangulation_on_sphere_2/triang_on_sphere.cpp Normal file
View File

@ -0,0 +1,46 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Delaunay_triangulation_on_sphere_traits_2.h>
#include <CGAL/Delaunay_triangulation_on_sphere_2.h>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Delaunay_triangulation_on_sphere_traits_2<K> Traits;
typedef CGAL::Delaunay_triangulation_on_sphere_2<Traits> DToS2;
typedef Traits::Point_3 Point_3;
int main(int, char**)
{
std::vector<Point_3> points;
points.emplace_back( 2, 1, 1);
points.emplace_back(-2, 1, 1); // not on the sphere
points.emplace_back( 0, 1, 1);
points.emplace_back( 1, 2, 1);
points.emplace_back( 0, 1, 1); // duplicate of #3
points.emplace_back( 1, 0, 1);
points.emplace_back( 1, 1, 2);
Traits traits(Point_3(1, 1, 1), 1); // sphere center on (1,1,1), with radius 1
DToS2 dtos(traits);
for(const Point_3& pt : points)
{
std::cout << "Inserting (" << pt
<< ") at squared distance " << CGAL::squared_distance(pt, traits.center())
<< " from the center of the sphere; is it on there sphere? "
<< (traits.is_on_sphere(pt) ? "yes" : "no") << std::endl;
dtos.insert(pt);
std::cout << "After insertion, the dimension of the triangulation is: " << dtos.dimension() << "\n";
std::cout << "It has:\n";
std::cout << dtos.number_of_vertices() << " vertices\n";
std::cout << dtos.number_of_edges() << " edges\n";
std::cout << dtos.number_of_faces() << " solid faces\n";
std::cout << dtos.number_of_ghost_faces() << " ghost faces\n" << std::endl;
}
CGAL::write_OFF("result.off", dtos, CGAL::parameters::stream_precision(17));
return EXIT_SUCCESS;
}

68
Triangulation_on_sphere_2/examples/Triangulation_on_sphere_2/triang_on_sphere_exact.cpp Normal file
View File

@ -0,0 +1,68 @@
#include <CGAL/Exact_predicates_exact_constructions_kernel_with_sqrt.h>
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Delaunay_triangulation_on_sphere_2.h>
#include <CGAL/Projection_on_sphere_traits_3.h>
#include <iostream>
#include <fstream>
template <typename Kernel>
void create_triangulation(const char* filename)
{
typedef typename Kernel::FT FT;
typedef CGAL::Projection_on_sphere_traits_3<Kernel> Traits;
typedef CGAL::Delaunay_triangulation_on_sphere_2<Traits> DToS2;
typedef typename Traits::Point_3 Point_3;
std::cout << "\n-- Constructing triangulation with Kernel: " << typeid(Kernel).name() << " --" << std::endl;
std::vector<Point_3> points;
double x, y, z;
std::ifstream in(filename);
if(!in)
{
std::cerr << "Invalid input file: " << filename << std::endl;
return;
}
while(in >> x >> y >> z)
points.emplace_back(x, y, z);
// Add an extra point that would be too close to 'p' with a basic kernel such as CGAL::EPICK,
const Point_3& p = points.back();
const FT tiny = 100 * std::numeric_limits<double>::epsilon();
points.emplace_back(p.x() + tiny, p.y() - tiny, p.z() + tiny);
std::cout << "Adding point " << points.back() << "\nvery close to " << p << std::endl;
std::cout << "Squared distance between points " << CGAL::squared_distance(points.back(), p) << std::endl;
std::cout << points.size() << " points in input" << std::endl;
Traits traits(Point_3(0, 0, 0), 100); // centered on (0,0,0), with radius 100
DToS2 dtos(points.begin(), points.end(), traits);
std::cout << dtos.number_of_vertices() << " vertices" << std::endl;
std::cout << dtos.number_of_faces() << " faces" << std::endl;
}
int main(int argc, char** argv)
{
std::cout.precision(17);
// This kernel CAN represent exactly all points of the sphere
typedef CGAL::Exact_predicates_exact_constructions_kernel_with_sqrt EPECK_w_SQRT;
// This kernel CANNOT represent exactly all points of the sphere
// and thus a separation mechanism is needed to ensure that no points are hidden
typedef CGAL::Exact_predicates_inexact_constructions_kernel EPICK;
const char* filename = (argc > 1) ? argv[1] : "data/poste_france.xyz";
create_triangulation<EPICK>(filename);
create_triangulation<EPECK_w_SQRT>(filename);
return EXIT_SUCCESS;
}

63
Triangulation_on_sphere_2/examples/Triangulation_on_sphere_2/triang_on_sphere_geo.cpp Normal file
View File

@ -0,0 +1,63 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Delaunay_triangulation_on_sphere_2.h>
#include <CGAL/Geographical_coordinates_traits_2.h>
#include <iostream>
#include <fstream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Geographical_coordinates_traits_2<K> Traits;
typedef CGAL::Delaunay_triangulation_on_sphere_2<Traits> DToS2;
typedef Traits::Point_3 Point_3;
typedef Traits::Point_on_sphere_2 Point;
// Note: Geographical_coordinates_traits_2 are currently undocumented because they use
// a conversion to the Cartesian domain R^3 to perform predicates and constructions,
// losing the benefit of the exact representation on the way.
// It can be useful on its own (especially since it contains code
// to handle generic convex domains), but hidden points are possible
// and are not handled properly.
int main(int argc, char** argv)
{
std::cout.precision(17);
const char* filename = (argc > 1) ? argv[1] : "data/poste_france.xyz";
Traits traits(Point_3(0, 0, 0), 100);
std::ifstream in(filename);
if(!in)
{
std::cerr << "Invalid input file: " << filename << std::endl;
return EXIT_FAILURE;
}
typename Traits::Construct_point_on_sphere_2 cps2 = traits.construct_point_on_sphere_2_object();
std::vector<Point> points;
double x, y, z;
while(in >> x >> y >> z)
{
Point_3 cp(x, y, z);
Point ps = cps2(cp);
std::cout << "Cartesian point: " << cp << " Coordinates on the sphere: " << ps << std::endl;
points.push_back(ps);
}
std::cout << points.size() << " points in input" << std::endl;
DToS2 dtos(points.begin(), points.end(), traits);
std::cout << dtos.number_of_vertices() << " vertices" << std::endl;
std::cout << dtos.number_of_faces() << " solid faces" << std::endl;
std::cout << dtos.number_of_ghost_faces() << " ghost faces" << std::endl;
CGAL::write_OFF("result.off", dtos, CGAL::parameters::stream_precision(17));
return EXIT_SUCCESS;
}

46
Triangulation_on_sphere_2/examples/Triangulation_on_sphere_2/triang_on_sphere_proj.cpp Normal file
View File

@ -0,0 +1,46 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Delaunay_triangulation_on_sphere_2.h>
#include <CGAL/Projection_on_sphere_traits_3.h>
#include <boost/iterator/transform_iterator.hpp>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Projection_on_sphere_traits_3<K> Traits;
typedef CGAL::Delaunay_triangulation_on_sphere_2<Traits> DToS2;
typedef Traits::Point_3 Point_3;
int main(int, char**)
{
std::vector<Point_3> points;
points.emplace_back( 3, 1, 1);
points.emplace_back(-8, 1, 1);
points.emplace_back( 1, 2, 1);
points.emplace_back( 1, -2, 1);
points.emplace_back( 1, 1, 10);
Traits traits(Point_3(1,1,1)); // radius is 1 by default
DToS2 dtos(traits);
Traits::Construct_point_on_sphere_2 cst = traits.construct_point_on_sphere_2_object();
for(const auto& pt : points)
{
std::cout << "----- Inserting (" << pt
<< ") at squared distance " << CGAL::squared_distance(pt, traits.center())
<< " from the center of the sphere" << std::endl;
dtos.insert(cst(pt));
std::cout << "The triangulation now has dimension: " << dtos.dimension() << " and\n";
std::cout << dtos.number_of_vertices() << " vertices" << std::endl;
std::cout << dtos.number_of_edges() << " edges" << std::endl;
std::cout << dtos.number_of_faces() << " solid faces" << std::endl;
std::cout << dtos.number_of_ghost_faces() << " ghost faces" << std::endl;
}
CGAL::write_OFF("result.off", dtos, CGAL::parameters::stream_precision(17));
return EXIT_SUCCESS;
}

46
Triangulation_on_sphere_2/examples/Triangulation_on_sphere_2/triang_on_sphere_range.cpp Normal file
View File

@ -0,0 +1,46 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Delaunay_triangulation_on_sphere_2.h>
#include <CGAL/Projection_on_sphere_traits_3.h>
#include <iostream>
#include <fstream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Projection_on_sphere_traits_3<K> Traits;
typedef CGAL::Delaunay_triangulation_on_sphere_2<Traits> DToS2;
typedef DToS2::Point_3 Point;
int main(int argc, char** argv)
{
std::cout.precision(17);
const char* filename = (argc > 1) ? argv[1] : "data/radar.xyz";
std::vector<Point> points;
double x, y, z;
std::ifstream in(filename);
if(!in)
{
std::cerr << "Invalid input file: " << filename << std::endl;
return EXIT_FAILURE;
}
while(in >> x >> y >> z)
points.emplace_back(x, y, z);
std::cout << points.size() << " points in input" << std::endl;
Traits traits(Point(0, 0, 0), 100);
DToS2 dtos(points.begin(), points.end(), traits);
std::cout << dtos.number_of_vertices() << " vertices" << std::endl;
std::cout << dtos.number_of_faces() << " solid faces" << std::endl;
CGAL::write_OFF("result.off", dtos, CGAL::parameters::stream_precision(17));
return EXIT_SUCCESS;
}

1195
Triangulation_on_sphere_2/include/CGAL/Delaunay_triangulation_on_sphere_2.h Normal file
View File

File diff suppressed because it is too large Load Diff

411
Triangulation_on_sphere_2/include/CGAL/Delaunay_triangulation_on_sphere_traits_2.h Normal file
View File

@ -0,0 +1,411 @@
// Copyright (c) 1997, 2012, 2019 INRIA Sophia-Antipolis (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org).
//
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial
//
// Author(s) : Mariette Yvinec,
// Claudia Werner,
// Mael Rouxel-Labbé
#ifndef CGAL_DELAUNAY_TRIANGULATION_SPHERE_TRAITS_2_H
#define CGAL_DELAUNAY_TRIANGULATION_SPHERE_TRAITS_2_H
#include <CGAL/license/Triangulation_on_sphere_2.h>
#include <CGAL/Triangulation_on_sphere_2/internal/get_precision_bounds.h>
#include <CGAL/Algebraic_kernel_for_spheres_2_3.h>
#include <CGAL/Has_conversion.h>
#include <CGAL/Spherical_kernel_3.h>
#include <CGAL/squared_distance_3.h>
#include <CGAL/triangulation_assertions.h>
#include <CGAL/enum.h>
namespace CGAL {
namespace internal {
template <typename LK>
class Orientation_on_sphere_2
{
public:
typedef typename LK::Point_3 Point_3;
typedef typename LK::Point_3 Point_on_sphere_2;
typedef Comparison_result result_type;
Orientation_on_sphere_2(const Point_3& center, const LK& lk)
: _center(center), _lk(lk)
{ }
Comparison_result operator()(const Point_on_sphere_2& p, const Point_on_sphere_2& q, const Point_on_sphere_2& r) const
{ return _lk.orientation_3_object()(_center, p, q, r); }
protected:
const Point_3& _center;
const LK& _lk;
};
template <typename LK>
class Equal_on_sphere_2
{
public:
typedef typename LK::Point_3 Point_3;
typedef typename LK::Point_3 Point_on_sphere_2;
typedef bool result_type;
Equal_on_sphere_2(const Point_3& center, const LK& lk)
: _center(center), _lk(lk)
{ }
bool operator()(const Point_on_sphere_2& p, const Point_on_sphere_2& q) const
{
return _lk.collinear_3_object()(_center, p, q) &&
!_lk.collinear_are_strictly_ordered_along_line_3_object()(p, _center, q);
}
protected:
const Point_3& _center;
const LK& _lk;
};
template <typename LK>
class Inside_cone_2
{
public:
typedef typename LK::Point_3 Point_3;
typedef typename LK::Point_3 Point_on_sphere_2;
typedef bool result_type;
Inside_cone_2(const Point_3& center, const LK& lk)
: _center(center), _lk(lk)
{ }
bool operator()(const Point_on_sphere_2& p, const Point_on_sphere_2& q, const Point_on_sphere_2& r) const
{
// @todo first two checks might be unnecessary because we should always have r != p|q (using equal_on_sphere)
if(_lk.collinear_3_object()(_center, p, r) ||
_lk.collinear_3_object()(_center, q, r) ||
_lk.orientation_3_object()(_center, p, q, r) != COLLINEAR)
{
return false;
}
// @todo also useless?
if(collinear(_center, p, q))
return true;
// Check if r is on the same side of the plane orthogonal to the plane (c,p,q)
// and going through c and p as q
const Orientation op = _lk.coplanar_orientation_3_object()(_center, p, q, r);
CGAL_triangulation_assertion(op != COLLINEAR);
if(op == NEGATIVE)
return false;
// Check if r is on the same side of the plane orthogonal to the plane (c,p,q)
// and going through c and q as p
const Orientation oq = _lk.coplanar_orientation_3_object()(_center, q, p, r);
CGAL_triangulation_assertion(oq != COLLINEAR);
if(oq == NEGATIVE)
return false;
return true; // both true
}
protected:
const Point_3& _center;
const LK& _lk;
};
template <typename LK>
typename LK::Point_3
ray_sphere_intersection(const typename LK::Point_3 o, // origin of the ray
const typename LK::Vector_3 dir, // direction of the ray
const typename LK::Point_3 sc, // center of the sphere
const typename LK::FT r, // radius of the sphere
const LK& lk)
{
typedef typename LK::FT FT;
typedef typename LK::Vector_3 Vector_3;
// Line-Sphere intersection is:
//
// || x - c ||^2 = r^2 // sphere
// x = o + t * d // ray
//
// thus || o + t*d - c ||^2 = r^2
//
// so a*t^2 + b*t + c = 0 with
//
// a := || d ||^2
// b := 2 * d * (o - c)
// c := || o - c||^2 - r^2
const Vector_3 sco = lk.construct_vector_3_object()(sc, o);
const FT a = lk.compute_squared_length_3_object()(dir);
const FT b = 2 * lk.compute_scalar_product_3_object()(dir, sco);
const FT c = lk.compute_squared_length_3_object()(sco) - r*r;
FT d = b*b - 4*a*c;
// by construction d > 0 is expected, but just to be safe
if(d < 0)
d = 0;
const FT dr = CGAL::sqrt(d);
const FT det = FT(0.5) / a;
// care only about the positive value, and we know det > 0 && dr > 0
const FT t0 = det * (-b + dr);
return lk.construct_translated_point_3_object()(o, lk.construct_scaled_vector_3_object()(dir, t0));
}
template <typename LK>
class Construct_circumcenter_on_sphere_2
{
public:
typedef typename LK::FT FT;
typedef typename LK::Point_3 Point_3;
typedef typename LK::Point_3 Point_on_sphere_2;
typedef typename LK::Vector_3 Vector_3;
typedef Point_on_sphere_2 result_type;
Construct_circumcenter_on_sphere_2(const Point_3& center, const FT radius, const LK& lk)
: _center(center), _radius(radius), _lk(lk)
{ }
result_type operator()(const Point_on_sphere_2& p, const Point_on_sphere_2& q, const Point_on_sphere_2& r) const
{
const Point_3 c = _lk.construct_circumcenter_3_object()(p, q, r);
const Vector_3 n = _lk.construct_normal_3_object()(p, q, r);
return ray_sphere_intersection(c, n, _center, _radius, _lk);
}
protected:
const Point_3& _center;
const FT _radius;
const LK& _lk;
};
template <typename LK, typename SK>
class Construct_arc_on_sphere_2
{
public:
typedef typename LK::FT FT;
typedef typename LK::Point_3 Point_3;
typedef typename LK::Point_3 Point_on_sphere_2;
typedef typename SK::Circular_arc_3 result_type;
Construct_arc_on_sphere_2(const Point_3& center, const FT radius, const LK& lk, const SK& sk)
: _center(center), _radius(radius), _lk(lk), _sk(sk)
{ }
result_type operator()(const Point_on_sphere_2& p, const Point_on_sphere_2& q) const
{
typedef typename CGAL::internal::Converter_selector<LK, SK>::type Converter;
Converter cv;
typename SK::Point_3 sc = cv(_center); // @tmp cache that
typename SK::Point_3 sp = cv(p);
typename SK::Point_3 sq = cv(q);
CGAL_precondition(!_sk.collinear_3_object()(sc, sp, sq));
typename SK::Plane_3 pl = _sk.construct_plane_3_object()(sc, sp, sq);
CGAL_assertion(is_valid(pl));
// CircleC3 has an annoying assertion to check pl.has_on(sc) that has no chance to
// pass with a non-exact kernel, so dodging it by passing the plane's normal,
//despite this creating internally another plane...
typename SK::Circle_3 c = _sk.construct_circle_3_object()(
sc, square(_radius), _sk.construct_orthogonal_vector_3_object()(pl));
typename SK::Circular_arc_point_3 cp = _sk.construct_circular_arc_point_3_object()(sp);
typename SK::Circular_arc_point_3 cq = _sk.construct_circular_arc_point_3_object()(sq);
return _sk.construct_circular_arc_3_object()(c, cp, cq);
}
protected:
const Point_3& _center;
const FT _radius;
const LK& _lk;
const SK& _sk;
};
} // namespace internal
template <typename LK_,
typename SK_ = CGAL::Spherical_kernel_3<
LK_, CGAL::Algebraic_kernel_for_spheres_2_3<typename LK_::FT> > >
class Delaunay_triangulation_on_sphere_traits_2
: public LK_, SK_
{
public:
typedef LK_ LK;
typedef SK_ SK;
private:
typedef Delaunay_triangulation_on_sphere_traits_2<LK> Self;
public:
typedef typename LK::FT FT;
typedef typename LK::Point_3 Point_on_sphere_2;
typedef typename LK::Point_3 Point_3;
typedef typename LK::Segment_3 Segment_3;
typedef typename LK::Triangle_3 Triangle_3;
typedef typename SK::Circular_arc_3 Arc_on_sphere_2;
// predicates
typedef internal::Inside_cone_2<LK> Collinear_are_strictly_ordered_on_great_circle_2;
typedef typename LK::Compare_xyz_3 Compare_on_sphere_2;
typedef internal::Equal_on_sphere_2<LK> Equal_on_sphere_2;
typedef internal::Orientation_on_sphere_2<LK> Orientation_on_sphere_2;
typedef typename LK::Orientation_3 Side_of_oriented_circle_on_sphere_2;
// constructions
typedef typename LK::Construct_point_3 Construct_point_on_sphere_2;
typedef typename LK::Construct_point_3 Construct_point_3;
typedef typename LK::Construct_segment_3 Construct_segment_3;
typedef typename LK::Construct_triangle_3 Construct_triangle_3;
// For the hilbert sort
typedef typename LK::Compute_x_3 Compute_x_3;
typedef typename LK::Compute_y_3 Compute_y_3;
typedef typename LK::Compute_z_3 Compute_z_3;
typedef internal::Construct_circumcenter_on_sphere_2<LK> Construct_circumcenter_on_sphere_2;
typedef internal::Construct_arc_on_sphere_2<LK, SK> Construct_arc_on_sphere_2;
typedef typename LK::Construct_circumcenter_3 Construct_circumcenter_3;
protected:
Point_3 _center;
FT _radius;
FT _minDistSquared; // minimal distance of two points to each other
FT _minRadiusSquared; // minimal distance of a point from center of the sphere
FT _maxRadiusSquared; // maximal distance of a point from center of the sphere
static constexpr bool _has_exact_rep =
is_same_or_derived<Field_with_sqrt_tag,
typename Algebraic_structure_traits<FT>::Algebraic_category>::value &&
!std::is_floating_point<FT>::value;
public:
Delaunay_triangulation_on_sphere_traits_2(const Point_3& center = CGAL::ORIGIN,
const FT radius = 1,
const LK& k = LK(),
const SK& sk = SK())
: LK(k), SK(sk), _center(center), _radius(radius)
{
initialize_bounds();
}
public:
const LK& lk() const { return static_cast<const LK&>(*this); }
const SK& sk() const { return static_cast<const SK&>(*this); }
void set_center(const Point_3& center) { _center = center; }
const Point_3& center() const { return _center; }
void set_radius(const FT radius) { _radius = radius; initialize_bounds(); }
FT radius() const { return _radius; }
private:
void initialize_bounds()
{
Triangulations_on_sphere_2::internal::ToS2_precision_bound<FT, _has_exact_rep> bound_getter(_radius);
_minDistSquared = bound_getter.get_squared_min_dist();
_minRadiusSquared = bound_getter.get_squared_min_radius();
_maxRadiusSquared = bound_getter.get_squared_max_radius();
}
public:
Collinear_are_strictly_ordered_on_great_circle_2
collinear_are_strictly_ordered_on_great_circle_2_object() const
{ return Collinear_are_strictly_ordered_on_great_circle_2(center(), lk()); }
Compare_on_sphere_2
compare_on_sphere_2_object() const
{ return lk().compare_xyz_3_object(); }
Equal_on_sphere_2
equal_on_sphere_2_object() const
{ return Equal_on_sphere_2(center(), lk()); }
Orientation_on_sphere_2
orientation_on_sphere_2_object() const
{ return Orientation_on_sphere_2(center(), lk()); }
Side_of_oriented_circle_on_sphere_2
side_of_oriented_circle_on_sphere_2_object() const
{ return lk().orientation_3_object(); }
public:
Compute_x_3
compute_x_3_object() const
{ return lk().compute_x_3_object(); }
Compute_y_3
compute_y_3_object() const
{ return lk().compute_y_3_object(); }
Compute_z_3
compute_z_3_object() const
{ return lk().compute_z_3_object(); }
Construct_point_on_sphere_2
construct_point_on_sphere_2_object() const
{ return lk().construct_point_on_sphere_2_object(); }
Construct_point_3
construct_point_3_object() const
{ return lk().construct_point_3_object(); }
Construct_segment_3
construct_segment_3_object() const
{ return lk().construct_segment_3_object(); }
Construct_triangle_3
construct_triangle_3_object() const
{ return lk().construct_triangle_3_object(); }
Construct_arc_on_sphere_2
construct_arc_on_sphere_2_object() const
{ return Construct_arc_on_sphere_2(center(), radius(), lk(), sk()); }
Construct_circumcenter_on_sphere_2
construct_circumcenter_on_sphere_2_object() const
{ return Construct_circumcenter_on_sphere_2(center(), radius(), lk()); }
Construct_circumcenter_3
construct_circumcenter_3_object() const
{ return lk().construct_circumcenter_3_object(); }
public:
bool is_on_sphere(const Point_on_sphere_2& p) const
{
const FT sq_dist = lk().compute_squared_distance_3_object()(p, _center);
if(_has_exact_rep)
return (sq_dist == _minRadiusSquared);
else
return (_minRadiusSquared <= sq_dist && sq_dist < _maxRadiusSquared);
}
bool are_points_too_close(const Point_on_sphere_2& p, const Point_on_sphere_2& q) const
{
if(_has_exact_rep)
return false;
else
return (CGAL::squared_distance(p, q) <= _minDistSquared);
}
};
} // namespace CGAL
#endif // CGAL_DELAUNAY_TRIANGULATION_SPHERE_TRAITS_2_H

307
Triangulation_on_sphere_2/include/CGAL/Geographical_coordinates_traits_2.h Normal file
View File

@ -0,0 +1,307 @@
// Copyright (c) 1997, 2012, 2019 INRIA Sophia-Antipolis (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org).
//
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial
//
// Author(s) : Mael Rouxel-Labbé
#ifndef CGAL_TRIANGULATION_ON_SPHERE_GEOGRAPHICAL_COORDINATES_TRAITS_2_H
#define CGAL_TRIANGULATION_ON_SPHERE_GEOGRAPHICAL_COORDINATES_TRAITS_2_H
#include <CGAL/license/Triangulation_on_sphere_2.h>
#include <CGAL/Delaunay_triangulation_on_sphere_traits_2.h>
#include <CGAL/triangulation_assertions.h>
#include <CGAL/number_utils_classes.h>
#include <cmath>
#include <fstream>
#include <iostream>
// @todo
// Note: this is currently undocumented because it uses a conversion to
// the Cartesian domain R^3 to perform predicates and constructions,
// losing the benefit of the exact representation on the way.
// It can be useful on its own (especially since it contains code
// to handle generic convex domains), but hidden points are possible
// and are not handled properly.
namespace CGAL {
template <typename LK>
struct Geographical_coordinates
{
typedef typename LK::FT FT;
typedef FT Latitude;
typedef FT Longitude;
Geographical_coordinates() : _la(360), _lo(360) { }
Geographical_coordinates(const FT la, const FT lo) : _la(la), _lo(lo) { }
Latitude latitude() const { return _la; }
Longitude longitude() const { return _lo; }
friend std::ostream& operator<<(std::ostream& os, const Geographical_coordinates& l)
{
os << l._la << " " << l._lo;
return os;
}
private:
FT _la;
FT _lo;
};
namespace internal {
template <typename LK>
struct Construct_geographical_coordinates
: public CGAL::cpp98::unary_function<typename LK::Point_3, CGAL::Geographical_coordinates<LK> >
{
typedef typename LK::FT FT;
typedef typename LK::Point_3 Point_3;
typedef CGAL::Geographical_coordinates<LK> result_type;
Construct_geographical_coordinates(const Point_3& center, const FT radius)
: _center(center), _radius(radius)
{ }
result_type operator()(const Point_3& pt) const
{
CGAL_assertion(pt != _center);
const FT la = 180. * std::asin((pt.z() - _center.z()) / _radius) / CGAL_PI;
const FT lo = 180. * std::atan2((pt.y() - _center.y()), (pt.x() - _center.x())) / CGAL_PI;
return result_type(la, lo);
}
private:
const Point_3& _center;
const FT _radius;
};
template <typename LK>
struct Construct_Cartesian_coordinates
: public CGAL::cpp98::unary_function<typename LK::Point_3, CGAL::Geographical_coordinates<LK> >
{
typedef typename LK::FT FT;
typedef typename LK::Point_3 Point_3;
typedef CGAL::Geographical_coordinates<LK> Coordinates;
Construct_Cartesian_coordinates(const Point_3& center, const FT radius)
: _center(center), _radius(radius)
{ }
Point_3 operator()(const Coordinates& l) const
{
#if 1
const FT rla = CGAL_PI * l.latitude() / 180.;
const FT rlo = CGAL_PI * l.longitude() / 180.;
const FT x = _center.x() + _radius * std::cos(rla) * std::cos(rlo);
const FT y = _center.y() + _radius * std::cos(rla) * std::sin(rlo);
const FT z = _center.z() + _radius * std::sin(rla);
#else
const FT rla = CGAL_PI * l.latitude() / 180.;
const FT rlo = CGAL_PI * l.longitude() / 180.;
// Quoting wikipedia:
// - The Earth's equatorial radius a, or semi-major axis, is the distance from its center
// to the equator and equals 6,378.1370 km (3,963.1906 mi)
// - The Earth's polar radius b, or semi-minor axis, is the distance from its center to the North
// and South Poles, and equals 6,356.7523 km (3,949.9028 mi).
const FT equ_radius = 6378.1370;
const FT pol_radius = 6356.7523;
FT radius_ratio = equ_radius / pol_radius;
const FT a2 = square(_radius * radius_ratio);
const FT b2 = square(_radius);
// Using https://en.wikipedia.org/wiki/Geographic_coordinate_conversion#From_geodetic_to_ECEF_coordinates
const FT N = a2 / sqrt(a2 * square(std::cos(rla)) + b2 * square(std::sin(rla)));
const FT h = 0;
const FT x = _center.x() + (N + h) * std::cos(rla) * std::cos(rlo);
const FT y = _center.y() + (N + h) * std::cos(rla) * std::sin(rlo);
const FT z = _center.z() + (b2 * N / a2 + h) * std::sin(rla);
#endif
return Point_3(x, y, z);
}
private:
const Point_3& _center;
const FT _radius;
};
template <typename LK, typename Base_functor_>
struct Construct_circumcenter_on_geo_sphere_2
{
public:
typedef typename LK::FT FT;
typedef typename LK::Point_3 Point_3;
typedef typename LK::Vector_3 Vector_3;
typedef Geographical_coordinates<LK> Point_on_sphere_2;
typedef Geographical_coordinates<LK> result_type;
Construct_circumcenter_on_geo_sphere_2( const Base_functor_ f, const Point_3& center, const FT radius)
: _center(center), _radius(radius), _f(f)
{ }
result_type operator()(const Point_on_sphere_2& p, const Point_on_sphere_2& q, const Point_on_sphere_2& r) const
{
Construct_Cartesian_coordinates<LK> to_P3(_center, _radius);
Construct_geographical_coordinates<LK> to_PoS2(_center, _radius);
return to_PoS2(_f(to_P3(p), to_P3(q), to_P3(r)));
}
private:
const Base_functor_ _f;
const Point_3& _center;
const FT _radius;
};
// Adaptor for calling the functors with the points projected on the sphere
template <typename LK, typename Functor_>
class Functor_geo_projection_adaptor
: public Functor_
{
typedef Functor_ Functor;
typedef Functor_ Base;
typedef typename LK::FT FT;
typedef typename LK::Point_3 Point_3;
typedef Geographical_coordinates<LK> Point;
typedef Construct_Cartesian_coordinates<LK> To_Point_3;
public:
Functor_geo_projection_adaptor(const Functor& f, const Point_3& center, const FT radius)
: Base(f), to_p3(center, radius)
{ }
public:
using Base::operator();
decltype(auto) operator()(const Point& p0, const Point& p1)
{ return Base::operator()(to_p3(p0), to_p3(p1)); }
decltype(auto) operator()(const Point& p0, const Point& p1, const Point& p2)
{ return Base::operator()(to_p3(p0), to_p3(p1), to_p3(p2)); }
decltype(auto) operator ()(const Point& p0, const Point& p1, const Point& p2, const Point& p3)
{ return Base::operator()(to_p3(p0), to_p3(p1), to_p3(p2), to_p3(p3)); }
decltype(auto) operator()(const Point& p0, const Point& p1, const Point& p2, const Point& p3, const Point& p4)
{ return Base::operator()(to_p3(p0), to_p3(p1), to_p3(p2), to_p3(p3), to_p3(p4)); }
private:
const To_Point_3 to_p3;
};
} // namespace internal
template <typename LK>
class Geographical_coordinates_traits_2
: public Delaunay_triangulation_on_sphere_traits_2<LK>
{
protected:
typedef Delaunay_triangulation_on_sphere_traits_2<LK> Base;
typedef Geographical_coordinates_traits_2<LK> Self;
public:
typedef typename LK::FT FT;
typedef typename LK::Point_3 Point_3;
typedef Geographical_coordinates<LK> Point_on_sphere_2;
// predicates
typedef internal::Functor_geo_projection_adaptor<LK, typename Base::Collinear_are_strictly_ordered_on_great_circle_2> Collinear_are_strictly_ordered_on_great_circle_2;
typedef internal::Functor_geo_projection_adaptor<LK, typename Base::Compare_on_sphere_2> Compare_on_sphere_2;
typedef internal::Functor_geo_projection_adaptor<LK, typename Base::Equal_on_sphere_2> Equal_on_sphere_2;
typedef internal::Functor_geo_projection_adaptor<LK, typename Base::Orientation_on_sphere_2> Orientation_on_sphere_2;
typedef internal::Functor_geo_projection_adaptor<LK, typename Base::Side_of_oriented_circle_on_sphere_2> Side_of_oriented_circle_on_sphere_2;
// constructions
typedef internal::Construct_geographical_coordinates<LK> Construct_point_on_sphere_2;
typedef internal::Construct_Cartesian_coordinates<LK> Construct_point_3;
typedef internal::Functor_geo_projection_adaptor<LK, typename Base::Construct_segment_3> Construct_segment_3;
typedef internal::Functor_geo_projection_adaptor<LK, typename Base::Construct_triangle_3> Construct_triangle_3;
typedef internal::Functor_geo_projection_adaptor<LK, typename Base::Construct_arc_on_sphere_2> Construct_arc_on_sphere_2;
typedef internal::Construct_circumcenter_on_geo_sphere_2<LK, typename Base::Construct_circumcenter_on_sphere_2> Construct_circumcenter_on_sphere_2;
public:
Geographical_coordinates_traits_2(const Point_3& sphere = CGAL::ORIGIN,
const FT radius = 1,
const LK& lk = LK())
: Base(sphere, radius, lk)
{ }
// predicates
public:
Collinear_are_strictly_ordered_on_great_circle_2
collinear_are_strictly_ordered_on_great_circle_2_object() const
{ return Collinear_are_strictly_ordered_on_great_circle_2(
Base::collinear_are_strictly_ordered_on_great_circle_2_object(),
Base::center(), Base::radius()); }
Compare_on_sphere_2
compare_on_sphere_2_object() const
{ return Compare_on_sphere_2(Base::compare_on_sphere_2_object(),
Base::center(), Base::radius()); }
Equal_on_sphere_2
equal_on_sphere_2_object() const
{ return Equal_on_sphere_2(Base::equal_on_sphere_2_object(),
Base::center(), Base::radius()); }
Orientation_on_sphere_2
orientation_on_sphere_2_object() const
{ return Orientation_on_sphere_2(Base::orientation_on_sphere_2_object(),
Base::center(), Base::radius()); }
Side_of_oriented_circle_on_sphere_2
side_of_oriented_circle_on_sphere_2_object() const
{ return Side_of_oriented_circle_on_sphere_2(Base::side_of_oriented_circle_on_sphere_2_object(),
Base::center(), Base::radius()); }
// constructions
public:
Construct_point_on_sphere_2
construct_point_on_sphere_2_object() const
{ return Construct_point_on_sphere_2(Base::center(), Base::radius()); }
Construct_point_3
construct_point_3_object() const
{ return Construct_point_3(Base::center(), Base::radius()); }
Construct_arc_on_sphere_2
construct_arc_on_sphere_2_object() const
{ return Construct_arc_on_sphere_2(Base::construct_arc_on_sphere_2_object(),
Base::center(), Base::radius()); }
Construct_circumcenter_on_sphere_2
construct_circumcenter_on_sphere_2_object() const
{ return Construct_circumcenter_on_sphere_2(Base::construct_circumcenter_on_sphere_2_object(),
Base::center(), Base::radius()); }
public:
bool is_on_sphere(const Point_on_sphere_2& /*p*/) const
{
return true;
}
bool are_points_too_close(const Point_on_sphere_2& p, const Point_on_sphere_2& q) const
{
return Base::are_points_too_close(construct_point_3_object()(p), construct_point_3_object()(q));
}
};
} // namespace CGAL
#endif // CGAL_TRIANGULATION_ON_SPHERE_GEOGRAPHICAL_COORDINATES_TRAITS_2_H

236
Triangulation_on_sphere_2/include/CGAL/Projection_on_sphere_traits_3.h Normal file
View File

@ -0,0 +1,236 @@
// Copyright (c) 1997, 2012, 2019 INRIA Sophia-Antipolis (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org).
//
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial
//
// Author(s) : Mariette Yvinec,
// Claudia Werner
// Mael Rouxel-Labbé
#ifndef CGAL_TRIANGULATION_ON_SPHERE_PROJECTION_SPHERE_TRAITS_3_H
#define CGAL_TRIANGULATION_ON_SPHERE_PROJECTION_SPHERE_TRAITS_3_H
#include <CGAL/license/Triangulation_on_sphere_2.h>
#include <CGAL/Delaunay_triangulation_on_sphere_traits_2.h>
#include <CGAL/triangulation_assertions.h>
#include <CGAL/number_utils_classes.h>
namespace CGAL {
namespace internal {
template <typename LK>
class Point_with_scale
: public LK::Point_3
{
public:
typedef typename LK::FT FT;
typedef typename LK::Point_3 Base_point;
typedef typename LK::Vector_3 Vector_3;
private:
mutable bool cached;
mutable Base_point proj_pt;
public:
Point_with_scale() : Base_point() { } // vertex base wants a default constructor
Point_with_scale(const Base_point& p) : Base_point(p), cached(false) { }
const Base_point& get_projection(const Base_point& center,
const FT radius,
const LK& lk) const
{
if(!cached)
compute_projection(center, radius, lk);
return proj_pt;
}
const Base_point& get_projection(const Base_point& center,
const FT radius) const
{
return get_projection(center, radius, LK());
}
private:
void compute_projection(const Base_point& center,
const FT radius,
const LK& lk) const
{
CGAL_precondition(static_cast<const Base_point&>(*this) != center);
Vector_3 v = lk.construct_vector_3_object()(center, static_cast<const Base_point&>(*this));
FT n = lk.compute_squared_length_3_object()(v);
CGAL_assertion(n > 0);
n = CGAL::sqrt(n);
v = lk.construct_scaled_vector_3_object()(v, radius / n);
proj_pt = lk.construct_translated_point_3_object()(center, v);
cached = true;
}
};
template <typename P3, typename PoS2>
struct Construct_point_with_scale
: public CGAL::cpp98::unary_function<P3, PoS2>
{
PoS2 operator()(const P3& pt) const { return PoS2(pt); }
};
// Adaptor for calling the functors with the points projected on the sphere
template <typename LK, typename Functor_>
class Functor_projection_adaptor
: public Functor_
{
typedef Functor_ Functor;
typedef Functor_ Base;
typedef typename LK::FT FT;
typedef typename LK::Point_3 Point_3;
typedef internal::Point_with_scale<LK> Point;
public:
Functor_projection_adaptor(const Functor& f, const Point_3& center, const FT radius)
: Base(f), _c(center), _r(radius)
{ }
public:
using Base::operator();
decltype(auto) operator()(const Point& p0, const Point& p1)
{
return Base::operator()(p0.get_projection(_c, _r), p1.get_projection(_c, _r));
}
decltype(auto) operator()(const Point& p0, const Point& p1, const Point& p2)
{
return Base::operator()(p0.get_projection(_c, _r), p1.get_projection(_c, _r),
p2.get_projection(_c, _r));
}
decltype(auto) operator ()(const Point& p0, const Point& p1, const Point& p2, const Point& p3)
{
return Base::operator()(p0.get_projection(_c, _r), p1.get_projection(_c, _r),
p2.get_projection(_c, _r), p3.get_projection(_c, _r));
}
decltype(auto) operator()(const Point& p0, const Point& p1, const Point& p2,
const Point& p3, const Point& p4)
{
return Base::operator()(p0.get_projection(_c, _r), p1.get_projection(_c, _r),
p2.get_projection(_c, _r), p3.get_projection(_c, _r),
p4.get_projection(_c, _r));
}
private:
const Point_3& _c;
const FT _r;
};
} // namespace internal
template <typename LK>
class Projection_on_sphere_traits_3
: public Delaunay_triangulation_on_sphere_traits_2<LK>
{
protected:
typedef Delaunay_triangulation_on_sphere_traits_2<LK> Base;
typedef Projection_on_sphere_traits_3<LK> Self;
public:
typedef typename LK::FT FT;
typedef typename LK::Point_3 Point_3;
typedef internal::Point_with_scale<LK> Point_on_sphere_2;
// predicates
typedef internal::Functor_projection_adaptor<LK, typename Base::Collinear_are_strictly_ordered_on_great_circle_2> Collinear_are_strictly_ordered_on_great_circle_2;
typedef internal::Functor_projection_adaptor<LK, typename Base::Compare_on_sphere_2> Compare_on_sphere_2;
typedef internal::Functor_projection_adaptor<LK, typename Base::Equal_on_sphere_2> Equal_on_sphere_2;
typedef internal::Functor_projection_adaptor<LK, typename Base::Orientation_on_sphere_2> Orientation_on_sphere_2;
typedef internal::Functor_projection_adaptor<LK, typename Base::Side_of_oriented_circle_on_sphere_2> Side_of_oriented_circle_on_sphere_2;
// constructions
typedef internal::Construct_point_with_scale<Point_3, Point_on_sphere_2> Construct_point_on_sphere_2;
typedef internal::Functor_projection_adaptor<LK, typename Base::Construct_point_3> Construct_point_3;
typedef internal::Functor_projection_adaptor<LK, typename Base::Construct_segment_3> Construct_segment_3;
typedef internal::Functor_projection_adaptor<LK, typename Base::Construct_triangle_3> Construct_triangle_3;
typedef internal::Functor_projection_adaptor<LK, typename Base::Construct_arc_on_sphere_2> Construct_arc_on_sphere_2;
typedef internal::Functor_projection_adaptor<LK, typename Base::Construct_circumcenter_on_sphere_2> Construct_circumcenter_on_sphere_2;
public:
Projection_on_sphere_traits_3(const Point_3& sphere = CGAL::ORIGIN,
const FT radius = 1,
const LK& lk = LK())
: Base(sphere, radius, lk)
{ }
// predicates
public:
Collinear_are_strictly_ordered_on_great_circle_2
collinear_are_strictly_ordered_on_great_circle_2_object() const
{ return Collinear_are_strictly_ordered_on_great_circle_2(
Base::collinear_are_strictly_ordered_on_great_circle_2_object(),
Base::center(), Base::radius()); }
Compare_on_sphere_2
compare_on_sphere_2_object() const
{ return Compare_on_sphere_2(Base::compare_on_sphere_2_object(),
Base::center(), Base::radius()); }
Equal_on_sphere_2
equal_on_sphere_2_object() const
{ return Equal_on_sphere_2(Base::equal_on_sphere_2_object(),
Base::center(), Base::radius()); }
Orientation_on_sphere_2
orientation_on_sphere_2_object() const
{ return Orientation_on_sphere_2(Base::orientation_on_sphere_2_object(),
Base::center(), Base::radius()); }
Side_of_oriented_circle_on_sphere_2
side_of_oriented_circle_on_sphere_2_object() const
{ return Side_of_oriented_circle_on_sphere_2(Base::side_of_oriented_circle_on_sphere_2_object(),
Base::center(), Base::radius()); }
// constructions
public:
Construct_point_on_sphere_2
construct_point_on_sphere_2_object() const
{ return Construct_point_on_sphere_2(); }
Construct_arc_on_sphere_2
construct_arc_on_sphere_2_object() const
{ return Construct_arc_on_sphere_2(Base::construct_arc_on_sphere_2_object(),
Base::center(), Base::radius()); }
Construct_circumcenter_on_sphere_2
construct_circumcenter_on_sphere_2_object() const
{ return Construct_circumcenter_on_sphere_2(Base::construct_circumcenter_on_sphere_2_object(),
Base::center(), Base::radius()); }
Construct_point_3
construct_point_3_object() const
{ return Construct_point_3(Base::construct_point_3_object(), Base::center(), Base::radius()); }
public:
bool is_on_sphere(const Point_on_sphere_2& /*p*/) const
{
return true;
}
bool are_points_too_close(const Point_on_sphere_2& p, const Point_on_sphere_2& q) const
{
return Base::are_points_too_close(p.get_projection(Base::center(), Base::radius(), Base::lk()),
q.get_projection(Base::center(), Base::radius(), Base::lk()));
}
};
} // namespace CGAL
#endif // CGAL_TRIANGULATION_ON_SPHERE_PROJECTION_SPHERE_TRAITS_3_H

1321
Triangulation_on_sphere_2/include/CGAL/Triangulation_on_sphere_2.h Normal file
View File

File diff suppressed because it is too large Load Diff

192
Triangulation_on_sphere_2/include/CGAL/Triangulation_on_sphere_2/IO/OFF.h Normal file
View File

@ -0,0 +1,192 @@
// Copyright (c) 1997-2021 INRIA Sophia-Antipolis (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org)
//
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial
//
// Author(s) : Mael Rouxel-Labbé
#ifndef CGAL_TOS2_IO_OFF_H
#define CGAL_TOS2_IO_OFF_H
#include <CGAL/license/Triangulation_on_sphere_2.h>
#include <CGAL/triangulation_assertions.h>
#include <CGAL/boost/graph/Named_function_parameters.h>
#include <CGAL/boost/graph/named_params_helper.h>
#include <CGAL/IO/helpers.h>
#include <fstream>
#include <string>
#include <unordered_map>
#ifdef DOXYGEN_RUNNING
#define CGAL_BGL_NP_TEMPLATE_PARAMETERS NamedParameters
#define CGAL_BGL_NP_CLASS NamedParameters
#define CGAL_DEPRECATED
#endif
namespace CGAL {
template <typename Gt, typename Tds>
class Triangulation_on_sphere_2;
////////////////////////////////////////////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////////////////////////////////////////////
// Read
////////////////////////////////////////////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////////////////////////////////////////////
// Write
/*!
\ingroup PkgPointSet3IOOFF
\brief writes the content of a triangulation on the sphere into an output stream in the \ref IOStreamOFF.
\tparam Gt the geometric traits type of the triangulation
\tparam Tds the triangulation data structure type of the triangulation
\tparam NamedParameters a sequence of \ref bgl_namedparameters "Named Parameters"
\param os the output stream
\param dt the triangulation
\param np an optional sequence of \ref bgl_namedparameters "Named Parameters" among the ones listed below
\cgalNamedParamsBegin
\cgalParamNBegin{stream_precision}
\cgalParamDescription{a parameter used to set the precision (i.e., how many digits are generated) of the output stream}
\cgalParamType{int}
\cgalParamDefault{the precision of the stream `os`}
\cgalParamNEnd
\cgalNamedParamsEnd
\return `true` if the writing was successful, `false` otherwise.
*/
template <typename Gt, typename Tds, typename CGAL_BGL_NP_TEMPLATE_PARAMETERS>
bool write_OFF(std::ostream& os,
const CGAL::Triangulation_on_sphere_2<Gt, Tds>& dt,
const CGAL_BGL_NP_CLASS& np)
{
typedef Triangulation_on_sphere_2<Gt,Tds> Tr;
typedef typename Tr::Vertex_handle Vertex_handle;
typedef typename Tr::Vertices_iterator Vertex_iterator;
typedef typename Tr::All_faces_iterator Face_iterator;
typedef typename Tr::size_type size_type;
set_stream_precision_from_NP(os, np);
const bool use_colors = parameters::choose_parameter(parameters::get_parameter(np, internal_np::output_color), true);
typename Gt::Construct_point_3 cp3 = dt.geom_traits().construct_point_3_object();
const size_type n = dt.number_of_vertices();
std::stringstream output;
set_stream_precision_from_NP(output, np);
// write the vertices
std::unordered_map<Vertex_handle, size_type> index_of_vertex;
size_type i = 0;
for(Vertex_iterator it = dt.vertices_begin(); it != dt.vertices_end(); ++it, ++i)
{
Vertex_handle vh = it;
output << cp3(vh->point()) << " 0 0 0 \n"; // '0 0 0' for colors
index_of_vertex[vh] = i;
}
CGAL_triangulation_assertion(i == n);
size_type number_of_triangles = 0;
for(Face_iterator fit = dt.all_faces_begin() ; fit != dt.all_faces_end() ; ++fit)
{
output << "3 "
<< index_of_vertex[fit->vertex(0)] << " "
<< index_of_vertex[fit->vertex(1)] << " "
<< index_of_vertex[fit->vertex(2)];
if(use_colors)
{
if(fit->is_ghost())
output << " 229 117 0\n";
else
output << " 0 117 229\n";
}
else
{
output << "\n";
}
++number_of_triangles;
}
if(use_colors)
os << "COFF \n";
else
os << "OFF \n";
os << n << " "
<< number_of_triangles << " 0\n"
<< output.str();
return !os.fail();
}
/// \cond SKIP_IN_MANUAL
template <typename Gt, typename Tds>
bool write_OFF(std::ostream& os, const CGAL::Triangulation_on_sphere_2<Gt, Tds>& dt)
{
return write_OFF(os, dt, parameters::all_default());
}
/// \endcond
/*!
\ingroup PkgPointSet3IOOFF
\brief writes the content of a triangulation on the sphere into an output file in the \ref IOStreamOFF.
\tparam Gt the geometric traits type of the triangulation
\tparam Tds the triangulation data structure type of the triangulation
\tparam NamedParameters a sequence of \ref bgl_namedparameters "Named Parameters"
\param fname the path to the output file
\param point_set the point set
\param np an optional sequence of \ref bgl_namedparameters "Named Parameters" among the ones listed below
\cgalNamedParamsBegin
\cgalParamNBegin{stream_precision}
\cgalParamDescription{a parameter used to set the precision (i.e., how many digits are generated) of the output stream}
\cgalParamType{int}
\cgalParamDefault{`6`}
\cgalParamNEnd
\cgalNamedParamsEnd
\return `true` if the writing was successful, `false` otherwise.
*/
template <typename Gt, typename Tds, typename CGAL_BGL_NP_TEMPLATE_PARAMETERS>
bool write_OFF(const std::string& fname,
const CGAL::Triangulation_on_sphere_2<Gt, Tds>& dt,
const CGAL_BGL_NP_CLASS& np)
{
std::ofstream os(fname); // stream precision will be set in the ostream overload
return write_OFF(os, dt, np);
}
/// \cond SKIP_IN_MANUAL
template <typename Gt, typename Tds>
bool write_OFF(const std::string& fname, const CGAL::Triangulation_on_sphere_2<Gt, Tds>& dt)
{
std::ofstream os(fname);
return write_OFF(os, dt, parameters::all_default());
}
/// \endcond
} // namespace CGAL
#endif // CGAL_TOS2_IO_OFF_H

159
Triangulation_on_sphere_2/include/CGAL/Triangulation_on_sphere_2/internal/arc_on_sphere_2_subsampling.h Normal file
View File

@ -0,0 +1,159 @@
// Copyright (c) 1997-2021 INRIA Sophia-Antipolis (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org).
//
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial
//
// Author(s) : Sébastien Loriot
// Mael Rouxel-Labbé
//
#ifndef CGAL_TOS2_INTERNAL_ARC_ON_SPHERE_SUBSAMPLING_H
#define CGAL_TOS2_INTERNAL_ARC_ON_SPHERE_SUBSAMPLING_H
#include <CGAL/license/Triangulation_on_sphere_2.h>
#ifdef CGAL_EIGEN3_ENABLED
#include <CGAL/Eigen_solver_traits.h>
#include <cmath>
#include <list>
namespace CGAL {
namespace Triangulations_on_sphere_2 {
namespace internal {
template <class Kernel>
double get_theta( typename Kernel::Point_3& pt,
typename Kernel::Vector_3& V1,
typename Kernel::Vector_3& V2,
typename Kernel::Vector_3& V3)
{
typedef typename Kernel::FT FT;
typedef Eigen::Matrix<FT, 3, 3, Eigen::DontAlign> Matrix;
typedef Eigen::Matrix<FT, 3, 1> Col;
auto V1c = V1.cartesian_begin(), V2c = V2.cartesian_begin(), V3c = V3.cartesian_begin();
Matrix M;
M(0, 0) = *V1c++; M(0, 1) = *V2c++; M(0, 2) = *V3c++;
M(1, 0) = *V1c++; M(1, 1) = *V2c++; M(1, 2) = *V3c++;
M(2, 0) = *V1c++; M(2, 1) = *V2c++; M(2, 2) = *V3c++;
Col pt_in_vect;
auto ptc = pt.cartesian_begin();
pt_in_vect(0) = *ptc++; pt_in_vect(1) = *ptc++; pt_in_vect(2) = *ptc++;
Matrix inverse = M.inverse();
FT X = inverse.row(0) * pt_in_vect;
FT Y = inverse.row(1) * pt_in_vect;
double angle = atan2(Y*sqrt(V2.squared_length()), X*sqrt(V1.squared_length()));
return angle;
}
template <class Kernel>
inline typename Kernel::Point_3 compute_point(const typename Kernel::Point_3& center,
double radius,
double theta,
const typename Kernel::Vector_3& b1,
const typename Kernel::Vector_3& b2)
{
return center +
radius * cos(theta) * b1 / sqrt(b1.squared_length()) +
radius * sin(theta) * b2 / sqrt(b2.squared_length());
}
template <class Kernel, class Output_iterator>
void subsample_arc_on_sphere_2(const typename Kernel::Circle_3& circle,
double source,
double target,
const typename Kernel::Vector_3& b1,
const typename Kernel::Vector_3& b2,
Output_iterator out_pts,
double min_edge_size)
{
if(source > target)
target += 2*CGAL_PI;
CGAL_triangulation_assertion(target > source);
const double radius = sqrt(circle.squared_radius());
const double edge_len = (target - source) * radius;
const int nb_of_segments = static_cast<int>(floor(edge_len / min_edge_size));
*out_pts++ = compute_point<Kernel>(circle.center(), radius, source, b1, b2);
const double step_size = (target - source) / static_cast<double>(nb_of_segments);
double current_theta = source;
for(int i=0; i<nb_of_segments-1; ++i)
{
current_theta += step_size;
CGAL_triangulation_assertion(current_theta <= target);
*out_pts++ = compute_point<Kernel>(circle.center(), radius, current_theta, b1, b2);
}
*out_pts++ = compute_point<Kernel>(circle.center(), radius, target, b1, b2);
}
//Subsample from source to target seen ccw from the side of plane pointed by its orthogonal_vector()
template <class Kernel,class Output_iterator>
void subsample_circle_3(const typename Kernel::Circle_3& circle,
Output_iterator out_pts,
double min_edge_size)
{
typename Kernel::Vector_3 b1 = circle.supporting_plane().base1();
typename Kernel::Vector_3 b2 = circle.supporting_plane().base2();
subsample_arc_on_sphere_2<Kernel>(circle, 0, 2*CGAL_PI, b1, b2, out_pts, min_edge_size);
}
//Subsample from source to target seen ccw from the side of plane pointed by its orthogonal_vector()
template <class Kernel,class Output_iterator>
void subsample_arc_on_sphere_2(const typename Kernel::Circle_3& circle,
const typename Kernel::Plane_3& plane,
const typename Kernel::Circular_arc_point_3& source,
const typename Kernel::Circular_arc_point_3& target,
Output_iterator out_pts,
double min_edge_size)
{
typename Kernel::Vector_3 b1 = plane.base1();
typename Kernel::Vector_3 b2 = plane.base2();
typename Kernel::Vector_3 b3 = plane.orthogonal_vector();
typename Kernel::Point_3 tmp = typename Kernel::Point_3(source.x() - circle.center().x(),
source.y() - circle.center().y(),
source.z() - circle.center().z());
const double theta_source = get_theta<Kernel>(tmp, b1, b2, b3);
tmp = typename Kernel::Point_3(target.x() - circle.center().x(),
target.y() - circle.center().y(),
target.z() - circle.center().z());
const double theta_target = get_theta<Kernel>(tmp, b1, b2, b3);
subsample_arc_on_sphere_2<Kernel>(circle, theta_source, theta_target, b1, b2, out_pts, min_edge_size);
}
template <class Kernel, class ArcOnSphere, class Output_iterator>
void subsample_arc_on_sphere_2(const ArcOnSphere& arc,
Output_iterator out_pts,
double min_edge_size = 0.01)
{
subsample_arc_on_sphere_2<Kernel>(arc.supporting_circle(),
arc.supporting_circle().supporting_plane(),
arc.source(),
arc.target(),
out_pts,
min_edge_size);
}
} // namespace internal
} // namespace Triangulations_on_sphere_2
} // namespace CGAL
#endif // CGAL_EIGEN3_ENABLED
#endif // CGAL_TOS2_INTERNAL_ARC_ON_SPHERE_SUBSAMPLING_H

76
Triangulation_on_sphere_2/include/CGAL/Triangulation_on_sphere_2/internal/get_precision_bounds.h Normal file
View File

@ -0,0 +1,76 @@
// Copyright (c) 1997-2021 INRIA Sophia-Antipolis (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org).
//
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial
//
// Author(s) : Mael Rouxel-Labbé
// Sébastien Loriot
#ifndef CGAL_TRIANGULATION_ON_SPHERE_GET_PRECISION_BOUNDS_H
#define CGAL_TRIANGULATION_ON_SPHERE_GET_PRECISION_BOUNDS_H
#include <CGAL/license/Triangulation_on_sphere_2.h>
#include <CGAL/double.h>
#include <CGAL/Algebraic_structure_traits.h>
#include <CGAL/number_utils.h>
#include <type_traits>
namespace CGAL {
namespace Triangulations_on_sphere_2 {
namespace internal {
// @todo could do something more suble than requiring exact SQRT representation (Root_of_2 etc.)
template <typename FT,
bool has_exact_rep =
is_same_or_derived<
Field_with_sqrt_tag,
typename Algebraic_structure_traits<FT>::Algebraic_category>::value &&
!std::is_floating_point<FT>::value>
struct ToS2_precision_bound
{
ToS2_precision_bound(const FT /*radius*/) { }
// If you get a link error here, your number type is not supported and needs to implement its own specific values.
// See the manual of the package for more information: https://doc.cgal.org/latest/Triangulation_on_sphere_2/
FT get_squared_min_dist() const;
FT get_squared_min_radius() const;
FT get_squared_max_radius() const;
};
template <typename FT>
struct ToS2_precision_bound<FT, true> // exact representation of points on the sphere
{
ToS2_precision_bound(const FT radius) : _sq_radius(CGAL::square(radius)) { }
FT get_squared_min_dist() const { return 0; }
FT get_squared_min_radius() const { return _sq_radius; }
FT get_squared_max_radius() const { return _sq_radius; }
private:
FT _sq_radius;
};
template <>
struct ToS2_precision_bound<double, false>
{
ToS2_precision_bound(const double radius) : _radius(radius) { }
double get_squared_min_dist() const { return CGAL::square(_radius * std::pow(2, -23)); }
double get_squared_min_radius() const { return CGAL::square(_radius * (1 - std::pow(2, -50))); }
double get_squared_max_radius() const { return CGAL::square(_radius * (1 + std::pow(2, -50))); }
private:
double _radius;
};
} // namespace internal
} // namespace Triangulations_on_sphere_2
} // namespace CGAL
#endif // CGAL_TRIANGULATION_ON_SPHERE_GET_PRECISION_BOUNDS_H

69
Triangulation_on_sphere_2/include/CGAL/Triangulation_on_sphere_face_base_2.h Normal file
View File

@ -0,0 +1,69 @@
// Copyright (c) 1997, 2012, 2019 INRIA Sophia-Antipolis (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org).
//
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial
//
// Author(s) : Mariette Yvinec,
// Claudia Werner
// Mael Rouxel-Labbé
#ifndef CGAL_TRIANGULATION_ON_SPHERE_FACE_BASE_2_H
#define CGAL_TRIANGULATION_ON_SPHERE_FACE_BASE_2_H
#include <CGAL/license/Triangulation_on_sphere_2.h>
#include <CGAL/Triangulation_ds_face_base_2.h>
namespace CGAL {
template <typename Gt,
typename Fb = Triangulation_ds_face_base_2<> >
class Triangulation_on_sphere_face_base_2
: public Fb
{
public:
typedef Gt Geom_traits;
typedef typename Fb::Vertex_handle Vertex_handle;
typedef typename Fb::Face_handle Face_handle;
template < typename TDS2 >
struct Rebind_TDS {
typedef typename Fb::template Rebind_TDS<TDS2>::Other Fb2;
typedef Triangulation_on_sphere_face_base_2<Gt, Fb2> Other;
};
public:
void set_ghost(const bool b) { _ghost_flag = b; }
bool is_ghost() const { return _ghost_flag; }
public:
Triangulation_on_sphere_face_base_2()
: Fb(), _ghost_flag(false)
{ }
Triangulation_on_sphere_face_base_2(Vertex_handle v0,
Vertex_handle v1,
Vertex_handle v2)
: Fb(v0, v1, v2), _ghost_flag(false)
{ }
Triangulation_on_sphere_face_base_2(Vertex_handle v0,
Vertex_handle v1,
Vertex_handle v2,
Face_handle n0,
Face_handle n1,
Face_handle n2)
: Fb(v0, v1, v2, n0, n1, n2), _ghost_flag(false)
{ }
protected:
bool _ghost_flag;
};
} // namespace CGAL
#endif //CGAL_TRIANGULATION_ON_SPHERE_FACE_BASE_2_H

80
Triangulation_on_sphere_2/include/CGAL/Triangulation_on_sphere_vertex_base_2.h Normal file
View File

@ -0,0 +1,80 @@
// Copyright (c) 1997 INRIA Sophia-Antipolis (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org).
//
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial
//
//
// Author(s) : Mariette Yvinec
#ifndef CGAL_TRIANGULATION_ON_SPHERE_VERTEX_BASE_2_H
#define CGAL_TRIANGULATION_ON_SPHERE_VERTEX_BASE_2_H
#include <CGAL/license/Triangulation_on_sphere_2.h>
#include <CGAL/Triangulation_ds_vertex_base_2.h>
namespace CGAL {
template <typename GT,
typename Vb = Triangulation_ds_vertex_base_2<> >
class Triangulation_on_sphere_vertex_base_2
: public Vb
{
typedef typename Vb::Triangulation_data_structure Tds;
public:
typedef GT Geom_traits;
typedef typename GT::Point_on_sphere_2 Point;
typedef Tds Triangulation_data_structure;
typedef typename Tds::Face_handle Face_handle;
typedef typename Tds::Vertex_handle Vertex_handle;
template < typename TDS2 >
struct Rebind_TDS {
typedef typename Vb::template Rebind_TDS<TDS2>::Other Vb2;
typedef Triangulation_on_sphere_vertex_base_2<GT, Vb2> Other;
};
private:
Point _p;
public:
Triangulation_on_sphere_vertex_base_2 () : Vb(), _p() { }
Triangulation_on_sphere_vertex_base_2(const Point & p) : Vb(), _p(p) { }
Triangulation_on_sphere_vertex_base_2(const Point & p, Face_handle f) : Vb(f), _p(p) { }
Triangulation_on_sphere_vertex_base_2(Face_handle f) : Vb(f) { }
void set_point(const Point & p) { _p = p; }
const Point& point() const { return _p; }
// the non-const version of point() is undocument but needed to make the point iterator works
// using Lutz projection scheme
Point& point() { return _p; }
// the following trivial is_valid to allow the user of derived face base classes
// to add their own purpose checking
bool is_valid(bool /* verbose */ = false, int /* level */ = 0) const {return true;}
};
template < class GT, class Vb >
std::istream&
operator>>(std::istream &is, Triangulation_on_sphere_vertex_base_2<GT, Vb> &v)
{
// non combinatorial information. Default = point
return is >> static_cast<Vb&>(v) >> v.point();
}
template < class GT, class Vb >
std::ostream&
operator<<(std::ostream &os, const Triangulation_on_sphere_vertex_base_2<GT, Vb> &v)
{
// non combinatorial information. Default = point
return os << static_cast<const Vb&>(v) << v.point();
}
} // namespace CGAL
#endif //CGAL_TRIANGULATION_ON_SPHERE_VERTEX_BASE_2_H

365
Triangulation_on_sphere_2/include/CGAL/Triangulation_sphere_line_face_circulator_2.h Normal file
View File

@ -0,0 +1,365 @@
// Copyright (c) 1997, 2O12, 2019 INRIA Sophia-Antipolis (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org).
//
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial
//
// Author(s) : Claudia Werner, Mariette Yvinec
#ifndef CGAL_TRIANGULATION_SPHERE_LINE_FACE_CIRCULATOR_2_H
#define CGAL_TRIANGULATION_SPHERE_LINE_FACE_CIRCULATOR_2_H
#include <CGAL/license/Triangulation_on_sphere_2.h>
#include <CGAL/triangulation_assertions.h>
namespace CGAL {
template <class Triangulation_>
class Triangulation_sphere_line_face_circulator_2
: public Bidirectional_circulator_base<typename Triangulation_::Triangulation_data_structure::Face,
std::size_t>,
public Triangulation_cw_ccw_2
{
public:
typedef Triangulation_sphere_line_face_circulator_2<Triangulation_> Line_face_circulator;
typedef Triangulation_ Triangulation;
typedef typename Triangulation::Geom_traits Gt;
typedef typename Triangulation_::Triangulation_data_structure Tds;
typedef typename Tds::Vertex Vertex;
typedef typename Tds::Edge Edge;
typedef typename Tds::Face Face;
typedef typename Tds::Vertex_handle Vertex_handle;
typedef typename Tds::Face_handle Face_handle;
typedef typename Tds::Face_circulator Face_circulator;
typedef typename Gt::Point_2 Point;
typedef typename Triangulation::Locate_type Locate_type;
enum State {undefined = -1,
vertex_vertex,
vertex_edge,
edge_vertex,
edge_edge};
private:
Face_handle pos;
const Triangulation* _tr;
State s;
int i;
Point p, q;
public:
Triangulation_sphere_line_face_circulator_2()
: pos(), _tr(NULL), s(undefined), i(-1)
{ }
Triangulation_sphere_line_face_circulator_2(Vertex_handle v,
const Triangulation* tr,
const Point& dir);
Line_face_circulator& operator++();
Line_face_circulator& operator--();
Line_face_circulator operator++(int);
Line_face_circulator operator--(int);
Face* operator->() { return &*pos; }
Face& operator*() { return *pos; }
Face_handle handle() { return pos; }
operator const Face_handle() const { return pos; }
bool operator==(const Line_face_circulator& lfc) const;
bool operator!=(const Line_face_circulator& lfc) const;
bool operator==(const Face_handle& fh) const { return fh == pos; }
bool operator!=(const Face_handle& fh) const { return fh != pos; }
bool operator==(Nullptr_t CGAL_assertion_code(n)) const;
bool operator!=(Nullptr_t n) const;
bool is_empty() const;
bool locate(const Point& t, Locate_type &lt, int &li);
//private:
Triangulation_sphere_line_face_circulator_2(const Face_handle& face,
int index,
State state,
const Triangulation * t,
const Point& pp,
const Point& qq);
private:
void increment();
void decrement();
};
template <typename Triangulation>
inline
bool
operator==(typename Triangulation::Triangulation_data_structure::Face_handle fh,
Triangulation_sphere_line_face_circulator_2<Triangulation> fc)
{
return (fc == fh);
}
template <typename Triangulation>
inline
bool
operator!=(typename Triangulation::Triangulation_data_structure::Face_handle fh,
Triangulation_sphere_line_face_circulator_2<Triangulation> fc)
{
return (fc != fh);
}
template <typename Triangulation>
Triangulation_sphere_line_face_circulator_2<Triangulation>::
Triangulation_sphere_line_face_circulator_2(Vertex_handle v,
const Triangulation* tr,
const Point& dir)
: pos(), _tr(tr), s(undefined)
// begin at the face incident to v, traversed by the ray from v to dir
// or null iterator
{
CGAL_precondition((_tr->dimension() == 2) &&
(! _tr->xy_equal(v->point(), dir)));
p = v->point();
q = dir;
// find a vertex to the left of pq
// if there is no, the line_face_circulator is null
Face_circulator fc = _tr->incident_faces(v);
Face_circulator done(fc);
int ic = fc->index(v);
Vertex_handle vt = fc->vertex(cw(ic));
//Orientation o = _tr->orientation(p, q, vt->point());
while( _tr->orientation(p, q, vt->point()) != LEFT_TURN)
{
++fc;
ic = fc->index(v);
vt= fc->vertex(cw(ic));
if(fc == done)
{
*this = Line_face_circulator();
return;
}
}
Vertex_handle vr = fc-> vertex(ccw(ic));
Orientation pqr = RIGHT_TURN; // warning "pqr might be used uninitialized"
while((pqr = _tr->orientation(p, q, vr->point())) == LEFT_TURN)
{
--fc;
ic = fc->index(v);
vr = fc-> vertex(ccw(ic));
}
// first intersected face found
// [pqr] is COLLINEAR or RIGHT_TURN
ic = fc->index(v);
vt= fc->vertex(cw(ic));
CGAL_assertion (_tr->orientation(p, q, vt->point()) == LEFT_TURN);
if(pqr == COLLINEAR)
{
pos = fc;
s = vertex_vertex;
i = ccw(ic);
}
else // pqr == RIGHT_TURN
{
pos = fc;
s = vertex_edge;
i = ic ;
}
}
template <typename Triangulation>
inline
void
Triangulation_sphere_line_face_circulator_2<Triangulation>::
increment()
{
CGAL_precondition(pos != Face_handle());
if(s == vertex_vertex || s == edge_vertex)
{
Orientation o;
do
{
Face_handle n = pos->neighbor(cw(i));
i = n->index(pos);
pos = n;
o = _tr->orientation(p, q, pos->vertex(i)->point());
i = cw(i);
} while(o == LEFT_TURN);
if(o == COLLINEAR)
{
s = vertex_vertex;
i = ccw(i);
}
else
{
s = vertex_edge;
}
}
else
{
Face_handle n = pos->neighbor(i);
int ni = n->index(pos);
pos = n ;
Orientation o = _tr->orientation(p, q, pos->vertex(ni)->point());
switch(o)
{
case LEFT_TURN:
s = edge_edge;
i = ccw(ni);
break;
case RIGHT_TURN:
s = edge_edge;
i = cw(ni);
break;
default:
s = edge_vertex;
i = ni;
}
}
}
template <typename Triangulation>
void
Triangulation_sphere_line_face_circulator_2<Triangulation>::
decrement()
{
CGAL_precondition(pos != Face_handle());
if(s == vertex_vertex || s == vertex_edge)
{
if(s == vertex_vertex)
i = cw(i);
Orientation o;
do
{
Face_handle n = pos->neighbor(ccw(i));
i = n->index(pos);
pos = n;
o = _tr->orientation(p, q, pos->vertex(i)->point());
i = ccw(i);
}
while(o == LEFT_TURN);
s = (o == COLLINEAR) ? vertex_vertex : edge_vertex;
}
else // s == edge_edge || s == edge_vertex
{
// the following is not nice. A better solution is to say
// that index i is at the vertex that is alone on one side of l(p, q)
if(s == edge_edge)
i = (_tr->orientation(p, q, pos->vertex(i)->point()) == LEFT_TURN) ? cw(i) : ccw(i);
Face_handle n = pos->neighbor(i);
i = n->index(pos);
pos = n;
Orientation o = _tr->orientation(p, q, pos->vertex(i)->point());
s = (o == COLLINEAR) ? vertex_edge : edge_edge;
}
}
template <typename Triangulation>
inline
Triangulation_sphere_line_face_circulator_2<Triangulation>&
Triangulation_sphere_line_face_circulator_2<Triangulation>::
operator++()
{
CGAL_precondition(pos != Face_handle()) ;
increment();
return *this;
}
template <typename Triangulation>
inline
Triangulation_sphere_line_face_circulator_2<Triangulation>&
Triangulation_sphere_line_face_circulator_2<Triangulation>::
operator--()
{
CGAL_precondition(pos != Face_handle()) ;
decrement();
return *this;
}
template <typename Triangulation>
inline
Triangulation_sphere_line_face_circulator_2<Triangulation>
Triangulation_sphere_line_face_circulator_2<Triangulation>::
operator++(int)
{
Line_face_circulator tmp(*this);
++(*this);
return tmp;
}
template <typename Triangulation>
inline
Triangulation_sphere_line_face_circulator_2<Triangulation>
Triangulation_sphere_line_face_circulator_2<Triangulation>::
operator--(int)
{
Line_face_circulator tmp(*this);
--(*this);
return tmp;
}
template <typename Triangulation>
inline bool
Triangulation_sphere_line_face_circulator_2<Triangulation>::
operator==(const Line_face_circulator& lfc) const
{
CGAL_precondition(pos != Face_handle() && lfc.pos != Face_handle());
return (pos == lfc.pos && _tr == lfc._tr &&
s == lfc.s && p==lfc.p && q==lfc.q);
}
template <typename Triangulation>
inline bool
Triangulation_sphere_line_face_circulator_2<Triangulation>::
operator!=(const Line_face_circulator& lfc) const
{
return !(*this == lfc);
}
template <typename Triangulation>
inline bool
Triangulation_sphere_line_face_circulator_2<Triangulation>::
is_empty() const
{
return pos == Face_handle();
}
template <typename Triangulation>
inline bool
Triangulation_sphere_line_face_circulator_2<Triangulation>::
operator==(Nullptr_t CGAL_assertion_code(n)) const
{
CGAL_triangulation_assertion(n == NULL);
return pos == Face_handle();
}
template <typename Triangulation>
inline bool
Triangulation_sphere_line_face_circulator_2<Triangulation>::
operator!=(Nullptr_t n) const
{
CGAL_triangulation_assertion(n == NULL);
return !(*this == n);
}
} // end namespace CGAL
#endif // CGAL_TRIANGULATION_SPHERE_LINE_FACE_CIRCULATOR_2_H

1
Triangulation_on_sphere_2/package_info/Triangulation_on_sphere_2/copyright Normal file
View File

@ -0,0 +1 @@
INRIA (France)

27
Triangulation_on_sphere_2/package_info/Triangulation_on_sphere_2/dependencies Normal file
View File

@ -0,0 +1,27 @@
Algebraic_foundations
Algebraic_kernel_for_spheres
Arithmetic_kernel
BGL
Cartesian_kernel
Circular_kernel_3
Circulator
Distance_2
Distance_3
Filtered_kernel
Hash_map
Homogeneous_kernel
Installation
Intersections_2
Interval_support
Kernel_23
Modular_arithmetic
Number_types
Polygon
Profiling_tools
Property_map
STL_Extension
Solver_interface
Spatial_sorting
Stream_support
TDS_2
Triangulation_on_sphere_2

1
Triangulation_on_sphere_2/package_info/Triangulation_on_sphere_2/license.txt Normal file
View File

@ -0,0 +1 @@
GPL (v3 or later)

0
Triangulation_on_sphere_2/package_info/Triangulation_on_sphere_2/maintainer Normal file
View File

28
Triangulation_on_sphere_2/test/Triangulation_on_sphere_2/CMakeLists.txt Normal file
View File

@ -0,0 +1,28 @@
cmake_minimum_required(VERSION 3.1...3.15)
project( Triangulation_on_sphere_2_Tests )
find_package(CGAL REQUIRED COMPONENTS Core)
find_package(Eigen3 3.1.0) #(requires 3.1.0 or greater)
include(CGAL_Eigen3_support)
if ( CGAL_FOUND )
create_single_source_cgal_program( "test_dtos.cpp" )
create_single_source_cgal_program( "test_dtos2_remove.cpp" )
create_single_source_cgal_program( "test_dtos_degenerate_cases.cpp" )
create_single_source_cgal_program( "test_dtos_illegal_points.cpp" )
create_single_source_cgal_program( "test_dtos_projection_traits.cpp" )
create_single_source_cgal_program( "test_dtos_traits.cpp" )
if(TARGET CGAL::Eigen3_support)
create_single_source_cgal_program( "test_dtos_dual.cpp" )
target_link_libraries(test_dtos_dual PUBLIC CGAL::Eigen3_support)
endif()
else()
message(STATUS "This program requires the CGAL library, and will not be compiled.")
endif()

189
Triangulation_on_sphere_2/test/Triangulation_on_sphere_2/test_dtos.cpp Normal file
View File

@ -0,0 +1,189 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Delaunay_triangulation_on_sphere_traits_2.h>
#include <CGAL/Delaunay_triangulation_on_sphere_2.h>
#include <CGAL/enum.h>
#include <iostream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef K::Point_3 Point;
typedef CGAL::Delaunay_triangulation_on_sphere_traits_2<K> Gt;
typedef CGAL::Delaunay_triangulation_on_sphere_2<Gt> Tr;
int main(int, char**)
{
Gt traits(Point(1, 2, 3), 5);
std::vector<Point> points;
points.emplace_back( 2, 1, 1);
points.emplace_back(-2, 1, 1); // not on the sphere
points.emplace_back( 0, 1, 1);
points.emplace_back( 1, 2, 1);
points.emplace_back( 0, 1, 1); // duplicate
points.emplace_back( 1, 0, 1);
points.emplace_back( 1, 1, 2);
// test constructors
Tr tr0;
assert(tr0.geom_traits().center() == Point(0, 0, 0));
assert(tr0.geom_traits().radius() == 1);
Tr tr(traits);
assert(tr.geom_traits().center() == Point(1, 2, 3));
assert(tr.geom_traits().radius() == 5);
Tr trbis(tr);
assert(trbis.geom_traits().center() == tr.geom_traits().center());
assert(trbis.geom_traits().radius() == tr.geom_traits().radius());
assert(trbis.number_of_vertices() == tr.number_of_vertices());
assert(trbis.number_of_edges() == tr.number_of_edges());
assert(trbis.number_of_solid_faces() == tr.number_of_solid_faces());
assert(trbis.number_of_ghost_faces() == tr.number_of_ghost_faces());
Tr tr2(tr);
assert(tr2.geom_traits().center() == tr.geom_traits().center());
assert(tr2.geom_traits().radius() == tr.geom_traits().radius());
Tr tr3 = tr2;
assert(tr3.geom_traits().center() == tr2.geom_traits().center());
assert(tr3.geom_traits().radius() == tr2.geom_traits().radius());
Tr tr4(points.begin(), points.end(), Point(1, 1, 1), 1);
assert(tr4.number_of_vertices() == 5);
tr4.set_radius(0.123);
assert(tr4.number_of_vertices() == 0);
Tr tr5(points.begin(), points.end());
assert(tr5.geom_traits().center() == Point(0, 0, 0));
assert(tr5.geom_traits().radius() == 1);
assert(tr5.number_of_vertices() == 0);
Gt traits2(Point(1, 1, 1), 1);
Tr tr6(points.begin(), points.end(), traits2);
assert(tr6.number_of_vertices() == 5);
tr6.set_center(Point(1,2,3));
assert(tr6.number_of_vertices() == 0);
// center / radius setting
tr.insert(Point(6, 2, 3));
assert(tr.number_of_vertices() == 1);
tr.set_center_and_radius(Point(1, 1, 1), 1);
assert(tr3.number_of_vertices() == 0);
assert(tr.geom_traits().center() == Point(1, 1, 1));
assert(tr.geom_traits().radius() == 1);
tr.insert(points.begin(), points.end());
write_OFF("test_dtos.off", tr);
assert(tr.is_valid());
//////////////////////////////////////////////////////////////////////////////////////////////////
// test ranges
Tr::Vertex_handles vs = tr.vertex_handles();
assert(vs.size() == tr.number_of_vertices() && vs.size() == 5);
Tr::All_edges es = tr.all_edges();
Tr::Solid_edges ses = tr.solid_edges();
assert(es.size() == 9 && ses.size() == 8);
Tr::All_face_handles afs = tr.all_face_handles();
assert(afs.size() == tr.number_of_faces() && afs.size() == 6);
Tr::Solid_face_handles sfs = tr.solid_faces();
assert(sfs.size() == 4);
assert(tr.number_of_solid_faces() + tr.number_of_ghost_faces() == tr.number_of_faces());
Tr::Points pts = tr.points();
assert(pts.size() == tr.number_of_vertices() && pts.size() == 5);
//////////////////////////////////////////////////////////////////////////////////////////////////
// test iterators
// All vertices iterator
Tr::Vertices_iterator vit = tr.vertices_begin(), vend = tr.vertices_end();
assert(std::distance(vit, vend) == 5);
const Tr::Point& p0 = tr.point(vit);
++vit;
const Tr::Point& p1 = tr.point(vit++);
assert(p0 != p1);
std::advance(vit, -2);
assert(p0 == vit->point());
// All edges iterator
Tr::All_edges_iterator eit = tr.all_edges_begin(), eend = tr.all_edges_end();
assert(std::distance(eit, eend) == 9);
Tr::Edge e0 = *eit;
++eit;
Tr::Edge e1 = *eit++;
assert(e0 != e1);
std::advance(eit, -2);
assert(e0 == *eit);
// Solid edges iterator
Tr::Solid_edges_iterator seit = tr.solid_edges_begin(), seend = tr.solid_edges_end();
assert(std::distance(seit, seend) == 8);
Tr::Edge se0 = *seit;
++seit;
Tr::Edge se1 = *seit++;
assert(se0 != se1);
std::advance(seit, -2);
assert(se0 == *seit);
// Contour edges iterator
Tr::Contour_edges_iterator ceit=tr.contour_edges_begin(), ceend=tr.contour_edges_end();
assert(std::distance(ceit, ceend) == 4);
std::set<Tr::Edge> unique_edges;
while(ceit != ceend)
{
const Tr::Edge& e = *ceit;
unique_edges.insert(e);
assert(tr.is_ghost(e.first) != tr.is_ghost(tr.mirror_edge(e).first)); // xor
++ceit;
}
assert(unique_edges.size() == 4);
std::advance(ceit, -4);
assert(ceit == tr.contour_edges_begin());
// All faces iterator
Tr::All_faces_iterator fit = tr.all_faces_begin(), fend = tr.all_faces_end();
assert(std::distance(fit, fend) == 6);
Tr::Face_handle f0 = fit;
++fit;
Tr::Face_handle f1 = fit++;
assert(f0 != f1);
std::advance(fit, -2);
assert(f0 == fit);
// Solid faces iterator
Tr::Solid_faces_iterator sfit = tr.solid_faces_begin(), sfend = tr.solid_faces_end();
assert(std::distance(sfit, sfend) == 4);
Tr::Face_handle sf0 = sfit;
++sfit;
Tr::Face_handle sf1 = sfit++;
assert(sf0 != sf1);
std::advance(sfit, -2);
Tr::Face_handle sf2 = sfit;
assert(sf0 == sf2);
// Face Circulator
vit = tr.vertices_begin();
std::advance(vit, 3); // top of the pyramid
Tr::Face_circulator fc = tr.incident_faces(vit), end = fc;
assert(std::distance(++fc, end) + 1 == 4);
std::cout << "Done!" << std::endl;
return EXIT_SUCCESS;
}

88
Triangulation_on_sphere_2/test/Triangulation_on_sphere_2/test_dtos2_remove.cpp Normal file
View File

@ -0,0 +1,88 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Delaunay_triangulation_on_sphere_traits_2.h>
#include <CGAL/Projection_on_sphere_traits_3.h>
#include <CGAL/Delaunay_triangulation_on_sphere_2.h>
#include <cmath>
#include <fstream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef K::Point_3 Point;
typedef CGAL::Delaunay_triangulation_on_sphere_traits_2<K> Gt;
typedef CGAL::Projection_on_sphere_traits_3<K> Gt2;
typedef CGAL::Delaunay_triangulation_on_sphere_2<Gt> DTOS;
typedef CGAL::Delaunay_triangulation_on_sphere_2<Gt2> DTOS2;
void test1()
{
DTOS dtos;
dtos.set_radius(10);
Point a( 0, 0, 10);
Point b(-10/sqrt(3), -10/sqrt(3), 10/sqrt(3));
Point c( 10/sqrt(3), 10/sqrt(3), 10/sqrt(3));
Point d( 10/sqrt(3), -10/sqrt(3), 10/sqrt(3));
Point e(-10/sqrt(3), 10/sqrt(3), 10/sqrt(3));
DTOS::Vertex_handle v1 = dtos.insert(a);
DTOS::Vertex_handle v2 = dtos.insert(b);
DTOS::Vertex_handle v3 = dtos.insert(c);
DTOS::Vertex_handle v4 = dtos.insert(d);
DTOS::Vertex_handle v5 = dtos.insert(e);
assert(dtos.number_of_ghost_faces() == 2);
assert(dtos.dimension() == 2);
dtos.remove(v1);
assert(dtos.dimension() == 1);
dtos.remove(v2);
assert(dtos.dimension() == 1);
dtos.remove(v3);
assert(dtos.dimension() == 0);
dtos.remove(v4);
assert(dtos.dimension() == -1);
dtos.remove(v5);
assert(dtos.dimension() == -2);
}
void test2()
{
DTOS2 dtos;
dtos.set_radius(10);
Point a( 0, 0, 10);
Point b(-10/sqrt(3), -10/sqrt(3), 10/sqrt(3));
Point c( 10/sqrt(3), 10/sqrt(3), 10/sqrt(3));
Point d( 10/sqrt(3), -10/sqrt(3), 10/sqrt(3));
Point e(-10/sqrt(3), 10/sqrt(3), 10/sqrt(3));
Point f( 0, 0, -10);
/*DTOS2::Vertex_handle v1 = */dtos.insert(a);
/*DTOS2::Vertex_handle v2 = */dtos.insert(b);
/*DTOS2::Vertex_handle v3 = */dtos.insert(c);
/*DTOS2::Vertex_handle v4 = */dtos.insert(d);
/*DTOS2::Vertex_handle v5 = */dtos.insert(e);
DTOS2::Vertex_handle v6 = dtos.insert(f);
assert(dtos.number_of_ghost_faces() == 0);
dtos.remove(v6);
assert(dtos.number_of_ghost_faces() == 2);
assert(dtos.dimension() == 2);
}
int main(int, char**)
{
test1();
test2();
std::cout << "Done" << std::endl;
return EXIT_SUCCESS;
}

225
Triangulation_on_sphere_2/test/Triangulation_on_sphere_2/test_dtos_degenerate_cases.cpp Normal file
View File

@ -0,0 +1,225 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Exact_predicates_exact_constructions_kernel_with_sqrt.h>
#include <CGAL/Delaunay_triangulation_on_sphere_traits_2.h>
#include <CGAL/Projection_on_sphere_traits_3.h>
#include <CGAL/Triangulation_on_sphere_2.h>
#include <CGAL/Delaunay_triangulation_on_sphere_2.h>
#include <algorithm>
#include <iostream>
#include <vector>
template <class Vertex_handle, class Face_handle>
bool has_face(const Face_handle fh,
const Vertex_handle v0,
const Vertex_handle v1,
const Vertex_handle v2)
{
bool test1, test2, test3;
for(int i=0; i<=2; ++i)
{
test1 = (v0->point() == fh->vertex(i)->point());
if(test1)
break;
}
if(!test1)
return false;
for(int i=0; i<=2; ++i)
{
test2 = (v1->point() == fh->vertex(i)->point());
if(test2)
break;
}
if(!test2)
return false;
for(int i=0; i<=2; ++i)
{
test3 = (v2->point() == fh->vertex(i)->point());
if(test3)
break;
}
if(!test3)
return false;
return true;
}
template <class Triangul>
bool are_equal(const Triangul& triA, const Triangul& triB)
{
typedef typename Triangul::Vertex_handle Vertex_handle;
typedef typename Triangul::Face_handle Face_handle;
typedef typename Triangul::All_faces_iterator Face_iterator;
std::cout << "First triangulation:" << std::endl;
std::cout << "dimension: " << triA.dimension() << std::endl;
std::cout << triA.number_of_vertices() << " nv" << std::endl;
std::cout << triA.number_of_faces() << " nf" << std::endl;
std::cout << triA.number_of_ghost_faces() << " gf" << std::endl;
std::cout << "Second triangulation:" << std::endl;
std::cout << "dimension: " << triB.dimension() << std::endl;
std::cout << triB.number_of_vertices() << " nv" << std::endl;
std::cout << triB.number_of_faces() << " nf" << std::endl;
std::cout << triB.number_of_ghost_faces() << " gf" << std::endl;
if (triA.number_of_vertices()!= triB.number_of_vertices())
return false;
if (triA.number_of_faces()!= triB.number_of_faces())
return false;
if(triA.number_of_ghost_faces()!=triB.number_of_ghost_faces())
return false;
bool found = true;
for(Face_iterator fiA = triA.all_faces_begin(); fiA != triA.all_faces_end(); ++fiA)
{
found = false;
for(Face_iterator fiB=triB.all_faces_begin(); fiB!=triB.all_faces_end(); ++fiB)
{
Face_handle fb = Face_handle(fiB);
Face_handle fa = Face_handle(fiA);
Vertex_handle v0 = fa->vertex(0);
Vertex_handle v1 = fa->vertex(1);
Vertex_handle v2 = fa->vertex(2);
if(has_face(fb, v0, v1, v2))
{
found = true;
break;
}
}
assert(found);
}
return found;
}
// tests whether it is possible to insert points in degenerated positions
// and whether the result is uniquely defined after this.
template <typename DTOS, typename FT, typename PointContainer>
void test(const FT radius,
PointContainer coplanar_points)
{
std::cout << coplanar_points.size() << " input points" << std::endl;
DTOS dtos(CGAL::ORIGIN, radius);
for(const auto& p : coplanar_points) // to avoid Hilbert sort
dtos.insert(p);
assert(dtos.is_valid());
std::cout << "Triangulation:" << std::endl;
std::cout << "dimension: " << dtos.dimension() << std::endl;
std::cout << dtos.number_of_vertices() << " nv" << std::endl;
std::cout << dtos.number_of_edges() << " ne" << std::endl;
std::cout << dtos.number_of_faces() << " nf" << std::endl;
std::cout << dtos.number_of_ghost_faces() << " gf" << std::endl;
for(int i=0; i<10; ++i)
{
CGAL::cpp98::random_shuffle(coplanar_points.begin(), coplanar_points.end());
DTOS dtos2(CGAL::ORIGIN, radius);
for(const auto& p : coplanar_points) // to avoid Hilbert sort
dtos2.insert(p);
assert(dtos2.is_valid());
assert(are_equal(dtos, dtos2));
std::cout << dtos2.number_of_vertices() << " nv and " << dtos2.number_of_faces() << " nf" << std::endl;
}
}
template <typename K>
void test_kernel()
{
std::cout << "Test: " << typeid(K).name() << std::endl;
typedef typename K::FT FT;
typedef typename K::Point_3 Point_3;
typedef CGAL::Delaunay_triangulation_on_sphere_traits_2<K> Gt;
typedef CGAL::Delaunay_triangulation_on_sphere_2<Gt> DTOS;
typedef CGAL::Projection_on_sphere_traits_3<K> PGt;
typedef CGAL::Delaunay_triangulation_on_sphere_2<PGt> PDTOS;
const FT radius = 100;
const FT radius2 = CGAL::square(radius);
const FT denom = FT(1) / CGAL::sqrt(FT(2));
const FT z = CGAL::sqrt(radius2 - 1);
std::vector<Point_3> coplanar_low_dim { Point_3(0,0,radius), Point_3(radius,0,0), Point_3(0,radius,0)
// , Point_3(0,0,-radius)
};
// Points are coplanar and coplanar with the center of the sphere
std::vector<Point_3> coplanar_points { Point_3( radius*denom, radius*denom, 0),
Point_3(- radius*denom, radius*denom, 0),
Point_3(- radius*denom, - radius*denom, 0),
Point_3( radius*denom, - radius*denom, 0),
Point_3( radius, 0, 0),
Point_3( 0, 0, radius) };
std::vector<Point_3> coplanar_points_on_great_circle { Point_3( 0, 0, radius),
Point_3( denom, denom, z),
Point_3(-denom, -denom, z),
Point_3( denom, denom, -z),
Point_3(-denom, -denom, -z),
Point_3( 0, 0,-radius) };
std::vector<Point_3> coplanar_points_on_circle { Point_3( denom, denom, z),
Point_3(-denom, -denom, z),
Point_3( 0, 1, z)/*,
Point_3( 1, 0, z),
Point_3(-denom, denom, z),
Point_3( denom, -denom, z)*/ };
// -----------------------------------------------------------------------------------------------
std::cout << "Testing with Delaunay_triangulation_sphere_traits:" << std::endl;
test<DTOS>(radius, coplanar_low_dim);
test<DTOS>(radius, coplanar_points);
test<DTOS>(radius, coplanar_points_on_great_circle);
test<DTOS>(radius, coplanar_points_on_circle);
return;
// -----------------------------------------------------------------------------------------------
std::cout << "Testing with projection: " << std::endl;
PGt traits(CGAL::ORIGIN, radius);
typename PGt::Construct_point_on_sphere_2 to_s2 = traits.construct_point_on_sphere_2_object();
std::vector<typename PGt::Point_on_sphere_2> coplanar_ppoints;
std::vector<typename PGt::Point_on_sphere_2> coplanar_ppoints_on_great_circle;
std::transform(coplanar_points.begin(), coplanar_points.end(), std::back_inserter(coplanar_ppoints), to_s2);
test<PDTOS>(radius, coplanar_ppoints);
// tests the convenience API (passing P3)
test<PDTOS>(radius, coplanar_low_dim);
test<PDTOS>(radius, coplanar_points_on_circle);
test<PDTOS>(radius, coplanar_points_on_great_circle);
test<PDTOS>(radius, coplanar_points_on_circle);
}
int main(int, char**)
{
typedef CGAL::Exact_predicates_exact_constructions_kernel_with_sqrt EPECK_w_SQRT;
typedef CGAL::Exact_predicates_inexact_constructions_kernel EPICK;
test_kernel<EPICK>();
test_kernel<EPECK_w_SQRT>();
std::cout << "Done" << std::endl;
return EXIT_SUCCESS;
}

164
Triangulation_on_sphere_2/test/Triangulation_on_sphere_2/test_dtos_dual.cpp Normal file
View File

@ -0,0 +1,164 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Geographical_coordinates_traits_2.h>
#include <CGAL/Delaunay_triangulation_on_sphere_traits_2.h>
#include <CGAL/Projection_on_sphere_traits_3.h>
#include <CGAL/Delaunay_triangulation_on_sphere_2.h>
#include <CGAL/enum.h>
#include <CGAL/point_generators_3.h>
#include <iostream>
#include <fstream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel EPICK;
template <typename K, typename PointRange>
void test(const PointRange& points)
{
// typedef CGAL::Geographical_coordinates_traits_2<K> Gt;
typedef CGAL::Projection_on_sphere_traits_3<K> Gt;
typedef CGAL::Delaunay_triangulation_on_sphere_2<Gt> Dt;
typedef typename Gt::SK SK;
typedef typename Dt::FT FT;
typedef typename Dt::Point_3 Point_3;
typedef typename Dt::Point Point_on_sphere_2;
typedef typename Dt::Segment_3 Segment_3;
typedef typename Dt::Arc_on_sphere_2 Arc_on_sphere_2;
typedef typename Dt::Face_handle Face_handle;
Dt dt;
dt.insert(points.begin(), points.end());
// The triangulation, using straight edges
write_OFF("result.off", dt, CGAL::parameters::stream_precision(17));
std::ofstream out_primal("edges_primal.polylines.cgal");
out_primal.precision(17);
std::ofstream out_dual("edges_dual.polylines.cgal");
out_dual.precision(17);
for(typename Dt::Solid_faces_iterator fit = dt.solid_faces_begin(); fit!=dt.solid_faces_end(); ++fit)
{
Face_handle fh = fit;
const Point_3 c = dt.dual(fh);
CGAL_USE(c);
const Point_on_sphere_2 cs = dt.dual_on_sphere(fh);
const FT r = dt.geom_traits().radius();
assert(CGAL::abs(CGAL::squared_distance(dt.construct_point(cs),
dt.geom_traits().center()) - r*r) <= 1e-10);
for(int i=0; i<3; ++i)
{
Face_handle nfh = fh->neighbor(i);
if(!dt.is_ghost(nfh) && fh < nfh)
continue;
typename Dt::Edge e(fh, i);
const bool diametral_edge = CGAL::collinear(dt.construct_point(dt.point(fh, (i+1)%3)),
dt.construct_point(dt.point(fh, (i+2)%3)),
dt.geom_traits().center());
// primal
if(!diametral_edge)
{
Arc_on_sphere_2 as = dt.segment_on_sphere(e);
std::vector<typename SK::Point_3> discretization_points;
CGAL::Triangulations_on_sphere_2::internal::subsample_arc_on_sphere_2<SK>(as, std::back_inserter(discretization_points));
assert(discretization_points.size() >= 2);
for(std::size_t i=0; i<discretization_points.size()-1; ++i)
out_primal << "2 " << discretization_points[i] << " " << discretization_points[i+1] << "\n";
}
else
{
Segment_3 s = dt.segment(e);
out_primal << "2 " << s.source() << " " << s.target() << "\n";
}
// Dual
if(dt.is_contour(e))
continue;
const Point_on_sphere_2 c1 = dt.circumcenter_on_sphere(e.first);
const Point_on_sphere_2 c2 = dt.circumcenter_on_sphere(dt.mirror_edge(e).first);
// That should never be possible, but with constructions...
const bool diametral_dual = CGAL::collinear(dt.construct_point(c1),
dt.construct_point(c2),
dt.geom_traits().center());
if(!diametral_dual)
{
Arc_on_sphere_2 ad = dt.dual_on_sphere(e);
std::vector<typename SK::Point_3> discretization_points;
CGAL::Triangulations_on_sphere_2::internal::subsample_arc_on_sphere_2<SK>(ad, std::back_inserter(discretization_points));
assert(discretization_points.size() >= 2);
for(std::size_t i=0; i<discretization_points.size()-1; ++i)
out_dual << "2 " << discretization_points[i] << " " << discretization_points[i+1] << "\n";
}
else
{
Segment_3 d = dt.dual(e);
out_dual << "2 " << d.source() << " " << d.target() << "\n";
};
}
}
}
template <typename K>
void test_with_ghost_faces()
{
std::vector<typename K::Point_3> points;
points.emplace_back( 1, 0, 0);
points.emplace_back( 0, 1, 0);
points.emplace_back( 1, 1, 1);
points.emplace_back( -1, -1, -1);
points.emplace_back( 0.5, 0.5, 0);
points.emplace_back( 1, 1, 0);
points.emplace_back(0.25, 0.25, 0);
// points.emplace_back( 0, 0, 1);
// points.emplace_back( 2, 0, 0);
// points.emplace_back( 0.5, 0, 0);
// points.emplace_back( 0.5, 0, 0.5);
// points.emplace_back( 0.5, 0, -0.5);
// points.emplace_back( 0, 0, 0.5);
// points.emplace_back( -1, 0, 0);
return test<K>(points);
}
template <typename K>
void test_random_data(const std::size_t num_of_pts)
{
typedef typename K::Point_3 Point_3;
typedef CGAL::Creator_uniform_3<double, Point_3> Creator;
CGAL::Random r;
std::cout << "Seed is " << r.get_seed() << std::endl;
CGAL::Random_points_on_sphere_3<Point_3, Creator> on_sphere(1 /*radius*/, r);
std::vector<Point_3> points;
points.reserve(num_of_pts);
for(std::size_t c=0; c<num_of_pts; ++c)
points.push_back(*on_sphere++);
return test<K>(points);
}
int main(int argc, char** argv)
{
const std::size_t num_pts = (argc > 1) ? std::atoi(argv[1]) : 100;
test_with_ghost_faces<EPICK>();
test_random_data<EPICK>(num_pts);
std::cout << "Done!" << std::endl;
return EXIT_SUCCESS;
}

141
Triangulation_on_sphere_2/test/Triangulation_on_sphere_2/test_dtos_illegal_points.cpp Normal file
View File

@ -0,0 +1,141 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Delaunay_triangulation_on_sphere_traits_2.h>
#include <CGAL/Projection_on_sphere_traits_3.h>
#include <CGAL/Delaunay_triangulation_on_sphere_2.h>
#include <iostream>
#include <vector>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Delaunay_triangulation_on_sphere_traits_2<K> Gt;
typedef CGAL::Projection_on_sphere_traits_3<K> Pgt;
typedef CGAL::Delaunay_triangulation_on_sphere_2<Gt> DTOS;
typedef CGAL::Delaunay_triangulation_on_sphere_2<Pgt> PDTOS;
typedef K::Point_3 Point;
void test_projection_traits()
{
const Point center = Point(0,0,0);
const double radius = 1000;
Pgt traits(center, radius);
PDTOS pdtos(traits);
std::vector<Point> points;
// legal points
Point p1( radius/sqrt(2), 0, radius/sqrt(2));
Point p2(-radius/sqrt(3), radius/sqrt(3), -radius/sqrt(3));
Point p3( radius, 0, 0);
Point p4(-radius/sqrt(2), -radius/sqrt(2), 0);
points.push_back(p1);
points.push_back(p2);
points.push_back(p3);
points.push_back(p4);
// points original not on sphere
Point p5(radius, 0, -0.5*radius);
points.push_back(p5);
// Points too close
Point p6( 5*radius, std::pow(2, -25), 0);
Point p7(-7*radius/sqrt(3), 7*radius/sqrt(3), -7*radius/sqrt(3) + std::pow(2,-35));
points.push_back(p6);
points.push_back(p7);
pdtos.insert(points.begin(), points.end());
assert(pdtos.number_of_vertices() == 5);
assert(pdtos.is_valid());
}
void test_Delaunay_traits()
{
typedef DTOS::Point Point;
const double radius = 1000;
DTOS dtos;
dtos.set_radius(radius);
std::vector<Point> points;
// legal points
Point p1( radius/sqrt(2), 0, radius/sqrt(2));
Point p2(-radius/sqrt(3), radius/sqrt(3), -radius/sqrt(3));
Point p3( radius, 0, 0);
Point p4(-radius/sqrt(2), -radius/sqrt(2), 0);
points.push_back(p1);
points.push_back(p2);
points.push_back(p3);
points.push_back(p4);
// points not on sphere
Point p5(-1.5*radius, 0, 0);
Point p6( radius, 0, -0.5*radius);
points.push_back(p5);
points.push_back(p6);
// points to close
Point p7(radius/sqrt(2), radius*std::pow(2, -15), radius/sqrt(2));
Point p8( radius, 0, radius*std::pow(2, -14));
points.push_back(p7);
points.push_back(p8);
dtos.insert(points.begin(), points.end());
assert(dtos.is_valid());
assert(dtos.number_of_vertices() == 4);
dtos.clear();
dtos.set_radius(radius);
dtos.insert(p1);
assert(dtos.number_of_vertices() == 1);
dtos.insert(p1);
assert(dtos.number_of_vertices() == 1);
dtos.insert(p2);
assert(dtos.number_of_vertices() == 2);
dtos.insert(p2);
assert(dtos.number_of_vertices() == 2);
dtos.insert(p3);
assert(dtos.number_of_vertices() == 3);
dtos.insert(p3);
assert(dtos.number_of_vertices() == 3);
dtos.insert(p4);
assert(dtos.number_of_vertices() == 4);
dtos.insert(p4);
assert(dtos.number_of_vertices() == 4);
dtos.insert(p5);
assert(dtos.number_of_vertices() == 4);
dtos.insert(p6);
assert(dtos.number_of_vertices() == 4);
dtos.insert(p7);
assert(dtos.number_of_vertices() == 4);
dtos.insert(p8);
assert(dtos.number_of_vertices() == 4);
dtos.insert(p8);
assert(dtos.number_of_vertices() == 4);
}
int main()
{
std::cout << "Test inserting illegal points with Delaunay_traits" << std::endl;
test_Delaunay_traits();
std::cout << "Test inserting illegal points with Projection_traits" << std::endl;
test_projection_traits();
std::cout << "Done" << std::endl;
return EXIT_SUCCESS;
}

121
Triangulation_on_sphere_2/test/Triangulation_on_sphere_2/test_dtos_projection_traits.cpp Normal file
View File

@ -0,0 +1,121 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Projection_on_sphere_traits_3.h>
#include <CGAL/enum.h>
#include <iostream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Projection_on_sphere_traits_3<K> Gt;
typedef Gt::Point_3 Point_3;
typedef Gt::Point_on_sphere_2 Point_on_sphere_2;
typedef Gt::Construct_point_on_sphere_2 Construct_point_on_sphere_2;
typedef Gt::Compare_on_sphere_2 Compare_on_sphere_2;
typedef Gt::Equal_on_sphere_2 Equal_on_sphere_2;
typedef Gt::Collinear_are_strictly_ordered_on_great_circle_2 Collinear_are_strictly_ordered_on_great_circle_2;
typedef Gt::Orientation_on_sphere_2 Orientation_on_sphere_2;
typedef Gt::Side_of_oriented_circle_on_sphere_2 Side_of_oriented_circle_on_sphere_2;
typedef Gt::Compare_on_sphere_2 Compare_on_sphere_2;
int main(int, char**)
{
Point_3 c(0, 0, 0);
Gt traits(c, 10);
Equal_on_sphere_2 equal_on_sphere = traits.equal_on_sphere_2_object();
Compare_on_sphere_2 compare_on_sphere = traits.compare_on_sphere_2_object();
Collinear_are_strictly_ordered_on_great_circle_2 inside_cone = traits.collinear_are_strictly_ordered_on_great_circle_2_object();
Orientation_on_sphere_2 orientation = traits.orientation_on_sphere_2_object();
Side_of_oriented_circle_on_sphere_2 side_of_oriented_circle = traits.side_of_oriented_circle_on_sphere_2_object();
Construct_point_on_sphere_2 cst = traits.construct_point_on_sphere_2_object();
// Testing with Points projected on unit sphere
Point_3 p11(1.5, 1.5, 1.5);
Point_3 p12( -1, 1, 1);
Point_3 p13( 0, 0.7, 0);
Point_3 p14( 0, -2, 0);
Point_3 p21( 2, 0, -1);
Point_3 p22(-1.8, 0, -0.7);
Point_3 p23( 0, 0, 2.6);
// Points with same coordinates
Point_3 p31(0.6, 0.3, -1);
Point_3 p32(0.6, 0.3, -1);
Point_3 p33(0.6, 0.3, 1);
// Equal points
Point_3 p41( 1, 1, 1);
Point_3 p42(0.7, 0.7, 0.7);
Point_3 p43(0.3, 0.6, -2.4);
Point_3 p44(0.6, 1.2, -4.8);
// Inside Cone
Point_3 p51( 1, 1, 1);
Point_3 p52( -1, 1, 1);
Point_3 p53( 1, 3, 3); // inside
Point_3 p54(0.8, 0.8, 0.8); // boundary
Point_3 p55( 0, -2, -2); // outside
//distance
Point_3 p61( 1, 1, 0);
Point_3 p62( -1, 1, 0);
// transform to projected points
Point_on_sphere_2 pp11 = cst(p11);
Point_on_sphere_2 pp12 = cst(p12);
Point_on_sphere_2 pp13 = cst(p13);
Point_on_sphere_2 pp14 = cst(p14);
Point_on_sphere_2 pp21 = cst(p21);
Point_on_sphere_2 pp22 = cst(p22);
Point_on_sphere_2 pp23 = cst(p23);
Point_on_sphere_2 pp31 = cst(p31);
Point_on_sphere_2 pp32 = cst(p32);
Point_on_sphere_2 pp33 = cst(p33);
Point_on_sphere_2 pp41 = cst(p41);
Point_on_sphere_2 pp42 = cst(p42);
Point_on_sphere_2 pp43 = cst(p43);
Point_on_sphere_2 pp44 = cst(p44);
Point_on_sphere_2 pp51 = cst(p51);
Point_on_sphere_2 pp52 = cst(p52);
Point_on_sphere_2 pp53 = cst(p53);
Point_on_sphere_2 pp54 = cst(p54);
Point_on_sphere_2 pp55 = cst(p55);
std::cout << "Test Orientation" << std::endl;
assert(orientation(pp11, pp12, pp13) == CGAL::NEGATIVE);
assert(orientation(pp13, pp12, pp11) == CGAL::POSITIVE);
assert(orientation(pp21, pp22, pp23) == CGAL::ON_ORIENTED_BOUNDARY);
std::cout << "Test Side_of_oriented_circle_on_sphere_2" << std::endl;
assert(side_of_oriented_circle(pp11, pp12, pp13, pp14) == CGAL::POSITIVE);
assert(side_of_oriented_circle(pp14, pp11, pp12, pp13) == CGAL::NEGATIVE);
assert(side_of_oriented_circle(pp21, pp22, pp23, pp11) == CGAL::POSITIVE);
std::cout << "Test Equal_on_sphere_2" << std::endl;
assert(equal_on_sphere(pp41, pp41));
assert(equal_on_sphere(pp41, pp42));
assert(equal_on_sphere(pp43, pp44));
std::cout << "Test Inside_cone_sphere_2" << std::endl;
assert(inside_cone(pp51, pp52, pp53));
assert(!inside_cone(pp51, pp52, pp54));
assert(!inside_cone(pp51, pp52, pp55));
std::cout << "Test Compare_on_sphere_2" << std::endl;
assert(compare_on_sphere(pp31, pp31) == CGAL::EQUAL);
assert(compare_on_sphere(pp31, pp32) == CGAL::EQUAL);
assert(compare_on_sphere(pp31, pp33) == CGAL::SMALLER);
assert(compare_on_sphere(pp33, pp31) == CGAL::LARGER);
std::cout << "Done" << std::endl;
return EXIT_SUCCESS;
}

162
Triangulation_on_sphere_2/test/Triangulation_on_sphere_2/test_dtos_traits.cpp Normal file
View File

@ -0,0 +1,162 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Delaunay_triangulation_on_sphere_traits_2.h>
#include <CGAL/enum.h>
#include <CGAL/Origin.h>
#include <iostream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef K::Point_3 Point;
typedef CGAL::Delaunay_triangulation_on_sphere_traits_2<K> Gt;
typedef Gt::Orientation_on_sphere_2 Orientation_on_sphere_2;
typedef Gt::Side_of_oriented_circle_on_sphere_2 Side_of_oriented_circle_on_sphere_2;
// This only tests non-trivial functions (a lot of functors are just calls to some functors
// of the 3D linear kernel, such as:
// typedef typename LK::Compare_xyz_3 Compare_on_sphere_2;
int main(int, char**)
{
Point P0 ( 0, 0, 0);
Point P1 ( 1, 0, 0);
Point P2 ( 0, 1, 0);
Point P3 ( 0, 0, 1);
Point P4 ( 1, 1, 1);
Point P5 ( -1, -1, -1);
Point P6 ( 0.5, 0.5, 0);
Point P7 ( 2, 0, 0);
Point P8 ( 0.5, 0, 0);
Point P9 ( 1, 1, 0);
Point P10(0.25, 0.25, 0);
Point P11( 0.5, 0, 0.5);
Point P12( 0.5, 0, -0.5);
Point P13( 0, 0, 0.5);
Point P14( -1, 0, 0);
Gt traits;
Gt traits_2(CGAL::ORIGIN);
Gt traits_3(Point(1, 2, 3), 15);
Gt traits_4(traits_3);
assert(traits_4.radius() == traits_3.radius());
assert(traits_4.center() == traits_3.center());
Gt traits_5 = traits_4;
assert(traits_4.radius() == traits_5.radius());
assert(traits_4.center() == traits_5.center());
std::cout << "Test Side_of_oriented_circle_on_sphere_2" << std::endl;
CGAL::Oriented_side result;
result = traits.side_of_oriented_circle_on_sphere_2_object()(P1, P2, P3, P4);
assert(result == CGAL::ON_POSITIVE_SIDE);
result = traits.side_of_oriented_circle_on_sphere_2_object()(P1, P2, P3, P5);
assert(result == CGAL::ON_NEGATIVE_SIDE);
result = traits.side_of_oriented_circle_on_sphere_2_object()(P1, P2, P3, P6);
assert(result == CGAL::ON_ORIENTED_BOUNDARY);
// Former tests for the degenerate cases,
// https://github.com/MaelRL/cgal/blob/Triangulation_on_sphere_2-GF-old/
// Triangulation_on_sphere_2/include/CGAL/Delaunay_triangulation_sphere_traits_2.h#L29
//
// power_test_2(p,q,r) where p,q and r are coplanar points with the center of the sphere
// result = traits.power_test_2_object()(P1, P2, P9);
// assert(result == CGAL::ON_POSITIVE_SIDE);
// result = traits.power_test_2_object()(P2, P1, P9);
// assert(result == CGAL::ON_POSITIVE_SIDE);
// result = traits.power_test_2_object()(P1, P2, P6);
// assert(result == CGAL::ON_ORIENTED_BOUNDARY);
// result = traits.power_test_2_object()(P1, P2, P10);
// assert(result == CGAL::ON_NEGATIVE_SIDE);
// result = traits.power_test_2_object()(P2, P1, P10);
// assert(result == CGAL::ON_NEGATIVE_SIDE);
// power_test_2(p,q) where p, q, and the center of the sphere are colinear
// result = traits.power_test_2_object()(P1, P7);
// assert(result == CGAL::ON_POSITIVE_SIDE);
// result = traits.power_test_2_object()(P1, P1);
// assert(result == CGAL::ON_ORIENTED_BOUNDARY);
// result = traits.power_test_2_object()(P1, P8);
// assert(result == CGAL::ON_NEGATIVE_SIDE);
// With a center different from the origin
Gt traits_6(Point(0,0,1));
// power_test_2(p,q,r,s)
result = traits_6.side_of_oriented_circle_on_sphere_2_object()(P0, P2, P1, P5);
assert(result == CGAL::ON_POSITIVE_SIDE);
result = traits_6.side_of_oriented_circle_on_sphere_2_object()(P0, P2, P1, P4);
assert(result == CGAL::ON_NEGATIVE_SIDE);
result = traits_6.side_of_oriented_circle_on_sphere_2_object()(P0, P1, P1, P10);
assert(result == CGAL::ON_ORIENTED_BOUNDARY);
// power_test_2(p,q,r) where p,q and r are coplanar points
// result = traits_6.power_test_2_object()(P14, P1, P12);
// assert(result == CGAL::ON_POSITIVE_SIDE);
// result = traits_6.power_test_2_object()(P1, P2, P6);
// assert(result == CGAL::ON_ORIENTED_BOUNDARY);
// result = traits_6.power_test_2_object()(P14, P1, P11);
// assert(result == CGAL::ON_NEGATIVE_SIDE);
// power_test_2(p,q) where p, q and sphere are colinear
// result = traits_6.power_test_2_object()(P13, P0);
// assert(result == CGAL::ON_POSITIVE_SIDE);
// result = traits_6.power_test_2_object()(P13, P13);
// assert(result == CGAL::ON_ORIENTED_BOUNDARY);
// result = traits_6.power_test_2_object()(P0, P13);
// assert(result == CGAL::ON_NEGATIVE_SIDE);
std::cout << "Test Orientation_on_sphere_2" << std::endl;
Point p21( 0.5, 0.5, sqrt(0.75));
Point p22(-0.5, 0.5, sqrt(0.75));
Point p23( 0, -1, sqrt(0.75));
Point p24( 0, 0, -1);
Point p25( 0, 0, -1.5);
Point p26( 1, 1, 2*sqrt(0.75));
result = traits.orientation_on_sphere_2_object()(p21, p22, p23);
assert(result == CGAL::ON_POSITIVE_SIDE);
result = traits.orientation_on_sphere_2_object()(p23, p22, p21);
assert(result == CGAL::ON_NEGATIVE_SIDE);
result = traits.side_of_oriented_circle_on_sphere_2_object()(p24, p21, p22, p23);
assert(result == CGAL::ON_POSITIVE_SIDE);
std::cout << "Test Equal_on_sphere_2" << std::endl;
bool coradial;
coradial = traits.equal_on_sphere_2_object()(p24, p25);
assert(coradial);
coradial = traits.equal_on_sphere_2_object()(p25, p24);
assert(coradial);
coradial = traits.equal_on_sphere_2_object()(p21, p26);
assert(coradial);
coradial = traits.equal_on_sphere_2_object()(p22, p26);
assert(!coradial);
std::cout << "Test Collinear_are_strictly_ordered_on_great_circle_2" << std::endl;
Point p31(1,1,1);
Point p32(-1,1,1);
Point p33(0, sqrt(1.5),sqrt(1.5));
Point p34(0, -sqrt(2),-1);
Point p35(0.9,0.9,0.9);
Point p36(0, sqrt(2.5),0);
bool inside;
// cone defined by 0, p31 and p32
inside = traits.collinear_are_strictly_ordered_on_great_circle_2_object()(p31, p32, p33);
assert(inside);
inside = traits.collinear_are_strictly_ordered_on_great_circle_2_object()(p31, p32, p34);
assert(!inside);
inside = traits.collinear_are_strictly_ordered_on_great_circle_2_object()(p31, p32, p35);
assert(!inside);
//not coplanar
inside = traits.collinear_are_strictly_ordered_on_great_circle_2_object()(p31, p32, p36);
assert(!inside);
std::cout << "Done" << std::endl;
return EXIT_SUCCESS;
}