Hervé Brönnimann, Andreas Fabri, Geert-Jan Giezeman, Susan Hert, Michael Hoffmann, Lutz Kettner, Sylvain Pion, and Stefan Schirra
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.
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.
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<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<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.
The use of kernel classes not only avoids problems, it also makes all Cgal classes very uniform. They always consist of:
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 [BBP01] 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.
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.
Point_2< Cartesian<double> > p(1.0, 1.0), q; Vector_2< Cartesian<double> > 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 p1 and p2, you can write2
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.
Lines (Line_2<Kernel>, Line_3<Kernel>) 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 (Ray_2<Kernel>, Ray_3<Kernel>) 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 (Segment_2<Kernel>, Segment_3<Kernel>) 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 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 3 and the embedding of 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 polygons3 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.
For these objects there is a 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()
Cgal provides predicates for the orientation of point sets (orientation, leftturn, rightturn, 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.
Affine transformations (Aff_transformation_2<Kernel>, Aff_transformation_3<Kernel>) 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:
{ 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; Object obj = intersection(segment_1, segment_2); if (const Point_2 *point = object_cast<Point_2>(&obj)) { /* do something with *point */ } else if (const Segment_2 *segment = object_cast<Segment_2>(&obj)) { /* do something with *segment*/ }
/* there was no intersection */ }
The
intersection
routine itself looks roughly as follows:
template < class Kernel > Object intersection(Segment_2<Kernel> s1, Segment_2<Kernel> s2) {
if (/* intersection in a point */ ) {
Point_2<Kernel> p = ... ; return make_object(p);
} else if (/* intersection in a segment */ ) {
Segment_2<Kernel> s = ... ; return make_object(s); } return Object(); }
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).
This manual section describe how users can plug user defined geometric classes in existing Cgal kernels. This is best illustrated by an example.
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).
Assume you 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: examples/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: examples/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: examples/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: examples/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()(CGAL::Origin o) const { return MyPointC2(0, 0, 0); } Point_2 operator()(const RT& x, const RT& y) const { return MyPointC2(x, y, 0); } 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: examples/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: examples/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::Filtered_kernel_adaptor<MK> K; 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); assert(s1.source().color() == RED); K::Intersect_2 intersection; CGAL::Object o = 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); squared_distance(r1, r2); squared_distance(r1, l2); squared_distance(r1, s2); 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; }
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.
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 for the Concepts'' sections of the concepts.
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.
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).
1 | Currently it requires having either LEDA or CORE installed |
2 | you might call midpoint(p_1,p_2) instead |
3 | Any 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. |