ttl Namespace Reference

Main interface to TTL. More...

Functions

Delaunay Triangulation



template<class TraitsType , class DartType , class PointType >
bool insertNode (DartType &dart, PointType &point)
 Inserts a new node in an existing Delaunay triangulation and swaps edges to obtain a new Delaunay triangulation.
template<class TraitsType , class DartType >
void removeRectangularBoundary (DartType &dart)
 Removes the rectangular boundary of a triangulation as a final step of an incremental Delaunay triangulation.
template<class TraitsType , class DartType >
void removeNode (DartType &dart)
 Removes the node associated with dart and updates the triangulation to be Delaunay.
template<class TraitsType , class DartType >
void removeBoundaryNode (DartType &dart)
 Removes the boundary node associated with dart and updates the triangulation to be Delaunay.
template<class TraitsType , class DartType >
void removeInteriorNode (DartType &dart)
 Removes the interior node associated with dart and updates the triangulation to be Delaunay.
template<class TraitsType , class ForwardIterator , class DartType >
void insertNodes (ForwardIterator first, ForwardIterator last, DartType &dart)
Topological and Geometric Queries



template<class TraitsType , class PointType , class DartType >
bool locateFaceSimplest (const PointType &point, DartType &dart)
 Locates the face containing a given point.
template<class TraitsType , class PointType , class DartType >
bool locateTriangle (const PointType &point, DartType &dart)
 Locates the triangle containing a given point.
template<class TraitsType , class PointType , class DartType >
bool inTriangleSimplest (const PointType &point, const DartType &dart)
 Checks if point is inside the triangle associated with dart.
template<class TraitsType , class PointType , class DartType >
bool inTriangle (const PointType &point, const DartType &dart)
 Checks if point is inside the triangle associated with dart.
template<class DartType , class DartListType >
void getBoundary (const DartType &dart, DartListType &boundary)
 Gets the boundary as sequence of darts, where the edges associated with the darts are boundary edges, given a dart with an associating edge at the boundary of a topology structure.
template<class DartType >
bool isBoundaryEdge (const DartType &dart)
 Checks if the edge associated with dart is at the boundary of the triangulation.
template<class DartType >
bool isBoundaryFace (const DartType &dart)
 Checks if the face associated with dart is at the boundary of the triangulation.
template<class DartType >
bool isBoundaryNode (const DartType &dart)
 Checks if the node associated with dart is at the boundary of the triangulation.
template<class DartType >
int getDegreeOfNode (const DartType &dart)
 Returns the degree of the node associated with dart.
template<class DartType , class DartListType >
void get_0_orbit_interior (const DartType &dart, DartListType &orbit)
 Gets the 0-orbit around an interior node.
template<class DartType , class DartListType >
void get_0_orbit_boundary (const DartType &dart, DartListType &orbit)
 Gets the 0-orbit around a node at the boundary.
template<class DartType >
bool same_0_orbit (const DartType &d1, const DartType &d2)
 Checks if the two darts belong to the same 0-orbit, i.e., if they share a node.
template<class DartType >
bool same_1_orbit (const DartType &d1, const DartType &d2)
 Checks if the two darts belong to the same 1-orbit, i.e., if they share an edge.
template<class DartType >
bool same_2_orbit (const DartType &d1, const DartType &d2)
 Checks if the two darts belong to the same 2-orbit, i.e., if they lie in the same triangle.
template<class TraitsType , class DartType >
bool swappableEdge (const DartType &dart, bool allowDegeneracy)
 Checks if the edge associated with dart is swappable, i.e., if the edge is a diagonal in a strictly convex (or convex) quadrilateral.
template<class DartType >
void positionAtNextBoundaryEdge (DartType &dart)
 Given a dart, CCW or CW, positioned in a 0-orbit at the boundary of a tessellation.
template<class TraitsType , class DartType >
bool convexBoundary (const DartType &dart)
 Checks if the boundary of a triangulation is convex.
template<class TopologyElementType , class DartType >
bool isMemberOfFace (const TopologyElementType &topologyElement, const DartType &dart)
template<class TraitsType , class NodeType , class DartType >
bool locateFaceWithNode (const NodeType &node, DartType &dart_iter)
template<class DartType >
void getAdjacentTriangles (const DartType &dart, DartType &t1, DartType &t2, DartType &t3)
template<class DartType >
void getNeighborNodes (const DartType &dart, list< DartType > &node_list, bool &boundary)
template<class TraitsType , class DartType >
bool degenerateTriangle (const DartType &dart)
Utilities for Delaunay Triangulation



template<class TraitsType , class DartType , class DartListType >
void optimizeDelaunay (DartListType &elist)
 Optimizes the edges in the given sequence according to the Delaunay criterion, i.e., such that the edge will fullfill the circumcircle criterion (or equivalently the MaxMin angle criterion) with respect to the quadrilaterals where they are diagonals.
template<class TraitsType , class DartType , class DartListType >
void optimizeDelaunay (DartListType &elist, const typename DartListType::iterator end)
template<class TraitsType , class DartType >
bool swapTestDelaunay (const DartType &dart, bool cycling_check)
 Checks if the edge associated with dart should be swapped according to the Delaunay criterion, i.e., the circumcircle criterion (or equivalently the MaxMin angle criterion).
template<class TraitsType , class DartType >
void recSwapDelaunay (DartType &diagonal)
 Recursively swaps edges in the triangulation according to the Delaunay criterion.
template<class TraitsType , class DartType , class ListType >
void swapEdgesAwayFromInteriorNode (DartType &dart, ListType &swapped_edges)
 Swaps edges away from the (interior) node associated with dart such that that exactly three edges remain incident with the node.
template<class TraitsType , class DartType , class ListType >
void swapEdgesAwayFromBoundaryNode (DartType &dart, ListType &swapped_edges)
 Swaps edges away from the (boundary) node associated with dart in such a way that when removing the edges that remain incident with the node, the boundary of the triangulation will be convex.
template<class TraitsType , class DartType , class DartListType >
void swapEdgeInList (const typename DartListType::iterator &it, DartListType &elist)
 Swap the the edge associated with iterator it and update affected darts in elist accordingly.
Constrained (Delaunay) Triangulation



template<class TraitsType , class DartType >
DartType insertConstraint (DartType &dstart, DartType &dend, bool optimize_delaunay)
 Inserts a constrained edge between two existing nodes in a triangulation.

Detailed Description

Main interface to TTL.

This namespace contains the basic generic algorithms for the TTL, the Triangulation Template Library.

Examples of functionality are:

General requirements and assumptions:
  • DartType and TraitsType should be implemented in accordance with the description in Application Programming Interface to TTL (API).
  • A "Requires:" section in the documentation of a function template shows which functionality is required in TraitsType to support that specific function.
    Functionalty required in DartType is the same (almost) for all function templates; see Application Programming Interface to TTL (API) and the example referred to.
  • When a reference to a dart object is passed to a function in TTL, it is assumed that it is oriented counterclockwise (CCW) in a triangle unless it is explicitly mentioned that it can also be clockwise (CW). The same applies for a dart that is passed from a function in TTL to the users TraitsType class (or struct).
  • When an edge (represented with a dart) is swapped, it is assumed that darts outside the quadrilateral where the edge is a diagonal are not affected by the swap. Thus, TraitsType::swapEdge must be implemented in accordance with this rule.
Glossary:
See also:
ttl_util and Application Programming Interface to TTL (API)
Author:
´┐Żyvind Hjelle, oyvindhj@ifi.uio.no

Function Documentation

template<class TraitsType , class DartType >
bool ttl::convexBoundary ( const DartType &  dart  )  [inline]

Checks if the boundary of a triangulation is convex.

Parameters:
dart A CCW dart at the boundary of the triangulation

Definition at line 1355 of file ttl.h.

References ttl_util::crossProduct2d(), and getBoundary().

01355                                               {
01356     
01357     list<DartType> blist;
01358     ttl::getBoundary(dart, blist);
01359     
01360     int no;
01361     no = blist.size();
01362     typename list<DartType>::const_iterator bit = blist.begin();
01363     DartType d1 = *bit;
01364     ++bit;
01365     DartType d2;
01366     bool convex = true;
01367     for (; bit != blist.end(); ++bit) {
01368       d2 = *bit;
01369       double crossProd = TraitsType::crossProduct2d(d1, d2);
01370       if (crossProd < 0.0) {
01371         //cout << "!!! Boundary is NOT convex: crossProd = " << crossProd << endl;
01372         convex = false;
01373         return convex;
01374       }
01375       d1 = d2;
01376     }
01377     
01378     // Check the last angle
01379     d2 = *blist.begin();
01380     double crossProd = TraitsType::crossProduct2d(d1, d2);
01381     if (crossProd < 0.0) {
01382       //cout << "!!! Boundary is NOT convex: crossProd = " << crossProd << endl;
01383       convex = false;
01384     }
01385     
01386     //if (convex)
01387     //  cout << "\n---> Boundary is convex\n" << endl;
01388     //cout << endl;
01389     return convex;
01390   }

template<class TraitsType , class DartType >
bool ttl::degenerateTriangle ( const DartType &  dart  )  [inline]

Definition at line 1249 of file ttl.h.

References ttl_util::crossProduct2d().

01249                                                   {
01250     
01251     // Check if triangle is degenerate
01252     // Assumes CCW dart
01253     
01254     DartType d1 = dart;
01255     DartType d2 = d1;
01256     d2.alpha1();
01257     if (TraitsType::crossProduct2d(d1,d2) == 0)
01258       return true;
01259     
01260     return false;
01261   }

template<class DartType , class DartListType >
void ttl::get_0_orbit_boundary ( const DartType &  dart,
DartListType &  orbit 
) [inline]

Gets the 0-orbit around a node at the boundary.

Parameters:
dart A dart (CCW or CW) positioned at a boundary node and at a boundary edge.
Return values:
orbit Sequence of darts with one orbit for each arc. All the darts, exept the last one, have the same orientation (CCW or CW) as dart, and dart is the first element in the sequence.
  • DartListType::push_back (DartType&)
Note:
  • The last dart in the sequence have opposite orientation compared to the others!
See also:
ttl::get_0_orbit_interior

Definition at line 1158 of file ttl.h.

01158                                                                          {
01159     
01160     DartType dart_prev;
01161     DartType d_iter = dart;
01162     do {
01163       orbit.push_back(d_iter);
01164       d_iter.alpha1();
01165       dart_prev = d_iter;
01166       d_iter.alpha2();
01167     } while (d_iter != dart_prev);
01168     
01169     orbit.push_back(d_iter); // the last one with opposite orientation
01170   }

template<class DartType , class DartListType >
void ttl::get_0_orbit_interior ( const DartType &  dart,
DartListType &  orbit 
) [inline]

Gets the 0-orbit around an interior node.

Parameters:
dart A dart (CCW or CW) positioned at an interior node.
Return values:
orbit Sequence of darts with one orbit for each arc. All the darts have the same orientation (CCW or CW) as dart, and dart is the first element in the sequence.
  • DartListType::push_back (DartType&)
See also:
ttl::get_0_orbit_boundary

Definition at line 1124 of file ttl.h.

01124                                                                          {
01125     
01126     DartType d_iter = dart;
01127     orbit.push_back(d_iter);
01128     d_iter.alpha1().alpha2();
01129     
01130     while (d_iter != dart) {
01131       orbit.push_back(d_iter);
01132       d_iter.alpha1().alpha2();
01133     }
01134   }

template<class DartType >
void ttl::getAdjacentTriangles ( const DartType &  dart,
DartType &  t1,
DartType &  t2,
DartType &  t3 
) [inline]

Definition at line 831 of file ttl.h.

00831                                                                                               {
00832     
00833     DartType dart_iter = dart;
00834     
00835     // add first
00836     if (dart_iter.alpha2() != dart) {
00837       t1 = dart_iter;
00838       dart_iter = dart;
00839     }
00840     
00841     // add second
00842     dart_iter.alpha0();
00843     dart_iter.alpha1();
00844     DartType dart_prev = dart_iter;
00845     if ((dart_iter.alpha2()) != dart_prev) {
00846       t2 = dart_iter;
00847       dart_iter = dart_prev;
00848     }
00849     
00850     // add third
00851     dart_iter.alpha0();
00852     dart_iter.alpha1();
00853     dart_prev = dart_iter;
00854     if ((dart_iter.alpha2()) != dart_prev)
00855       t3 = dart_iter;
00856   }

template<class DartType , class DartListType >
void ttl::getBoundary ( const DartType &  dart,
DartListType &  boundary 
) [inline]

Gets the boundary as sequence of darts, where the edges associated with the darts are boundary edges, given a dart with an associating edge at the boundary of a topology structure.

The first dart in the sequence will be the given one, and the others will have the same orientation (CCW or CW) as the first. Assumes that the given dart is at the boundary.

Parameters:
dart A dart at the boundary (CCW or CW)
boundary A sequence of darts, where the associated edges are the boundary edges
  • DartListType::push_back (DartType&)

Definition at line 876 of file ttl.h.

References positionAtNextBoundaryEdge().

Referenced by convexBoundary().

00876                                                                    {
00877     // assumes the given dart is at the boundary (by edge)
00878     
00879     DartType dart_iter(dart);
00880     boundary.push_back(dart_iter); // Given dart as first element
00881     dart_iter.alpha0();
00882     positionAtNextBoundaryEdge(dart_iter);
00883     
00884     while (dart_iter != dart) {
00885       boundary.push_back(dart_iter);
00886       dart_iter.alpha0();
00887       positionAtNextBoundaryEdge(dart_iter);
00888     }
00889   }

template<class DartType >
int ttl::getDegreeOfNode ( const DartType &  dart  )  [inline]

Returns the degree of the node associated with dart.

Definition:
The degree (or valency) of a node V in a triangulation, is defined as the number of edges incident with V, i.e., the number of edges joining V with another node in the triangulation.

Definition at line 999 of file ttl.h.

Referenced by swapEdgesAwayFromInteriorNode().

00999                                               {
01000     
01001     DartType dart_iter(dart);
01002     DartType dart_prev;
01003     
01004     // If input dart is reached again, then interior node
01005     // If alpha2(d)=d, then boundary
01006     
01007     int degree = 0;
01008     bool boundaryVisited = false;
01009     do {
01010       dart_iter.alpha1();
01011       degree++;
01012       dart_prev = dart_iter;
01013       
01014       dart_iter.alpha2();
01015       
01016       if (dart_iter == dart_prev) {
01017         if (!boundaryVisited) {
01018           boundaryVisited = true;
01019           // boundary is reached first time, count in the reversed direction
01020           degree++; // count the start since it is not done above
01021           dart_iter = dart;
01022           dart_iter.alpha2();
01023         }
01024         else
01025           return degree;
01026       }
01027       
01028     } while (dart_iter != dart);
01029     
01030     return degree;
01031   }

template<class DartType >
void ttl::getNeighborNodes ( const DartType &  dart,
list< DartType > &  node_list,
bool &  boundary 
) [inline]

Definition at line 1059 of file ttl.h.

01059                                                                                            {
01060     
01061     DartType dart_iter(dart);
01062     
01063     dart_iter.alpha0(); // position the dart at an opposite node
01064     
01065     DartType dart_prev = dart_iter;
01066     
01067     bool start_at_boundary = false;
01068     dart_iter.alpha2();
01069     if (dart_iter == dart_prev)
01070       start_at_boundary = true;
01071     else
01072       dart_iter = dart_prev; // back again 
01073     
01074     DartType dart_start = dart_iter;
01075     
01076     do {
01077       node_list.push_back(dart_iter);
01078       dart_iter.alpha1();
01079       dart_iter.alpha0();
01080       dart_iter.alpha1();
01081       dart_prev = dart_iter;
01082       dart_iter.alpha2();
01083       if (dart_iter == dart_prev) {
01084         // boundary reached
01085         boundary = true;
01086         if (start_at_boundary == true) {
01087           // add the dart which now is positioned at the opposite boundary
01088           node_list.push_back(dart_iter);
01089           return;
01090         }
01091         else {
01092           // call the function again such that we start at the boundary
01093           // first clear the list and reposition to the initial node
01094           dart_iter.alpha0();
01095           node_list.clear();
01096           getNeighborNodes(dart_iter, node_list, boundary);
01097           return; // after one recursive step
01098         }
01099       }
01100     } while (dart_iter != dart_start);
01101     
01102     boundary = false;
01103   }

template<class TraitsType , class DartType >
DartType ttl::insertConstraint ( DartType &  dstart,
DartType &  dend,
bool  optimize_delaunay 
) [inline]

Inserts a constrained edge between two existing nodes in a triangulation.

If the constraint falls on one or more existing nodes and this is detected by the predicate TraitsType::orient2d, which should return zero in this case, the constraint is split. Otherwise a degenerate triangle will be made along the constraint.

Parameters:
dstart A CCW dart with the start node of the constraint as the source node
dend A CCW dart with the end node of the constraint as the source node
optimize_delaunay If set to true, the resulting triangulation will be a constrained Delaunay triangulation. If set to false, the resulting triangulation will not necessarily be of constrained Delaunay type.
Return values:
DartType A dart representing the constrained edge.
Assumes:
  • The constrained edge must be inside the existing triangulation (and it cannot cross the boundary of the triangulation).

Definition at line 543 of file ttl_constr.h.

References same_0_orbit().

00543                                                                                         {
00544     
00545     // Assumes:
00546     // - It is the users responsibility to avoid crossing constraints
00547     // - The constraint cannot cross the boundary, i.e., the boundary must be
00548     //   convex in the area of crossing edges.
00549     // - dtart and dend are preserved (same node associated.)
00550     
00551 
00552     // Find edges crossing the constraint and put them in elist.
00553     // If findCrossingEdges reaches a Node lying