\( \newcommand{\E}{\mathrm{E}} \) \( \newcommand{\A}{\mathrm{A}} \) \( \newcommand{\R}{\mathrm{R}} \) \( \newcommand{\N}{\mathrm{N}} \) \( \newcommand{\Q}{\mathrm{Q}} \) \( \newcommand{\Z}{\mathrm{Z}} \) \( \def\ccSum #1#2#3{ \sum_{#1}^{#2}{#3} } \def\ccProd #1#2#3{ \sum_{#1}^{#2}{#3} }\)
CGAL 4.12 - 2D and 3D Linear Geometry Kernel
User Manual

Authors
Hervé Brönnimann, Andreas Fabri, Geert-Jan Giezeman, Susan Hert, Michael Hoffmann, Lutz Kettner, Sylvain Pion, and Stefan Schirra

Introduction

CGAL, the Computational Geometry Algorithms Library, is written in C++ and consists of three major parts. The first part is the kernel, which consists of constant-size non-modifiable geometric primitive objects and operations on these objects. The objects are represented both as stand-alone classes that are parameterized by a representation class, which specifies the underlying number types used for calculations and as members of the kernel classes, which allows for more flexibility and adaptability of the kernel. The second part is a collection of basic geometric data structures and algorithms, which are parameterized by traits classes that define the interface between the data structure or algorithm and the primitives they use. In many cases, the kernel classes provided in CGAL can be used as traits classes for these data structures and algorithms. The third part of the library consists of non-geometric support facilities, such as circulators, random sources, I/O support for debugging and for interfacing CGAL to various visualization tools.

This part of the reference manual covers the kernel. The kernel contains objects of constant size, such as point, vector, direction, line, ray, segment, triangle, iso-oriented rectangle and tetrahedron. With each type comes a set of functions which can be applied to an object of this type. You will typically find access functions (e.g. to the coordinates of a point), tests of the position of a point relative to the object, a function returning the bounding box, the length, or the area of an object, and so on. The CGAL kernel further contains basic operations such as affine transformations, detection and computation of intersections, and distance computations.

Robustness

The correctness proof of nearly all geometric algorithms presented in theory papers assumes exact computation with real numbers. This leads to a fundamental problem with the implementation of geometric algorithms. Naively, often the exact real arithmetic is replaced by inexact floating-point arithmetic in the implementation. This often leads to acceptable results for many input data. However, even for the implementation of the simplest geometric algorithms this simplification occasionally does not work. Rounding errors introduced by an inaccurate arithmetic may lead to inconsistent decisions, causing unexpected failures for some correct input data. There are many approaches to this problem, one of them is to compute exactly (compute so accurate that all decisions made by the algorithm are exact) which is possible in many cases but more expensive than standard floating-point arithmetic. C. M. Hoffmann [3], [2] illustrates some of the problems arising in the implementation of geometric algorithms and discusses some approaches to solve them. A more recent overview is given in [5]. The exact computation paradigm is discussed by Yap and Dubé [6] and Yap [7].

In CGAL you can choose the underlying number types and arithmetic. You can use different types of arithmetic simultaneously and the choice can be easily changed, e.g. for testing. So you can choose between implementations with fast but occasionally inexact arithmetic and implementations guaranteeing exact computation and exact results. Of course you have to pay for the exactness in terms of execution time and storage space. See the dedicated chapter for more details on number types and their capabilities and performance.

Kernel Representations

Our object of study is the \( d\)-dimensional affine Euclidean space. Here we are mainly concerned with cases \( d=2\) and \( d=3\). Objects in that space are sets of points. A common way to represent the points is the use of Cartesian coordinates, which assumes a reference frame (an origin and \( d\) orthogonal axes). In that framework, a point is represented by a \( d\)-tuple \( (c_0,c_1,\ldots,c_{d-1})\), and so are vectors in the underlying linear space. Each point is represented uniquely by such Cartesian coordinates. Another way to represent points is by homogeneous coordinates. In that framework, a point is represented by a \( (d+1)\)-tuple \( (h_0,h_1,\ldots,h_d)\). Via the formulae \( c_i = h_i/h_d\), the corresponding point with Cartesian coordinates \( (c_0,c_1,\ldots,c_{d-1})\) can be computed. Note that homogeneous coordinates are not unique. For \( \lambda\ne 0\), the tuples \( (h_0,h_1,\ldots,h_d)\) and \( (\lambda\cdot h_0,\lambda\cdot h_1,\ldots,\lambda\cdot h_d)\) represent the same point. For a point with Cartesian coordinates \( (c_0,c_1,\ldots,c_{d-1})\) a possible homogeneous representation is \( (c_0,c_1,\ldots,c_{d-1},1)\). Homogeneous coordinates in fact allow to represent objects in a more general space, the projective space \( \mathbb{P}^d\). In CGAL we do not compute in projective geometry. Rather, we use homogeneous coordinates to avoid division operations, since the additional coordinate can serve as a common denominator.

Genericity Through Parameterization

Almost all the kernel objects (and the corresponding functions) are templates with a parameter that allows the user to choose the representation of the kernel objects. A type that is used as an argument for this parameter must fulfill certain requirements on syntax and semantics. The list of requirements defines an abstract kernel concept. For all kernel objects types, the types CGAL::Type<Kernel> and Kernel::Type are identical.

CGAL offers four families of concrete models for the concept Kernel, two based on the Cartesian representation of points and two based on the homogeneous representation of points. The interface of the kernel objects is designed such that it works well with both Cartesian and homogeneous representation. For example, points in 2D have a constructor with three arguments as well (the three homogeneous coordinates of the point). The common interfaces parameterized with a kernel class allow one to develop code independent of the chosen representation. We said "families" of models, because both families are parameterized too. A user can choose the number type used to represent the coordinates.

For reasons that will become evident later, a kernel class provides two typenames for number types, namely Kernel::FT and Kernel::RT. The type Kernel::FT must fulfill the requirements on what is called a FieldNumberType in CGAL. This roughly means that Kernel::FT is a type for which operations \( +\), \( -\), \( *\) and \( /\) are defined with semantics (approximately) corresponding to those of a field in a mathematical sense. Note that, strictly speaking, the built-in type int does not fulfill the requirements on a field type, since ints correspond to elements of a ring rather than a field, especially operation \( /\) is not the inverse of \( *\). The requirements on the type Kernel::RT are weaker. This type must fulfill the requirements on what is called a RingNumberType in CGAL. This roughly means that Kernel::RT is a type for which operations \( +\), \( -\), \( *\) are defined with semantics (approximately) corresponding to those of a ring in a mathematical sense.

Cartesian Kernels

With Cartesian<FieldNumberType> you can choose a Cartesian representation of coordinates. When you choose Cartesian representation you have to declare at the same time the type of the coordinates. A number type used with the Cartesian representation class should be a FieldNumberType as described above. As mentioned above, the built-in type int is not a FieldNumberType. However, for some computations with Cartesian representation, no division operation is needed, i.e., a RingNumberType is sufficient in this case. With Cartesian<FieldNumberType>, both Cartesian<FieldNumberType>::FT and Cartesian<FieldNumberType>::RT are mapped to FieldNumberType.

Cartesian<FieldNumberType> uses reference counting internally to save copying costs. CGAL also provides Simple_cartesian<FieldNumberType>, a kernel that uses Cartesian representation but no reference counting. Debugging is easier with Simple_cartesian<FieldNumberType>, since the coordinates are stored within the class and hence direct access to the coordinates is possible. Depending on the algorithm, it can also be slightly more or less efficient than Cartesian<FieldNumberType>. Again, in Simple_cartesian<FieldNumberType> both Simple_cartesian<FieldNumberType>::FT and Simple_cartesian<FieldNumberType>::RT are mapped to FieldNumberType.

Homogeneous Kernels

Homogeneous coordinates permit to avoid division operations in numerical computations, since the additional coordinate can serve as a common denominator. Avoiding divisions can be useful for exact geometric computation. With Homogeneous<RingNumberType> you can choose a homogeneous representation for the coordinates of the kernel objects. As for the Cartesian representation, one has to declare the type used to store the coordinates. Since the homogeneous representation does not use divisions, the number type associated with a homogeneous representation class must be a model for the weaker concept RingNumberType only. However, some operations provided by this kernel involve divisions, for example computing squared distances or Cartesian coordinates. To keep the requirements on the number type parameter of Homogeneous low, the number type Quotient<RingNumberType> is used for operations that require divisions. This number type can be viewed as an adaptor which turns a RingNumberType into a FieldNumberType. It maintains numbers as quotients, i.e., a numerator and a denominator. With Homogeneous<RingNumberType>, Homogeneous<RingNumberType>::FT is equal to Quotient<RingNumberType>, while Homogeneous<RingNumberType>::RT is equal to RingNumberType.

Homogeneous<RingNumberType> uses reference counting internally to save copying costs. CGAL also provides Simple_homogeneous<RingNumberType>, a kernel that uses homogeneous representation but no reference counting. Debugging is easier with Simple_homogeneous<RingNumberType>, since the coordinates are stored within the class and hence direct access to the coordinates is possible. Depending on the algorithm, it can also be slightly more or less efficient than Homogeneous<RingNumberType>. Again, in Simple_homogeneous<RingNumberType> the type Simple_homogeneous<RingNumberType>::FT is equal to Quotient<RingNumberType> while Simple_homogeneous<RingNumberType>::RT is equal to RingNumberType.

Naming Conventions

The use of kernel classes not only avoids problems, it also makes all CGAL classes very uniform. They always consist of:

  1. The capitalized base name of the geometric object, such as Point, Segment, or Triangle.

  2. An underscore followed by the dimension of the object, for example \( \_2\), \( \_3\), or \( \_d\).

  3. A kernel class as parameter, which itself is parameterized with a number type, such as Cartesian<double> or Homogeneous<leda_integer>.

Kernel as a Traits Class

Algorithms and data structures in the basic library of CGAL are parameterized by a traits class that subsumes the objects on which the algorithm or data structure operates as well as the operations to do so. For most of the algorithms and data structures in the basic library you can use a kernel as a traits class. For some algorithms you even do not have to specify the kernel; it is detected automatically using the types of the geometric objects passed to the algorithm. In some other cases, the algorithms or data structures needs more than is provided by the kernel concept. In these cases, a kernel can not be used as a traits class.

Choosing a Kernel and Predefined Kernels

If you start with integral Cartesian coordinates, many geometric computations will involve integral numerical values only. Especially, this is true for geometric computations that evaluate only predicates, which are tantamount to determinant computations. Examples are triangulation of point sets and convex hull computation. In this case, the Cartesian representation is probably the first choice, even with a ring type. You might use limited precision integer types like int or long, use double to present your integers (they have more bits in their mantissa than an int and overflow nicely), or an arbitrary precision integer type like the wrapper Gmpz for the GMP integers, leda_integer, or MP_Float. Note, that unless you use an arbitrary precision ring type, incorrect results might arise due to overflow.

If new points are to be constructed, for example the intersection point of two lines, computation of Cartesian coordinates usually involves divisions. Hence, one needs to use a FieldNumberType with Cartesian representation, or alternatively, switch to homogeneous representation. The type double is a - though imprecise - model for FieldNumberType. You can also put any RingNumberType into the Quotient adaptor to get a field type which then can be put into Cartesian. But using homogeneous representation on the RingNumberType is usually the better option. Other valid FieldNumberTypes are leda_rational and leda_real.

If it is crucial for you that the computation is reliable, the right choice is probably a number type that guarantees exact computation. The Filtered_kernel provides a way to apply filtering techniques [1] to achieve a kernel with exact and efficient predicates. Still other people will prefer the built-in type double, because they need speed and can live with approximate results, or even algorithms that, from time to time, crash or compute incorrect results due to accumulated rounding errors.

Predefined Kernels

For the user's convenience, CGAL provides 3 typedefs to generally useful kernels.

Kernel Geometry

Points and Vectors

In CGAL we strictly distinguish between points, vectors and directions. A point is a point in the Euclidean space \( \E^d\), a vector is the difference of two points \( p_2\), \( p_1\) and denotes the direction and the distance from \( p_1\) to \( p_2\) in the vector space \( \mathbb{R}^d\), and a direction is a vector where we forget about its length. They are different mathematical concepts. For example, they behave different under affine transformations and an addition of two points is meaningless in affine geometry. By putting them in different classes we not only get cleaner code, but also type checking by the compiler which avoids ambiguous expressions. Hence, it pays twice to make this distinction.

CGAL defines a symbolic constant ORIGIN of type Origin which denotes the point at the origin. This constant is used in the conversion between points and vectors. Subtracting it from a point \( p\) results in the locus vector of \( p\).

Cartesian<double>::Point_2 p(1.0, 1.0), q;
Cartesian<double>::Vector_2 v;
v = p - ORIGIN;
q = ORIGIN + v;
assert( p == q );

In order to obtain the point corresponding to a vector \( v\) you simply have to add \( v\) to ORIGIN. If you want to determine the point \( q\) in the middle between two points \( p_1\) and \( p_2\), you can writeyou might call midpoint(p_1,p_2) instead.

q = p_1 + (p_2 - p_1) / 2.0;

Note that these constructions do not involve any performance overhead for the conversion with the currently available representation classes.

Kernel Objects

Besides points (Kernel::Point_2, Kernel::Point_3), vectors (Kernel::Vector_2, Kernel::Vector_3), and directions (Kernel::Direction_2, Kernel::Direction_3), CGAL provides lines, rays, segments, planes, triangles, tetrahedra, iso-rectangles, iso-cuboids, circles and spheres.

Lines (Kernel::Line_2, Kernel::Line_3) in CGAL are oriented. In two-dimensional space, they induce a partition of the plane into a positive side and a negative side. Any two points on a line induce an orientation of this line. A ray (Kernel::Ray_2, Kernel::Ray_3) is semi-infinite interval on a line, and this line is oriented from the finite endpoint of this interval towards any other point in this interval. A segment (Kernel::Segment_2, Kernel::Segment_3) is a bounded interval on a directed line, and the endpoints are ordered so that they induce the same direction as that of the line.

Planes are affine subspaces of dimension two in \( \E^3\), passing through three points, or a point and a line, ray, or segment. CGAL provides a correspondence between any plane in the ambient space \( \E^3\) and the embedding of \( \E^2\) in that space. Just like lines, planes are oriented and partition space into a positive side and a negative side. In CGAL, there are no special classes for half-spaces. Half-spaces in 2D and 3D are supposed to be represented by oriented lines and planes, respectively.

Concerning polygons and polyhedra, the kernel provides triangles, iso-oriented rectangles, iso-oriented cuboids and tetrahedra. More complex polygonsAny sequence of points can be seen as a (not necessary simple) polygon or polyline. This view is used frequently in the basic library as well. and polyhedra or polyhedral surfaces can be obtained from the basic library (Polygon_2, Polyhedron_3), so they are not part of the kernel. As with any Jordan curves, triangles, iso-oriented rectangles and circles separate the plane into two regions, one bounded and one unbounded.

Orientation and Relative Position

Geometric objects in CGAL have member functions that test the position of a point relative to the object. Full dimensional objects and their boundaries are represented by the same type, e.g. half-spaces and hyperplanes are not distinguished, neither are balls and spheres and discs and circles. Such objects split the ambient space into two full-dimensional parts, a bounded part and an unbounded part (e.g. circles), or two unbounded parts (e.g. hyperplanes). By default these objects are oriented, i.e., one of the resulting parts is called the positive side, the other one is called the negative side. Both of these may be unbounded.

These objects have a member function oriented_side() that determines whether a test point is on the positive side, the negative side, or on the oriented boundary. These function returns a value of type Oriented_side.

Those objects that split the space in a bounded and an unbounded part, have a member function bounded_side() with return type Bounded_side.

If an object is lower dimensional, e.g. a triangle in three-dimensional space or a segment in two-dimensional space, there is only a test whether a point belongs to the object or not. This member function, which takes a point as an argument and returns a Boolean value, is called has_on().

Predicates and Constructions

Predicates

Predicates are at the heart of a geometry kernel. They are basic units for the composition of geometric algorithms and encapsulate decisions. Hence their correctness is crucial for the control flow and hence for the correctness of an implementation of a geometric algorithm. CGAL uses the term predicate in a generalized sense. Not only components returning a Boolean value are called predicates but also components returning an enumeration type like a Comparison_result or an Orientation. We say components, because predicates are implemented both as functions and function objects (provided by a kernel class).

CGAL provides predicates for the orientation of point sets (orientation(), left_turn(), right_turn(), collinear(), coplanar()), for comparing points according to some given order, especially for comparing Cartesian coordinates (e.g. lexicographically_xy_smaller()), in-circle and in-sphere tests, and predicates to compare distances.

Constructions

Functions and function objects that generate objects that are neither of type bool nor enum types are called constructions. Constructions involve computation of new numerical values and may be imprecise due to rounding errors unless a kernel with an exact number type is used.

Affine transformations (Kernel::Aff_transformation_2, Kernel::Aff_transformation_3) allow to generate new object instances under arbitrary affine transformations. These transformations include translations, rotations (in 2D only) and scaling. Most of the geometric objects in a kernel have a member function transform(Aff_transformation t) which applies the transformation to the object instance.

CGAL also provides a set of functions that detect or compute the intersection between objects of the 2D kernel, and many objects in the 3D kernel, and functions to calculate their squared distance. Moreover, some member functions of kernel objects are constructions.

So there are routines that compute the square of the Euclidean distance, but no routines that compute the distance itself. Why? First of all, the two values can be derived from each other quite easily (by taking the square root or taking the square). So, supplying only the one and not the other is only a minor inconvenience for the user. Second, often either value can be used. This is for example the case when (squared) distances are compared. Third, the library wants to stimulate the use of the squared distance instead of the distance. The squared distance can be computed in more cases and the computation is cheaper. We do this by not providing the perhaps more natural routine, The problem of a distance routine is that it needs the sqrt operation. This has two drawbacks:

  • The sqrt operation can be costly. Even if it is not very costly for a specific number type and platform, avoiding it is always cheaper.
  • There are number types on which no sqrt operation is defined, especially integer types and rationals.

Intersections and Variant Return Types

Some functions, for example intersection(), can return different types of objects. To achieve this in a type-safe way CGAL uses return values of type boost::optional< boost::variant< T... > > were T... is a list of all possible resulting geometric objects. The exact result type of an intersection can be determined through the metafunction cpp11::result_of<Kernel::Intersect_2(Type1, Type2)> or cpp11::result_of<Kernel::Intersect_3(Type1, Type2)>, where Type1 and Type2 are the types of the objects used in the intersection computation.

Example

In the following example, result_of is used to query the type of the return value for the intersection computation:

typedef Cartesian<double> K;
typedef K::Point_2 Point_2;
typedef K::Segment_2 Segment_2;
Segment_2 segment_1, segment_2;
std::cin >> segment_1 >> segment_2;
/* C++11 */
// auto v = intersection(segment_1, segment_2);
/* C++03 */
v = intersection(segment_1, segment_2);
if(v) {
/* not empty */
if (const Point_2 *p = boost::get<Point_2>(&*v) ) {
/* do something with *p */
} else {
const Segment_2 *s = boost::get<Segment_2>(&*v);
/* do something with *s */
}
} else {
/* empty intersection */
}

Constructive Predicates

For testing where a point p lies with respect to a plane defined by three points q, r and s, one may be tempted to construct the plane Kernel::Plane_3(q,r,s) and use the method oriented_side(p). This may pay off if many tests with respect to the plane are made. Nevertheless, unless the number type is exact, the constructed plane is only approximated, and round-off errors may lead oriented_side(p) to return an orientation which is different from the real orientation of p, q, r, and s.

In CGAL, we provide predicates in which such geometric decisions are made directly with a reference to the input points p, q, r, s, without an intermediary object like a plane. For the above test, the recommended way to get the result is to use orientation(p,q,r,s). For exact number types, the situation is different. If several tests are to be made with the same plane, it pays off to construct the plane and to use oriented_side(p).

Extensible Kernel

This manual section describe how users can plug user defined geometric classes in existing CGAL kernels. This is best illustrated by an example.

Introduction

CGAL defines the concept of a geometry kernel. Such a kernel provides types, construction objects and generalized predicates. Most implementations of Computational Geometry algorithms and data structures in the basic library of CGAL were done in a way that classes or functions can be parametrized with a geometric traits class.

In most cases this geometric traits class must be a model of the CGAL geometry kernel concept (but there are some exceptions).

An Extensive Example

Assume we have the following point class, where the coordinates are stored in an array of doubles, where we have another data member color, which shows up in the constructor.


File Kernel_23/MyPointC2.h

#ifndef MY_POINTC2_H
#define MY_POINTC2_H
#include <CGAL/Origin.h>
#include <CGAL/Bbox_2.h>
class MyPointC2 {
private:
double vec[2];
int col;
public:
MyPointC2()
: col(0)
{
*vec = 0;
*(vec+1) = 0;
}
MyPointC2(const double x, const double y, int c = 0)
: col(c)
{
*vec = x;
*(vec+1) = y;
}
const double& x() const { return *vec; }
const double& y() const { return *(vec+1); }
double & x() { return *vec; }
double& y() { return *(vec+1); }
int color() const { return col; }
int& color() { return col; }
bool operator==(const MyPointC2 &p) const
{
return ( *vec == *(p.vec) ) && ( *(vec+1) == *(p.vec + 1) && ( col == p.col) );
}
bool operator!=(const MyPointC2 &p) const
{
return !(*this == p);
}
};
#endif // MY_POINTC2_H

As said earlier the class is pretty minimalistic, for example it has no bbox() method. One might assume that a basic library algorithm which computes a bounding box (e.g, to compute the bounding box of a polygon), will not compile. Luckily it will, because it does not use of member functions of geometric objects, but it makes use of the functor Kernel::Construct_bbox_2.

To make the right thing happen with MyPointC2 we have to provide the following functor.


File Kernel_23/MyConstruct_bbox_2.h

#ifndef MYCONSTRUCT_BBOX_2_H
#define MYCONSTRUCT_BBOX_2_H
template <class ConstructBbox_2>
class MyConstruct_bbox_2 : public ConstructBbox_2 {
public:
using ConstructBbox_2::operator();
CGAL::Bbox_2 operator()(const MyPointC2& p) const {
return CGAL::Bbox_2(p.x(), p.y(), p.x(), p.y());
}
};
#endif //MYCONSTRUCT_BBOX_2_H

Things are similar for random access to the Cartesian coordinates of a point. As the coordinates are stored in an array of doubles we can use double* as random access iterator.


File Kernel_23/MyConstruct_coord_iterator.h

#ifndef MYCONSTRUCT_COORD_ITERATOR_H
#define MYCONSTRUCT_COORD_ITERATOR_H
class MyConstruct_coord_iterator {
public:
const double* operator()(const MyPointC2& p)
{
return &p.x();
}
const double* operator()(const MyPointC2& p, int)
{
const double* pyptr = &p.y();
pyptr++;
return pyptr;
}
};
#endif //MYCONSTRUCT_COORD_ITERATOR_H

The last functor we have to provide is the one which constructs points. That is you are not forced to add the constructor with the Origin as parameter to your class, nor the constructor with homogeneous coordinates. The functor is a kind of glue layer between the CGAL algorithms and your class.


File Kernel_23/MyConstruct_point_2.h

#ifndef MYCONSTRUCT_POINT_2_H
#define MYCONSTRUCT_POINT_2_H
template <typename K, typename OldK>
class MyConstruct_point_2
{
typedef typename K::RT RT;
typedef typename K::Point_2 Point_2;
typedef typename K::Line_2 Line_2;
typedef typename Point_2::Rep Rep;
public:
typedef Point_2 result_type;
// Note : the CGAL::Return_base_tag is really internal CGAL stuff.
// Unfortunately it is needed for optimizing away copy-constructions,
// due to current lack of delegating constructors in the C++ standard.
Rep // Point_2
operator()(CGAL::Return_base_tag, CGAL::Origin o) const
{ return Rep(o); }
Rep // Point_2
operator()(CGAL::Return_base_tag, const RT& x, const RT& y) const
{ return Rep(x, y); }
Rep // Point_2
operator()(CGAL::Return_base_tag, const RT& x, const RT& y, const RT& w) const
{ return Rep(x, y, w); }
Point_2
operator()(const CGAL::Origin&) const
{ return MyPointC2(0, 0, 0); }
Point_2
operator()(const RT& x, const RT& y) const
{
return MyPointC2(x, y, 0);
}
const Point_2&
operator()(const Point_2 & p) const
{
return p;
}
Point_2
operator()(const Line_2& l) const
{
typename OldK::Construct_point_2 base_operator;
Point_2 p = base_operator(l);
return p;
}
Point_2
operator()(const Line_2& l, int i) const
{
typename OldK::Construct_point_2 base_operator;
return base_operator(l, i);
}
// We need this one, as such a functor is in the Filtered_kernel
Point_2
operator()(const RT& x, const RT& y, const RT& w) const
{
if(w != 1){
return MyPointC2(x/w, y/w, 0);
} else {
return MyPointC2(x,y, 0);
}
}
};
#endif //MYCONSTRUCT_POINT_2_H

Now we are ready to put the puzzle together. We won't explain it in detail, but you see that there are typedefs to the new point class and the functors. All the other types are inherited.


File Kernel_23/MyKernel.h

#ifndef MYKERNEL_H
#define MYKERNEL_H
#include <CGAL/Cartesian.h>
#include "MyPointC2.h"
#include "MySegmentC2.h"
#include "MyConstruct_bbox_2.h"
#include "MyConstruct_coord_iterator.h"
#include "MyConstruct_point_2.h"
// K_ is the new kernel, and K_Base is the old kernel
template < typename K_, typename K_Base >
class MyCartesian_base
: public K_Base::template Base<K_>::Type
{
typedef typename K_Base::template Base<K_>::Type OldK;
public:
typedef K_ Kernel;
typedef MyPointC2 Point_2;
typedef MySegmentC2<Kernel> Segment_2;
typedef MyConstruct_point_2<Kernel, OldK> Construct_point_2;
typedef const double* Cartesian_const_iterator_2;
typedef MyConstruct_coord_iterator Construct_cartesian_const_iterator_2;
typedef MyConstruct_bbox_2<typename OldK::Construct_bbox_2>
Construct_bbox_2;
Construct_point_2
construct_point_2_object() const
{ return Construct_point_2(); }
Construct_bbox_2
construct_bbox_2_object() const
{ return Construct_bbox_2(); }
Construct_cartesian_const_iterator_2
construct_cartesian_const_iterator_2_object() const
{ return Construct_cartesian_const_iterator_2(); }
template < typename Kernel2 >
struct Base { typedef MyCartesian_base<Kernel2, K_Base> Type; };
};
template < typename FT_ >
struct MyKernel
: public CGAL::Type_equality_wrapper<
MyCartesian_base<MyKernel<FT_>, CGAL::Cartesian<FT_> >,
MyKernel<FT_> >
{};
#endif // MYKERNEL_H

Finally, we give an example how this new kernel can be used. Predicates and constructions work with the new point, they can be a used to construct segments and triangles with, and data structures from the Basic Library, as the Delaunay triangulation work with them.

The kernel itself can be made robust by plugging it in the Filtered_kernel.


File Kernel_23/MyKernel.cpp

#include <CGAL/basic.h>
#include <CGAL/Filtered_kernel.h>
#include <CGAL/Delaunay_triangulation_2.h>
#include <CGAL/squared_distance_2.h>
#include <cassert>
#include "MyKernel.h"
#include "MyPointC2_iostream.h"
typedef MyKernel<double> MK;
typedef CGAL::Delaunay_triangulation_2<K> Delaunay_triangulation_2;
typedef K::Point_2 Point;
typedef K::Segment_2 Segment;
typedef K::Ray_2 Ray;
typedef K::Line_2 Line;
typedef K::Triangle_2 Triangle;
typedef K::Iso_rectangle_2 Iso_rectangle;
const int RED= 1;
const int BLACK=2;
int main()
{
Point a(0,0), b(1,0), c(1,1), d(0,1);
a.color()=RED;
b.color()=BLACK;
d.color()=RED;
Delaunay_triangulation_2 dt;
dt.insert(a);
K::Orientation_2 orientation;
orientation(a,b,c);
Point p(1,2), q;
p.color() = RED;
q.color() = BLACK;
std::cout << p << std::endl;
K::Compute_squared_distance_2 squared_distance;
std::cout << "squared_distance(a, b) == "
<< squared_distance(a, b) << std::endl;
Segment s1(p,q), s2(a, c);
K::Construct_midpoint_2 construct_midpoint_2;
Point mp = construct_midpoint_2(p,q);
std::cout << "midpoint(" << p << " , " << q << ") == " << mp << std::endl;
assert(s1.source().color() == RED);
K::Intersect_2 intersection;
intersect = intersection(s1, s2);
K::Construct_cartesian_const_iterator_2 construct_it;
K::Cartesian_const_iterator_2 cit = construct_it(a);
assert(*cit == a.x());
cit = construct_it(a,0);
cit--;
assert(*cit == a.y());
Line l1(a,b), l2(p, q);
intersection(l1, l2);
intersection(s1, l1);
Ray r1(d,b), r2(d,c);
intersection(r1, r2);
intersection(r1, l1);
Triangle t1(a,b,c), t2(a,c,d);
intersection(t1, t2);
intersection(t1, l1);
intersection(t1, s1);
intersection(t1, r1);
Iso_rectangle i1(a,c), i2(d,p);
intersection(i1, i2);
intersection(i1, s1);
intersection(i1, r1);
intersection(i1, l1);
t1.orientation();
std::cout << s1.source() << std::endl;
std::cout << t1.bbox() << std::endl;
std::cout << "done" << std::endl;
return 0;
}

Limitations

The point class must have member functions x() and y() (and z() for the 3d point). We will probably introduce function objects that take care of coordinate access.

As we enforce type equality between MyKernel::Point_2 and Point_2<MyKernel>, the constructor with the color as third argument is not available.

Projection Traits Classes

It is sometimes useful to apply 2D algorithms to the projection of 3D points on a plane. Examples are triangulated terrains, which are points with elevation, or surface reconstruction from parallel slices, where one wants to check the simplicity or orientation of polygons.

For this purpose CGAL provides several projection traits classes, which are a model of traits class concepts of 2D triangulations, 2D polygon and 2D convex hull traits classes. The projection traits classes are listed in the "Is Model Of" sections of the concepts.

Design and Implementation History

At a meeting at Utrecht University in January 1995, Olivier Devillers, Andreas Fabri, Wolfgang Freiseisen, Geert-Jan Giezeman, Mark Overmars, Stefan Schirra, Otfried Schwarzkopf (now Otfried Cheong), and Sven Schönherr discussed the foundations of the CGAL kernel. Many design and software engineering issues were addressed, e.g. naming conventions, coupling of classes (flat versus deep class hierarchy), memory allocation, programming conventions, mutability of atomic objects, points and vectors, storing additional information, orthogonality of operations on the kernel objects, viewing non-constant-size objects like polygons as dynamic data structures (and hence not as part of the (innermost) kernel).

The people attending the meeting delegated the compilation of a draft specification to Stefan Schirra. The resulting draft specification was intentionally modeled on CGAL's precursors C++gal and Plageo as well as on the geometric part of LEDA. The specification already featured coexistence of Cartesian and homogeneous representation of point/vector data and parameterization by number type(s). During the discussion of the draft a kernel design group was formed. The members of this group were Andreas Fabri, Geert-Jan Giezeman, Lutz Kettner, Stefan Schirra, and Sven Schönherr. The work of the kernel design group led to significant changes and improvements of the original design, e.g. the strong separation between points and vectors. Probably the most important enhancement was the design of a common superstructure for the previously uncoupled Cartesian and homogeneous representations. One can say, that the kernel was designed by this group. The kernel was later revised based on suggestions by Hervé Brönnimann, Bernd Gärtner, Michael Hoffmann, and Lutz Kettner.

A first version of the kernel was internally made available at the beginning of the CGAL-project (esprit ltr iv project number 21957). Since then many more people contributed to the evolution of the kernel through discussions on the CGAL mailing lists. The implementation based on Cartesian representation was (initially) provided by Andreas Fabri, the homogeneous representation (initially) by Stefan Schirra. Intersection and distance computations were implemented by Geert-Jan Giezeman. Further work has been done by Susan Hert on the overall maintenance of the kernel. Philippe Guigue has provided efficient intersection tests for 3D triangles. Andreas Fabri, Michael Hoffmann and Sylvain Pion have improved the support for the extensibility and adaptability of the kernel. Pedro Machado Manhães de Castro and Monique Teillaud introduced 3D circles. In 2010, Pierre Alliez, Stéphane Tayeb and Camille Wormser added intersection constructions for 3D triangles and efficient intersection tests for bounding boxes.

Acknowledgment

This work was supported by the Graduiertenkolleg 'Algorithmische Diskrete Mathematik', under grant DFG We 1265/2-1, and by ESPRIT IV Long Term Research Projects No. 21957 (CGAL) and No. 28155 (GALIA).