Graphite Version 3
An experimental 3D geometry processing program
|
Constrained Delaunay Triangulation with vertices that are exact points. Can be used to implement 2D CSG. More...
#include <geogram/delaunay/CDT_2d.h>
Public Types | |
typedef exact::vec2h | ExactPoint |
Public Member Functions | |
ExactCDT2d () | |
ExactCDT2d constructor. | |
~ExactCDT2d () override | |
ExactCDT2d destructor. | |
void | clear () override |
Removes everything from this triangulation. | |
index_t | insert (const ExactPoint &p, index_t id=0, index_t hint=index_t(-1)) |
Inserts a point. | |
void | insert_constraint (index_t v1, index_t v2, index_t operand_bits=0) |
Inserts a constraint. | |
void | create_enclosing_quad (const ExactPoint &p1, const ExactPoint &p2, const ExactPoint &p3, const ExactPoint &p4) |
Creates a first large enclosing quad. | |
void | create_enclosing_rectangle (double x1, double y1, double x2, double y2) |
Creates a first large enclosing rectangle. | |
const ExactPoint & | point (index_t v) const |
Gets a point by vertex index. | |
index_t | vertex_id (index_t v) const |
Gets a vertex id by vertex index. | |
void | set_vertex_id (index_t v, index_t id) |
Sets a vertex id by vertex index. | |
void | classify_triangles (const std::string &boolean_expression, bool mark_only=false) |
Used by 2D CSG operations, discards triangles according to a boolean operation. | |
void | save (const std::string &filename) const override |
![]() | |
CDTBase2d () | |
CDTBase2d constructor. | |
virtual | ~CDTBase2d () |
CDTBase2d destructor. | |
void | insert_constraint (index_t i, index_t j) |
Inserts a constraint. | |
void | remove_external_triangles (bool remove_internal_holes=false) |
Recursively removes all the triangles adjacent to the border, and keeps what's surrounded by constraints. | |
void | set_delaunay (bool delaunay) |
Specifies whether a constrained Delaunay triangulation should be constructed, or just a plain constrained triangulation. | |
index_t | nT () const |
Gets the number of triangles. | |
index_t | nv () const |
Gets the number of vertices. | |
index_t | ncnstr () const |
Gets the number of constraints. | |
index_t | Tv (index_t t, index_t lv) const |
Gets a vertex of a triangle. | |
index_t | Tv_find (index_t t, index_t v) const |
Finds the local index of a vertex in a triangle. | |
index_t | Tadj (index_t t, index_t le) const |
Gets a triangle adjacent to a triangle. | |
index_t | Tadj_find (index_t t1, index_t t2) const |
Finds the edge accross which a triangle is adjacent to another one. | |
index_t | vT (index_t v) const |
Gets a triangle incident to a given vertex. | |
index_t | Tedge_cnstr_first (index_t t, index_t le) const |
Gets the constraint associated with an edge. | |
index_t | edge_cnstr_next (index_t ecit) const |
Gets the successor of an edge constraint iterator. | |
index_t | edge_cnstr (index_t ecit) const |
Gets an edge constraint from an edge constraint iterator. | |
index_t | Tedge_cnstr_nb (index_t t, index_t le) const |
Gets the number of constraints associated with a triange edge. | |
bool | Tedge_is_Delaunay (index_t t, index_t le) const |
Tests whether a triangle edge is Delaunay. | |
void | check_consistency () const |
Checks both combinatorics and geometry, aborts on unconsistency. | |
Protected Member Functions | |
void | add_point (const ExactPoint &p, index_t id=index_t(-1)) |
void | begin_insert_transaction () override |
void | commit_insert_transaction () override |
void | rollback_insert_transaction () override |
Sign | orient2d (index_t i, index_t j, index_t k) const override |
Sign | incircle (index_t i, index_t j, index_t k, index_t l) const override |
Incircle predicate. | |
index_t | create_intersection (index_t E1, index_t i, index_t j, index_t E2, index_t k, index_t l) override |
Given two segments that have an intersection, create the intersection. | |
![]() | |
index_t | insert (index_t v, index_t hint=index_t(-1)) |
Inserts a new point. | |
void | create_enclosing_triangle (index_t v1, index_t v2, index_t v3) |
Creates the combinatorics for a first large enclosing triangle. | |
void | create_enclosing_quad (index_t v1, index_t v2, index_t v3, index_t v4) |
Creates the combinatorics for a first large enclosing quad. | |
void | Tset_flag (index_t t, index_t flag) |
Sets a triangle flag. | |
void | Treset_flag (index_t t, index_t flag) |
Resets a triangle flag. | |
bool | Tflag_is_set (index_t t, index_t flag) |
Tests a triangle flag. | |
bool | Tis_in_list (index_t t) const |
Tests whether a triangle is in a DList. | |
void | insert_vertex_in_edge (index_t v, index_t t, index_t le, DList &S) |
Inserts a vertex in an edge. | |
void | insert_vertex_in_edge (index_t v, index_t t, index_t le) |
Inserts a vertex in an edge. | |
void | insert_vertex_in_triangle (index_t v, index_t t, DList &S) |
Inserts a vertex in a triangle. | |
index_t | find_intersected_edges (index_t i, index_t j, DList &Q) |
Finds the edges intersected by a constraint. | |
void | walk_constraint_v (CDT2d_ConstraintWalker &W) |
Used by find_intersected_edges() | |
void | walk_constraint_t (CDT2d_ConstraintWalker &W, DList &Q) |
Used by find_intersected_edges() | |
void | constrain_edges (index_t i, index_t j, DList &Q, DList &N) |
Constrains an edge by iteratively flipping the intersected edges. | |
void | Delaunayize_vertex_neighbors (index_t from_v) |
Restores Delaunay condition starting from the triangles incident to a given vertex. | |
void | Delaunayize_vertex_neighbors (index_t v, DList &S) |
Restores Delaunay condition starting from the triangles incident to a given vertex. | |
void | Delaunayize_new_edges (DList &N) |
Restores Delaunay condition for a set of edges after inserting a constrained edge. | |
void | Tset (index_t t, index_t v1, index_t v2, index_t v3, index_t adj1, index_t adj2, index_t adj3, index_t e1cnstr=index_t(-1), index_t e2cnstr=index_t(-1), index_t e3cnstr=index_t(-1)) |
Sets all the combinatorial information of a triangle and edge flags. | |
void | Trot (index_t t, index_t lv) |
Rotates indices in triangle t in such a way that a given vertex becomes vertex 0. | |
void | swap_edge (index_t t1, bool swap_t1_t2=false) |
Swaps an edge. | |
void | Tadj_set (index_t t, index_t le, index_t adj) |
Sets a triangle adjacency relation. | |
index_t | Topp (index_t t, index_t e=0) const |
Gets the neighboring triangle vertex opposite to a given vertex. | |
void | Tadj_back_connect (index_t t1, index_t le1, index_t prev_t2_adj_e2) |
After having changed connections from triangle to a neighbor, creates connections from neighbor to triangle. | |
index_t | Tnew () |
Creates a new triangle. | |
void | Tset_edge_cnstr_first (index_t t, index_t le, index_t ecit) |
Sets the constraints list associated with an edge. | |
void | Tadd_edge_cnstr (index_t t, index_t le, index_t cnstr_id) |
Adds a constraint to a triangle edge. | |
void | Tadd_edge_cnstr_with_neighbor (index_t t, index_t le, index_t cnstr_id) |
Adds a constraint to a triangle edge and to the neighboring edge if it exists. | |
bool | Tedge_is_constrained (index_t t, index_t le) const |
Tests whether an edge is constrained. | |
void | for_each_T_around_v (index_t v, std::function< bool(index_t t, index_t lv)> doit) |
Calls a user-defined function for each triangle around a vertex. | |
index_t | locate (index_t v, index_t hint=index_t(-1), Sign *orient=nullptr) const |
Locates a vertex. | |
bool | is_convex_quad (index_t t) const |
Tests whether triange t and its neighbor accross edge 0 form a strictly convex quad. | |
void | remove_marked_triangles () |
Removes all the triangles that have the flag T_MARKED_FLAG set. | |
void | Tcheck (index_t t) const |
Consistency check for a triangle. | |
void | debug_Tcheck (index_t t) const |
Consistency check for a triangle in debug mode, ignored in release mode. | |
void | check_combinatorics () const |
Consistency combinatorial check for all the triangles. | |
void | debug_check_combinatorics () const |
Consistency combinatorial check for all the triangles in debug mode, ignored in release mode. | |
virtual void | check_geometry () const |
Consistency geometrical check for all the triangles. | |
void | debug_check_geometry () const |
Consistency geometrical check for all the triangles in debug mode, ignored in release mode. | |
void | debug_check_consistency () const |
Checks both combinatorics and geometry in debug mode, ignored in release mode, aborts on unconsistency. | |
bool | segment_segment_intersect (index_t u1, index_t u2, index_t v1, index_t v2) const |
Tests whether two segments have a frank intersection. | |
bool | segment_edge_intersect (index_t v1, index_t v2, index_t t, index_t le) const |
Tests whether an edge triangle and a segment have a frank intersection. | |
void | check_edge_intersections (index_t v1, index_t v2, const DList &Q) |
Checks that the edges stored in a DList exactly correspond to all edge intersections between a segment and the triangle edges. | |
index_t | eT (Edge E) |
Gets a triangle incident a a given edge. | |
index_t | locate_naive (index_t v, index_t hint=index_t(-1), Sign *orient=nullptr) const |
Simpler version of locate() kept for reference. | |
void | constrain_edges_naive (index_t i, index_t j, DList &Q, vector< Edge > &N) |
Simpler version of constrain_edges() kept for reference. | |
void | Delaunayize_new_edges_naive (vector< Edge > &N) |
Simpler version of Delaunayize_new_edges() that uses a vector instead of a DList, kept for reference. | |
Protected Attributes | |
vector< ExactPoint > | point_ |
vector< double > | length_ |
vector< index_t > | id_ |
vector< index_t > | cnstr_operand_bits_ |
vector< index_t > | facet_inclusion_bits_ |
std::map< trindex, Sign > | pred_cache_ |
bool | use_pred_cache_insert_buffer_ |
std::vector< std::pair< trindex, Sign > > | pred_cache_insert_buffer_ |
vector< bindex > | constraints_ |
![]() | |
index_t | nv_ |
index_t | ncnstr_ |
vector< index_t > | T_ |
vector< index_t > | Tadj_ |
vector< index_t > | v2T_ |
vector< uint8_t > | Tflags_ |
vector< index_t > | Tecnstr_first_ |
vector< index_t > | ecnstr_val_ |
vector< index_t > | ecnstr_next_ |
vector< index_t > | Tnext_ |
vector< index_t > | Tprev_ |
bool | delaunay_ |
Sign | orient_012_ |
bool | exact_incircle_ |
bool | exact_intersections_ |
Additional Inherited Members | |
![]() | |
enum | { DLIST_S_ID =0 , DLIST_Q_ID =1 , DLIST_N_ID =2 , DLIST_NB =3 } |
Constants for list_id. More... | |
enum | { T_MARKED_FLAG = DLIST_NB , T_VISITED_FLAG = DLIST_NB+1 } |
Constants for triangle flags. More... | |
typedef std::pair< index_t, index_t > | Edge |
![]() | |
static index_t | find_3 (const index_t *T, index_t v) |
Finds the index of an integer in an array of three integers. | |
Constrained Delaunay Triangulation with vertices that are exact points. Can be used to implement 2D CSG.
Points are represented using exact 2d homogeneous coordinates. Unlike CDT2d, this ensures exact representation of constraints intersections with guaranteed behavior. Under the hood, it inherits CDTBase2d (constrained Delaunay triangulation), and redefines orient2d(), incircle2d() and create_intersection() using vectors with homogeneous coordinates stored as arithmetic expansions (vec2HE) or arbitrary-precision floating point numbers (vec2HEx) if compiled with Tessael's geogramplus extension package.
|
overrideprotectedvirtual |
Reimplemented from GEO::CDTBase2d.
void GEO::ExactCDT2d::classify_triangles | ( | const std::string & | boolean_expression, |
bool | mark_only = false |
||
) |
Used by 2D CSG operations, discards triangles according to a boolean operation.
Discards all the triangles that are outside the object defined by the boolean expression. It uses the operand bits associated with the constraints.
[in] | boolean_expression | a string with the boolean expression, as defined by BooleanExpression constructor. Each variable corresponds to an operand bit associated with the constraints. There can be up to 32 operands. |
mark_only | if set, triangles to be discarded are marked (but not discarded) |
|
overridevirtual |
Removes everything from this triangulation.
Reimplemented from GEO::CDTBase2d.
|
overrideprotectedvirtual |
Reimplemented from GEO::CDTBase2d.
void GEO::ExactCDT2d::create_enclosing_quad | ( | const ExactPoint & | p1, |
const ExactPoint & | p2, | ||
const ExactPoint & | p3, | ||
const ExactPoint & | p4 | ||
) |
Creates a first large enclosing quad.
[in] | p1,p2,p3,p4 | the four vertices of the quad |
The quad needs to be convex. create_enclosing_rectangle() or create_enclosing_quad() need to be called before anything else
|
inline |
Creates a first large enclosing rectangle.
[in] | x1,y1,x2,y2 | rectangle bounds |
create_enclosing_triangle(), create_enclosing_rectangle() or create_enclosing_quad() need to be called before anything else
|
overrideprotectedvirtual |
Given two segments that have an intersection, create the intersection.
The intersection is given both as the indices of segment extremities (i,j) and (k,l), that one can use to retreive the points in derived classes, and constraint indices E1 and E2, that derived classes may use to retreive symbolic information attached to the constraint
[in] | E1 | the index of the first edge, corresponding to the value of ncnstr() when insert_constraint() was called for that edge |
[in] | i,j | the vertices of the first segment |
[in] | E2 | the index of the second edge, corresponding to the value of ncnstr() when insert_constraint() was called for that edge |
[in] | k,l | the vertices of the second segment |
i
, j
] and [k
, l
] Implements GEO::CDTBase2d.
|
overrideprotectedvirtual |
Incircle predicate.
[in] | i,j,k | the three vertices of a triangle |
[in] | l | another vertex |
POSITIVE | if l is inside the circumscribed circle of the triangle |
ZERO | if l is on the circumscribed circle of the triangle |
NEGATIVE | if l is outside the circumscribed circle of the triangle |
Implements GEO::CDTBase2d.
index_t GEO::ExactCDT2d::insert | ( | const ExactPoint & | p, |
index_t | id = 0 , |
||
index_t | hint = index_t(-1) |
||
) |
Inserts a point.
[in] | p | the point to be inserted |
[in] | hint | a triangle not too far away from the point to be inserted |
[in] | id | an opaque identifier attached to the vertex that can be used by client code for instance to keep relations with a mesh. It can be queried using the vertex_id() function. |
Inserts a constraint.
[in] | v1,v2 | the two extremities of the constraint, as returned by insert() |
[in] | operand_bits | optional bitfield used by 2D CSG, indicating on which primitive boundaries the constraint is. Each bit set corresponds to a primitive. It is used by the classify() function. |
Implements GEO::CDTBase2d.
|
inline |
|
overrideprotectedvirtual |
Reimplemented from GEO::CDTBase2d.
|
overridevirtual |
Implements GEO::CDTBase2d.
|
protected |
|
protected |