From 5640b84b2be76f01ac08ef5b6d41bc8c6fd495f1 Mon Sep 17 00:00:00 2001 From: Thijs van Lankveld Date: Fri, 13 Jun 2014 14:42:51 +0200 Subject: [PATCH] First setup for the segment traverser small feature. --- .../conforming_gabriel_3.cpp | 105 -- .../Conforming_Delaunay_triangulation_3.h | 1268 ----------------- .../Conforming_triangulation_cell_base_3.h | 179 --- .../Conforming_triangulation_vertex_base_3.h | 69 - .../CGAL/Delaunay_triangulation_utils_3.h | 406 ------ .../include/CGAL/Encroaching_collecter_3.h | 135 -- .../CGAL/Triangulation_segment_traverser_3.h | 695 --------- .../Conforming_triangulation_3/copyright | 2 - .../description.txt | 3 - .../Conforming_triangulation_3/license.txt | 1 - .../Conforming_triangulation_3/maintainer | 1 - .../doc/Triangulation_3/examples.txt | 1 + .../examples/Triangulation_3/README | 5 + .../Triangulation_3/segment_traverser_3.cpp | 45 + .../include/CGAL/Triangulation_3.h | 81 +- .../CGAL/Triangulation_segment_traverser_3.h | 356 +++++ .../Triangulation_segment_traverser_3_impl.h | 665 +++++++++ .../include/CGAL/Triangulation_utils_3.h | 90 +- .../include/CGAL/internal/Dummy_gt.h | 37 + 19 files changed, 1244 insertions(+), 2900 deletions(-) delete mode 100644 Conforming_triangulation_3/example/Conforming_triangulation_3/conforming_gabriel_3.cpp delete mode 100644 Conforming_triangulation_3/include/CGAL/Conforming_Delaunay_triangulation_3.h delete mode 100644 Conforming_triangulation_3/include/CGAL/Conforming_triangulation_cell_base_3.h delete mode 100644 Conforming_triangulation_3/include/CGAL/Conforming_triangulation_vertex_base_3.h delete mode 100644 Conforming_triangulation_3/include/CGAL/Delaunay_triangulation_utils_3.h delete mode 100644 Conforming_triangulation_3/include/CGAL/Encroaching_collecter_3.h delete mode 100644 Conforming_triangulation_3/include/CGAL/Triangulation_segment_traverser_3.h delete mode 100644 Conforming_triangulation_3/package_info/Conforming_triangulation_3/copyright delete mode 100644 Conforming_triangulation_3/package_info/Conforming_triangulation_3/description.txt delete mode 100644 Conforming_triangulation_3/package_info/Conforming_triangulation_3/license.txt delete mode 100644 Conforming_triangulation_3/package_info/Conforming_triangulation_3/maintainer create mode 100644 Triangulation_3/examples/Triangulation_3/segment_traverser_3.cpp create mode 100644 Triangulation_3/include/CGAL/Triangulation_segment_traverser_3.h create mode 100644 Triangulation_3/include/CGAL/Triangulation_segment_traverser_3_impl.h create mode 100644 Triangulation_3/include/CGAL/internal/Dummy_gt.h diff --git a/Conforming_triangulation_3/example/Conforming_triangulation_3/conforming_gabriel_3.cpp b/Conforming_triangulation_3/example/Conforming_triangulation_3/conforming_gabriel_3.cpp deleted file mode 100644 index 103d12e13e9..00000000000 --- a/Conforming_triangulation_3/example/Conforming_triangulation_3/conforming_gabriel_3.cpp +++ /dev/null @@ -1,105 +0,0 @@ -#include -#include -#include - - -typedef CGAL::Exact_predicates_exact_constructions_kernel K; -typedef K::Point_3 Point_3; - -typedef CGAL::Triangulation_vertex_base_with_info_3 > Vb; -typedef CGAL::Conforming_triangulation_cell_base_3 Cb; -typedef CGAL::Triangulation_data_structure_3 < Vb, Cb > Tds; - -typedef CGAL::Conforming_Delaunay_triangulation_3< K, Tds > CDT3; - - -template -std::ptrdiff_t insert_with_info(Triangulation& t, InputIterator first,InputIterator last) -{ - typedef typename Triangulation::Geom_traits Geom_traits; - typedef typename Triangulation::Vertex_handle Vertex_handle; - - std::size_t n = t.number_of_vertices(); - std::vector indices; - std::vector points; - std::vector infos; - std::ptrdiff_t index=0; - for (InputIterator it=first;it!=last;++it){ - std::pair value=*it; - points.push_back( value.first ); - infos.push_back ( value.second ); - indices.push_back(index++); - } - - typedef CGAL::Spatial_sort_traits_adapter_3 Search_traits; - - CGAL::spatial_sort(indices.begin(),indices.end(),Search_traits(&(points[0]),Geom_traits())); - - Vertex_handle hint; - for (typename std::vector::const_iterator - it = indices.begin(), end = indices.end(); - it != end; ++it){ - hint = t.insert(points[*it], hint); - if (hint!=Vertex_handle()) hint->info()=infos[*it]; - } - - return t.number_of_vertices() - n; -} - -int main() -{ - - std::vector< std::pair > points; - points.push_back( std::make_pair(Point_3(0,0,0),0) ); - points.push_back( std::make_pair(Point_3(1,0,0),1) ); - points.push_back( std::make_pair(Point_3(0,1,0),2) ); - points.push_back( std::make_pair(Point_3(0.2,0.2,100),3) ); - points.push_back( std::make_pair(Point_3(0.2,0.2,-100),4) ); - - - std::vector< std::pair > constraints; - constraints.push_back( std::make_pair(0,1) ); - constraints.push_back( std::make_pair(0,2) ); - constraints.push_back( std::make_pair(1,2) ); - constraints.push_back( std::make_pair(3,4) ); - - CDT3 cdt3; - insert_with_info( cdt3, points.begin(), points.end() ); - - - std::vector< CDT3::Vertex_handle > vertices( points.size() ); - - for (CDT3::Finite_vertices_iterator vit=cdt3.finite_vertices_begin(), - vit_end=cdt3.finite_vertices_end();vit!=vit_end;++vit) - { - vertices[ vit->info() ]=vit; - } - - std::cout << cdt3.number_of_vertices() << std::endl; - - std::cout << "Making Delaunay conform" << std::endl; - for (std::vector< std::pair >::iterator cst_it=constraints.begin(), - cst_end=constraints.end(); - cst_it!=cst_end; ++cst_it) - { - cdt3.insert_conforming( std::make_pair(vertices[cst_it->first], vertices[cst_it->second]) ); - } - - - for(CDT3::Finite_edges_iterator it=cdt3.finite_edges_begin(); it!=cdt3.finite_edges_end(); ++it) - if ( it->first->vertex(it->second)->steiner() || it->first->vertex(it->third)->steiner() ) - std::cout << "edge with steiner endpoint" << std::endl; - - std::cout << "Making Gabriel conform" << std::endl; - for (std::vector< std::pair >::iterator cst_it=constraints.begin(), - cst_end=constraints.end(); - cst_it!=cst_end; ++cst_it) - { - cdt3.insert_conforming_Gabriel( vertices[cst_it->first], vertices[cst_it->second] ); - } - - for(CDT3::Finite_edges_iterator it=cdt3.finite_edges_begin(); it!=cdt3.finite_edges_end(); ++it) - if ( it->first->vertex(it->second)->steiner() || it->first->vertex(it->third)->steiner() ) - std::cout << "edge with steiner endpoint" << std::endl; -} - diff --git a/Conforming_triangulation_3/include/CGAL/Conforming_Delaunay_triangulation_3.h b/Conforming_triangulation_3/include/CGAL/Conforming_Delaunay_triangulation_3.h deleted file mode 100644 index 764c6e6b431..00000000000 --- a/Conforming_triangulation_3/include/CGAL/Conforming_Delaunay_triangulation_3.h +++ /dev/null @@ -1,1268 +0,0 @@ -//The conforming Delaunay triangulation structure. -//Copyright (C) 2012 Utrecht University -// -//This program is free software: you can redistribute it and/or modify -//it under the terms of the GNU General Public License as published by -//the Free Software Foundation, either version 3 of the License, or -//(at your option) any later version. -// -//This program is distributed in the hope that it will be useful, -//but WITHOUT ANY WARRANTY; without even the implied warranty of -//MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -//GNU General Public License for more details. -// -//You should have received a copy of the GNU General Public License -//along with this program. If not, see . -// -// Author(s): Thijs van Lankveld - -// Based on CGAL/Delaunay_triangulation_3.h -// -// For the love of something, why is the 3D triangulation package so badly written? - -#ifndef CGAL_CONFORMING_DELAUNAY_TRIANGULATION_3_H -#define CGAL_CONFORMING_DELAUNAY_TRIANGULATION_3_H - -#include - -#include -#include - -#include "Delaunay_triangulation_utils_3.h" -#include "Conforming_triangulation_cell_base_3.h" -#include "Conforming_triangulation_vertex_base_3.h" -#include "Encroaching_collecter_3.h" - -namespace CGAL { - -template < class Tr > class Natural_neighbors_3; - -template < class Gt, - class Tds_ = Triangulation_data_structure_3 < Conforming_triangulation_vertex_base_3, - Conforming_triangulation_cell_base_3 >, - class Itag = No_intersection_tag > -class Conforming_Delaunay_triangulation_3: public Delaunay_triangulation_utils_3 { - //typedef Triangulation_data_structure Tds; - typedef Tds_ Tds; - typedef Conforming_Delaunay_triangulation_3 cDT; - typedef Delaunay_triangulation_utils_3 DT; - typedef Triangulation_3 Tr; - - friend class Natural_neighbors_3; - typedef typename Gt::Line_3 Line; -public: - typedef Tds Triangulation_data_structure; - typedef Gt Geom_traits; - typedef Itag Intersection_tag; - - typedef typename Gt::FT FT; - - typedef typename Gt::Point_3 Point; - typedef typename Gt::Vector_3 Vector; - typedef typename Gt::Segment_3 Segment; - typedef typename Gt::Triangle_3 Triangle; - typedef typename Gt::Tetrahedron_3 Tetrahedron; - typedef typename Gt::Plane_3 Plane; - - typedef typename DT::Cell_handle Cell_handle; - typedef typename DT::Vertex_handle Vertex_handle; - - typedef typename DT::Cell Cell; - typedef typename DT::Vertex Vertex; - typedef typename DT::Facet Facet; - typedef typename DT::Edge Edge; - - typedef typename DT::Bi_vertex Bi_vertex; - - typedef std::pair Constraint_1; - - typedef typename DT::Cell_circulator Cell_circulator; - typedef typename DT::Facet_circulator Facet_circulator; - typedef typename DT::Cell_iterator Cell_iterator; - typedef typename DT::Facet_iterator Facet_iterator; - typedef typename DT::Edge_iterator Edge_iterator; - typedef typename DT::Vertex_iterator Vertex_iterator; - - typedef typename DT::Finite_vertices_iterator Finite_vertices_iterator; - typedef typename DT::Finite_cells_iterator Finite_cells_iterator; - typedef typename DT::Finite_facets_iterator Finite_facets_iterator; - typedef typename DT::Finite_edges_iterator Finite_edges_iterator; - - typedef typename DT::All_cells_iterator All_cells_iterator; - - typedef typename DT::Cell_traverser Cell_traverser; - - typedef typename DT::Locate_type Locate_type; - -private: - typedef std::list Point_list; - typedef typename Point_list::iterator Points_iterator; - typedef std::back_insert_iterator Points_output; - -public: - typedef Encroaching_collecter_3 Encroaching_collecter; - -#ifndef CGAL_CFG_USING_BASE_MEMBER_BUG_2 - using DT::cw; - using DT::ccw; - using DT::geom_traits; - using DT::number_of_vertices; - using DT::dimension; - using DT::finite_facets_begin; - using DT::finite_facets_end; - using DT::finite_vertices_begin; - using DT::finite_vertices_end; - using DT::finite_cells_begin; - using DT::finite_cells_end; - using DT::finite_edges_begin; - using DT::finite_edges_end; - using DT::tds; - using DT::infinite_vertex; - using DT::insert; - using DT::next_around_edge; - using DT::vertex_triple_index; - using DT::mirror_vertex; - using DT::move_point; - using DT::coplanar; - using DT::coplanar_orientation; - using DT::orientation; - using DT::is_infinite; - using DT::finite_incident_edges; - using DT::incident_cells; - using DT::insert_in_conflict; - using DT::is_valid_finite; -#endif - -protected: - // Visitor that should be run before and after inserting a point - // to maintain the conforming edges in the triangulation. - template < class Tr > - class Conforming_visitor { - typedef std::pair CBV; - - mutable std::list _conforming; - - protected: - Tr *_tr; - - public: - Conforming_visitor(Tr* tr): _tr(tr) {} - - // Store and remove the conforming edges from a collection of cells. - template - void process_cells_in_conflict(InputIterator start, InputIterator end) const { - for (InputIterator it = start; it != end; ++it) { - for (int i = 0; i < _tr->dimension(); ++i) { - for (int j = i+1; j < _tr->dimension()+1; ++j) { - if ((*it)->is_conforming(i, j)) { - _conforming.push_back(CBV(_tr->to_bi_vertex(*it, i, j), (*it)->is_marked(i, j))); - _tr->set_conforming(*it, i, j, false); - } - } - } - } - } - - // Reinsert the conforming edges into the triangulation. - void reinsert_vertices(Vertex_handle) { - while (!_conforming.empty()) { - CBV b = _conforming.front(); - _conforming.pop_front(); - if (b.second) - _tr->insert_marked(b.first); - else - _tr->insert_conforming(b.first); - } - } - - // Reinsert the conforming edges into the triangulation - // under the assumption that no Steiner points need to be added. - void reinsert_no_steiner() { - Cell_handle c; int li, lj; - while (!_conforming.empty()) { - CBV b = _conforming.front(); - _conforming.pop_front(); - CGAL_triangulation_assertion_code(bool ok =) - _tr->is_edge(b.first.first, b.first.second, c, li, lj); - CGAL_triangulation_assertion(ok); - - if (b.second) - _tr->mark_edge(c, li, lj); - else - _tr->set_conforming(c, li, lj, true); - } - } - - Vertex_handle replace_vertex(Cell_handle c, int index, const Point&) {return c->vertex(index);} - void hide_point(Cell_handle, const Point&) {} - }; // class Conforming_visitor - - friend class Conflict_tester_3; - friend class Conflict_tester_2; - friend class Encroaching_collecter_3; - - Conforming_visitor conforming_visitor; - -protected: - // PREDICATES - - // Checks if point c is encroaching upon the segment ab. - // A point is encroaching if it is inside the ball with diameter ab. - bool is_encroaching(const Point& a, const Point& b, const Point& c) const { - CGAL_triangulation_assertion(a != b); - CGAL_triangulation_assertion(a != c); - CGAL_triangulation_assertion(b != c); - return side_of_bounded_sphere(a, b, c) != ON_UNBOUNDED_SIDE; - } - -public: - // CONSTRUCTORS - Conforming_Delaunay_triangulation_3(const Gt& gt = Gt()): DT(gt), conforming_visitor(this) {} - - // Copy constructor duplicates vertices and cells. - Conforming_Delaunay_triangulation_3(const Conforming_Delaunay_triangulation_3& tr): DT(tr), conforming_visitor(this) { - CGAL_triangulation_postcondition(is_valid()); - } - - // Create a conforming 3D Delaunay triangulation from a number of constraints. - template < class InputIterator > - Conforming_Delaunay_triangulation_3(InputIterator begin, InputIterator end, const Gt& gt = Gt()): - DT(gt), conforming_visitor(this) { - insert_conforming(begin, end); - CGAL_triangulation_postcondition(is_valid()); - } - -protected: - // Compute the intersection of ab with a line or plane. - bool compute_intersection(const Line& ab, const Line& l, Point& p) const { - Object result = Gt().intersect_3_object()(ab, l); - return assign(p, result); - } - bool compute_intersection(const Line& ab, const Plane& pl, Point& p) const { - Object result = Gt().intersect_3_object()(ab, pl); - return assign(p, result); - } - - // Compute the intersection of ab with either a line (edge) or plane (facet). - bool compute_intersection(Locate_type lt, Cell_handle c, int li, int lj, Vertex_handle va, Vertex_handle vb, Point& p) const { - if (lt == DT::EDGE) - return compute_intersection(Line(va->point(), vb->point()), - Line(c->vertex(li)->point(), c->vertex(lj)->point()), p); - else - return compute_intersection(Line(va->point(), vb->point()), plane(c, li), p); - } - - // Insert the intersection of segment ab with either a facet or an edge. - Vertex_handle intersect(Cell_handle c, Locate_type lt, int li, int lj, Vertex_handle va, Vertex_handle vb) { - return intersect(c, lt, li, lj, va, vb, Itag());} - Vertex_handle intersect(Cell_handle c, Locate_type lt, int li, int lj, Vertex_handle va, Vertex_handle vb, No_intersection_tag); - Vertex_handle intersect(Cell_handle c, Locate_type lt, int li, int lj, Vertex_handle va, Vertex_handle vb, Exact_intersections_tag); - Vertex_handle intersect(Cell_handle c, Locate_type lt, int li, int lj, Vertex_handle va, Vertex_handle vb, Exact_predicates_tag); - - // Split segment ab based on a reference point. - Point construct_split_point(Vertex_handle va, Vertex_handle vb, const Point& ref); - -protected: - // Conform an edge to important encroaching points. - // The boolean indicates whether the edge should be made conforming Gabriel - template - void conform_segment(Vertex_handle va, Vertex_handle vb, Points_iterator encr_begin, Points_iterator encr_end); - - // Mark a segment that it is assumed to already be in the triangulation. - void mark_segment(Vertex_handle va, Vertex_handle vb); - // Insert a marked segment; - // this includes conforming the segment and marking its subsegments as conforming. - void insert_marked(Vertex_handle va, Vertex_handle vb); - void insert_marked(const Bi_vertex& m) {insert_marked(m.first, m.second);} - - // Mark an existing edge or facet as (un)constrained. - void set_conforming(Cell_handle c, int li, int lj, bool C); - void set_conforming(const Edge& e, bool C) {set_conforming(e.first, e.second, e.third, C);} - void mark_edge(Cell_handle c, int li, int lj); - void mark_edge(const Edge& e) {mark_edge(e.first, e.second, e.third);} - -protected: - // Restore the conforming edges of a facet based on the neighboring cell. - void restore_conforming(Cell_handle c, int li); - void restore_conforming(const Facet& f) {restore_conforming(f.first, f.second);} - -public: - // INSERT POINT - - virtual Vertex_handle insert(const Point& p, Locate_type lt, Cell_handle c, int li, int lj); - -public: - // INSERT CONFORMING - - // Insert a conforming segment. - void insert_conforming(Vertex_handle va, Vertex_handle vb) { - Point_list encr; conform_segment(va, vb, encr.begin(), encr.end());} - void insert_conforming_Gabriel(Vertex_handle va, Vertex_handle vb) { - Point_list encr; conform_segment(va, vb, encr.begin(), encr.end());} - void insert_conforming(const Bi_vertex& c) { - insert_conforming(c.first, c.second);} - void insert_conforming(const Point& a, const Point& b, Cell_handle hint = Cell_handle()); - void insert_conforming(const Constraint_1& c, Cell_handle hint = Cell_handle()) { - insert_conforming(c.first, c.second, hint);} - - // Insert a number of conforming segments, presented as a collection of pairs of Point_3. - template < class InputIterator > - int insert_conforming(InputIterator begin, InputIterator end, Cell_handle hint = Cell_handle()); - - // Insert a loop of conforming segments, presented as a collection of Point_3. - template < class InputIterator > - int insert_conforming_loop(InputIterator begin, InputIterator end, Cell_handle hint = Cell_handle()); - template < class InputIterator > - int insert_marked_loop(InputIterator begin, InputIterator end, Cell_handle hint = Cell_handle()); - template < class InputIterator > - int remove_conforming_loop(InputIterator begin, InputIterator end, Cell_handle hint = Cell_handle()); - -protected: - template < class InputIterator, class OutputIterator > - void insert_loop_points(InputIterator begin, InputIterator end, OutputIterator out, Cell_handle hint = Cell_handle()); - -public: - // REMOVE - - // Vertices that are incident to a conforming edge cannot be removed. - bool remove(Vertex_handle v); - - // Remove a conforming edge or segment. - bool remove_conforming(Cell_handle c, int li, int lj); - bool remove_conforming(const Edge& e) { - return remove_conforming(e.first, e.second, e.third);} - bool remove_conforming(Vertex_handle va, Vertex_handle vb); - bool remove_conforming(Bi_vertex c) { - return remove_conforming(c.first, c.second);} - -public: - // QUERY CONFORMING - - // Check if an edge is conforming. - bool is_conforming(Cell_handle c, int li, int lj) const {return c->is_conforming(li, lj);} - bool is_conforming(const Edge& e) const {return is_conforming(e.first, e.second, e.third);} - - // Check if an edge is marked. - bool is_marked(Cell_handle c, int li, int lj) const {return c->is_marked(li, lj);} - bool is_marked(const Edge& e) const {return is_marked(e.first, e.second, e.third);} - - // Check if a vertex is incident to a conforming edge. - bool are_there_incident_conforming(Vertex_handle v) const; - - // Check if the triangulation contains an edge ac on the segment ab, where c(li, lj) is that edge. - bool includes_edge(Vertex_handle va, Vertex_handle vb, Vertex_handle& vc, Cell_handle& c, int& li, int& lj) const; - - // Check if the triangulation is valid. - virtual bool is_valid(bool verbose = false, int level = 0) const; - virtual bool is_valid(Cell_handle c, bool verbose = false, int level = 0) const; - -private: - // Collect the points encroaching upon the segment ac, where c is the first point on the ray ab. - // encroaching must point to a Point. - // Returns true if a conforming edge or constrained facet was encountered and the intersection inserted. - template - bool collect_encroaching(Vertex_handle va, Vertex_handle vb, OutputIterator encroaching, Vertex_handle& vc); - - // Sort the points encroaching upon the segment ab to the from of the collection. - // begin and end must point to Point. - Points_iterator sort_encroaching(Points_iterator begin, Points_iterator end, const Point& a, const Point& b) const; -}; // Conforming_Delaunay_triangulation_3 - -template < class Gt, class Tds, class Itag > -typename Conforming_Delaunay_triangulation_3::Vertex_handle -Conforming_Delaunay_triangulation_3:: -intersect(Cell_handle /* c */, Locate_type /* lt */, int /* li */, int /* lj */, Vertex_handle /* va */, Vertex_handle /* vb */, No_intersection_tag) { - std::cerr << " sorry, this triangulation does not deal with" << std::endl - << " intersecting conforming edges" << std::endl; - CGAL_triangulation_assertion(false); - return Vertex_handle(); -} - -template < class Gt, class Tds, class Itag > -typename Conforming_Delaunay_triangulation_3::Vertex_handle -Conforming_Delaunay_triangulation_3:: -intersect(Cell_handle c, Locate_type lt, int li, int lj, Vertex_handle va, Vertex_handle vb, Exact_intersections_tag) { - // Compute the intersection of the segment ab and a facet or edge of cell c. - // Note that no further check is made whether this intersection actually occurs inside the segment(s) or triangle. - Point p; - CGAL_triangulation_assertion_code( bool ok = ) - compute_intersection(lt, c, li, lj, va, vb, p); - CGAL_triangulation_assertion(ok); - Vertex_handle v = insert(p, c); - v->steiner() = true; - return v; -} - -template < class Gt, class Tds, class Itag > -typename Conforming_Delaunay_triangulation_3::Vertex_handle -Conforming_Delaunay_triangulation_3:: -intersect(Cell_handle c, Locate_type lt, int li, int lj, Vertex_handle va, Vertex_handle vb, Exact_predicates_tag) { - // Compute the intersection of the segment ab and a facet or edge of cell c. - // Note that no further check is made whether this intersection actually occurs inside the segment(s) or triangle. - Point p; - bool ok = compute_intersection(lt, c, li, lj, va, vb, p); - if (ok) { - Vertex_handle v = insert(p, c); - v->steiner() = true; - return v; - } - - // Return the vertex closest to the intersection. - Vertex_handle closest; - switch (lt) { - case DT::EDGE: { - Line ab(va->point(), vb->point()); - Line ij(c->vertex(li)->point(), c->vertex(lj)->point()); - FT dist2 = squared_distance(va->point(), ij); - FT db = squared_distance(vb->point(), ij); - FT di = squared_distance(c->vertex(li)->point(), ab); - FT dj = squared_distance(c->vertex(lj)->point(), ab); - closest = va; - if (db < dist2) {dist2 = db; closest = vb;} - if (di < dist2) {dist2 = di; closest = c->vertex(li);} - if (dj < dist2) {dist2 = dj; closest = c->vertex(lj);} - } - case DT::FACET: { - Line ab(va->point(), vb->point()); - Plane p = plane(c, li); - FT dist2 = squared_distance(va->point(), p); - FT db = squared_distance(vb->point(), p); - FT d1 = squared_distance(c->vertex((li+1)&3)->point(), ab); - FT d2 = squared_distance(c->vertex((li+2)&3)->point(), ab); - FT d3 = squared_distance(c->vertex((li+3)&3)->point(), ab); - closest = va; - if (db < dist2) {dist2 = db; closest = vb;} - if (d1 < dist2) {dist2 = d1; closest = c->vertex((li+1)&3);} - if (d2 < dist2) {dist2 = d2; closest = c->vertex((li+2)&3);} - if (d3 < dist2) {dist2 = d3; closest = c->vertex((li+3)&3);} - } - default: - CGAL_triangulation_assertion(false); - } - - closest->steiner() = true; - return closest; -} - -template < class Gt, class Tds, class Itag > -typename Conforming_Delaunay_triangulation_3::Point -Conforming_Delaunay_triangulation_3:: -construct_split_point(Vertex_handle va, Vertex_handle vb, const Point& ref) { - // Blend of Si's and Shewchuk's segment recovery. - // We are actually interested in constructing a layer of protecting Steiner points around each acute vertex, - // such that no new Steiner point is inserted inside this layer. - // Shewchuk and Si place the protecting Steiner points on a ball around the acute vertex. - // However, determining if a vertex is acute is expensive. Additionally, the main purpose of the protecting points - // is to place new points outside the diametral ball with their acute vertex. - // This is always done if we place new Steiner points on the projection of the reference point onto the segment. - // To guarantee a minimum edge length, we never place a Steiner point closer to an endpoint than to the reference. - -#ifdef SPLIT_P_PROJ - // Place the split point at the midpoint, unless the reference encroaches upon the resulting edges. - // In that case, the split point is taken such that its projection onto ap or bp is p. - Point s = midpoint(va->point(), vb->point()); - - typename Geom_traits::Side_of_bounded_sphere_3 side_of_bounded_sphere = geom_traits().side_of_bounded_sphere_3_object(); - - if (side_of_bounded_sphere(va->point(), mid, ref) == ON_BOUNDED_SIDE) { - Object result = Gt().intersect_3_object()(Line(va->point(), vb->point()), Plane(ref, va->point() - ref)); - CGAL_triangulation_assertion(assign(s, result)); - } - else if (side_of_bounded_sphere(vb->point(), mid, ref) == ON_BOUNDED_SIDE) { - Object result = Gt().intersect_3_object()(Line(va->point(), vb->point()), Plane(ref, vb->point() - ref)); - CGAL_triangulation_assertion(assign(s, result)); - } - - return s; -#endif - - // Place the split point at the projection of the reference onto ab, unless that is too close to a or b, or it is not on the edge. - // In that case, the split point is the midpoint between a and b. - Line l(va->point(), vb->point()); - Point s = l.projection(ref); - - // The following will enforce a minimum edge length. - if (has_larger_distance_to_point(s, ref, va->point()) || has_larger_distance_to_point(s, ref, vb->point())) - s = midpoint(va->point(), vb->point()); - - // Note that for inexact construction kernels, this assertion may fail for both projections and midpoints. - //CGAL_triangulation_assertion(l.has_on(s)); - - // This method doesn't actually produce a correct bounding ball, because the Steiner points may be inserted closer to the endpoint than the reference point.. - - return s; - -#ifdef SPLIT_SI - // Si's point segment recovery. - typedef Gt::Sphere_3 Sphere; - Point s; - Object result; - if (c == Vertex_handle()) { - // Rule 1. - FT sd, r = squared_distance(a->point(), b->point()); - if (4*(sd = squared_distance(a->point(), p->point())) < r) { - result = Gt().intersect_3_object()(Sphere(a->point(), sd), Segment(a->point(), b->point())); - } - else if (4*(sd = squared_distance(b->point(), p->point())) < r) { - result = Gt().intersect_3_object()(Sphere(b->point(), sd), Segment(a->point(), b->point())); - } - else { - Vertex_handle v = insert(midpoint(a->point(), b->point()), a->cell()); - return v; - } - if (!assign(s, result)) return Vertex_handle(); - } - else { - // Rule 2. - FT r = squared_distance(c->point(), p->point()); - result = Gt().intersect_3_object()(Sphere(c->point(), r), Segment(a->point(), b->point())); - if (!assign(s, result)) return Vertex_handle(); - - FT sd = squared_distance(s, p->point()); - if (sd > squared_distance(s, b->point())) { - //Rule 3. - if (4*sd < squared_distance(s, a->point())) { - result = Gt().intersect_3_object()(Sphere(s, sd), Segment(a->point(), b->point())); - if (!assign(s, result)) - return Vertex_handle(); - } - else { - s = midpoint(a->point(), s); - } - } - } - - return s; -#endif -} - -template < class Gt, class Tds, class Itag > -template -void Conforming_Delaunay_triangulation_3:: -conform_segment(Vertex_handle va, Vertex_handle vb, Points_iterator encr_begin, Points_iterator encr_end) { - CGAL_triangulation_precondition(va != vb); - // This insertion method is an adjustment of the segment insertion - // in the paper "Constrained Delaunay tetrahedral mesh generation and refinement" by Si (2010). - // The points of the tetrahedra intersected by ab (encroaching points) are collected; - // from these points, the point making the largest circumsphere with ab is the reference point. - // The reference point indicates where ab is split and the subsegments each get a subset - // of the encroaching points and then the process is repeated. - - // Set the initial part of the segment that is already in the triangulation to conforming. - Vertex_handle vi; - Cell_handle c; - int li, lj; - while (includes_edge(va, vb, vi, c, li, lj) && (!make_gabriel || this->is_Gabriel(c, li, lj) ) ) { - set_conforming(c, li, lj, true); - va = vi; - if (va == vb) - return; - } - - // If the edge is not in the triangulation, there must be encroaching points. - // Collect the points encroaching upon and close to ab. - if (encr_begin == encr_end) { - Point_list encroaching; - bool intersect = collect_encroaching(va, vb, std::back_inserter(encroaching), vi); - - if (!intersect && encroaching.empty()) { - // There are occurrences where none of the vertices of the cells intersected by ab encroach upon ab. - // In this case, we insert the midpoint of ab. - vi = insert(midpoint(va->point(), vb->point()), va->cell()); - vi->steiner() = true; - intersect = true; - } - - // Insert ai, and if needed ib. - conform_segment(va, vi, encroaching.begin(), encroaching.end()); - if (vi != vb) { - conform_segment(vi, vb, encr_begin, encr_end); - } - return; - } - - // Find the reference point. - // This is the encroaching point closest to the segment. - FT sr, radius2 = 0; - Point ref; - for (Points_iterator it = encr_begin; it != encr_end; ++it) { - CGAL_triangulation_assertion(!collinear(va->point(), *it, vb->point())); - sr = squared_radius(va->point(), *it, vb->point()); - if (radius2 < sr) { - radius2 = sr; - ref = *it; - } - } - CGAL_triangulation_assertion(radius2 > 0); - - // Split the edge based on its reference point. - Point s = construct_split_point(va, vb, ref); - Vertex_handle vs = insert(s, va->cell()); - vs->steiner() = true; - - // Sort the encroaching points to encroaching upon as, sb, and only ab. - Points_iterator after_as = sort_encroaching(encr_begin, encr_end, va->point(), vs->point()); - Points_iterator after_sb = sort_encroaching(after_as, encr_end, vs->point(), vb->point()); - - // Insert the segments on both sides of s as constraints. - conform_segment(va, vs, encr_begin, after_as); - conform_segment(vs, vb, after_as, after_sb); -} - -template < class Gt, class Tds, class Itag > -void Conforming_Delaunay_triangulation_3:: -mark_segment(Vertex_handle va, Vertex_handle vb) { - Vertex_handle vi; - Cell_handle c; - int li, lj; - while (va != vb) { - // All segments must be in the triangulation. - CGAL_triangulation_assertion_code(bool ok = ) - includes_edge(va, vb, vi, c, li, lj); - CGAL_triangulation_assertion(ok); - - mark_edge(c, li, lj); - - va = vi; - } -} - -template < class Gt, class Tds, class Itag > -void Conforming_Delaunay_triangulation_3:: -insert_marked(Vertex_handle va, Vertex_handle vb) { - // Insert the segment and then mark all parts of the segment. - insert_conforming(va, vb); - mark_segment(va, vb); -} - -template < class Gt, class Tds, class Itag > -void Conforming_Delaunay_triangulation_3:: -set_conforming(Cell_handle c, int li, int lj, bool C) { - // The cell is already in the triangulation, so the edge is already Delaunay. - Vertex_handle vi = c->vertex(li), vj = c->vertex(lj); - switch (dimension()) { - case 3: { - // Constrain the edge in all its incident cells. - Cell_circulator it = incident_cells(c, li, lj), start(it); - do { - it->set_conforming(it->index(vi), it->index(vj), C); - ++it; - } while (it != start); - break; - } - case 2: { - // Constrain the edge in both cells that share it. - Cell_handle n = c->neighbor(3-li-lj); - n->set_conforming(n->index(vi), n->index(vj), C); - } - case 1: { - // Constrain the edge in the cell itself. - c->set_conforming(li, lj, C); - break; - } - } -} - -template < class Gt, class Tds, class Itag > -void Conforming_Delaunay_triangulation_3:: -mark_edge(Cell_handle c, int li, int lj) { - // The cell is already in the triangulation, so the edge is already Delaunay. - Vertex_handle vi = c->vertex(li), vj = c->vertex(lj); - switch (dimension()) { - case 3: { - // Constrain the edge in all its incident cells. - Cell_circulator it = incident_cells(c, li, lj), start(it); - do { - it->mark(it->index(vi), it->index(vj)); - it++; - } while (it != start); - break; - } - case 2: { - // Constrain the edge in both cells that share it. - Cell_handle n = c->neighbor(3-li-lj); - n->mark(n->index(vi), n->index(vj)); - } - case 1: { - // Constrain the edge in the cell itself. - c->mark(li, lj); - break; - } - } -} - -template < class Gt, class Tds, class Itag > -template < class InputIterator > -int Conforming_Delaunay_triangulation_3:: -insert_conforming_loop(InputIterator begin, InputIterator end, Cell_handle hint) { - // Insert the points in the loop. - std::list segments; - insert_loop_points(begin, end, std::back_inserter(segments), hint); - - // Insert the segments in the loop. - for (typename std::list::iterator it = segments.begin(); it != segments.end(); ++it) - insert_conforming(*it); - - return segments.size(); -} - -template < class Gt, class Tds, class Itag > -template < class InputIterator > -int Conforming_Delaunay_triangulation_3:: -insert_marked_loop(InputIterator begin, InputIterator end, Cell_handle hint) { - // Insert the points in the loop. - std::list segments; - insert_loop_points(begin, end, std::back_inserter(segments), hint); - - // Insert the segments in the loop. - for (typename std::list::iterator it = segments.begin(); it != segments.end(); ++it) - insert_marked(*it); - - return segments.size(); -} - -template < class Gt, class Tds, class Itag > -template < class InputIterator > -int Conforming_Delaunay_triangulation_3:: -remove_conforming_loop(InputIterator begin, InputIterator end, Cell_handle hint) { - // Insert the points in the loop. - std::list segments; - insert_loop_points(begin, end, std::back_inserter(segments), hint); - - // Insert the segments in the loop. - for (typename std::list::iterator it = segments.begin(); it != segments.end(); ++it) - remove_conforming(*it); - - return segments.size(); -} - -template < class Gt, class Tds, class Itag > -template < class InputIterator, class OutputIterator > -void Conforming_Delaunay_triangulation_3:: -insert_loop_points(InputIterator begin, InputIterator end, OutputIterator out, Cell_handle hint) { - if (begin == end) - return; - - // Insert the points in the loop. - Vertex_handle prev = insert(*begin, hint), first(prev); - hint = prev->cell(); - for (InputIterator it = ++begin; it != end; ++it) { - Vertex_handle cur = insert(*it, hint); - if (prev != cur) - *out++ = Bi_vertex(prev, cur); - hint = cur->cell(); - prev = cur; - } - if (first != prev) - *out++ = Bi_vertex(prev, first); -} - -template < class Gt, class Tds, class Itag > -void Conforming_Delaunay_triangulation_3:: -restore_conforming(Cell_handle c, int li) { - // Restore the conforming edges of facet i of c, based on its neighboring cell. - Cell_handle n = c->neighbor(li); - Vertex_handle v1, v2; Edge e; - for (int j = 0; j < dimension()+1; ++j) { - if (j != li) { - e = opposite_edge(c, li, (li+j)&3); - v1 = c->vertex(e.second); - v2 = c->vertex(e.third); - if (n->is_marked(n->index(v1), n->index(v2))) - e.first->mark(e.second, e.third); - else - e.first->set_conforming(e.second, e.third, n->is_conforming(n->index(v1), n->index(v2))); - } - } -} - -template < class Gt, class Tds, class Itag > -typename Conforming_Delaunay_triangulation_3::Vertex_handle -Conforming_Delaunay_triangulation_3:: -insert(const Point& p, Locate_type lt, Cell_handle c, int li, int lj) { - switch (dimension()) { - case 3: { - typename DT::Conflict_tester_3 tester(p, this); - Vertex_handle v = insert_in_conflict(p, lt, c, li, lj, tester, conforming_visitor); - return v; - } // dim 3 - case 2: { - typename DT::Conflict_tester_2 tester(p, this); - Vertex_handle v = insert_in_conflict(p, lt, c, li, lj, tester, conforming_visitor); - return v; - } // dim 2 - default : - // dimension <= 1 - // Do not use the generic insert. - return Tr::insert(p, c); - } -} - -template < class Gt, class Tds, class Itag > -void Conforming_Delaunay_triangulation_3:: -insert_conforming(const Point& a, const Point& b, Cell_handle hint) { - // Insert the two points and then conform the segment between them. - Vertex_handle va = insert(a, hint); - Vertex_handle vb = insert(b, va->cell()); - if (va != vb) - insert_conforming(va, vb); -} - -template < class Gt, class Tds, class Itag > -template < class InputIterator > -int Conforming_Delaunay_triangulation_3:: -insert_conforming(InputIterator begin, InputIterator end, Cell_handle hint) { - // Insert the points. - std::list segments; - for (InputIterator it = begin; it != end; ++it) { - Vertex_handle va = insert(it->first, hint); - hint = va->cell(); - if (it->first != it->second) { - Vertex_handle vb = insert(it->second, hint); - hint = vb->cell(); - segments.push_back(Bi_vertex(va, vb)); - } - } - - // Insert the segments - for (typename std::list::iterator it = segments.begin(); it != segments.end(); ++it) - insert_conforming(*it); - - return segments.size(); -} - -template < class Gt, class Tds, class Itag > -bool Conforming_Delaunay_triangulation_3:: -remove(Vertex_handle v) { - // The vertex cannot be removed if it is incident to a conforming edge. - if (are_there_incident_conforming(v)) - return false; - - // Collect the opposite facets as seen from the outside. - std::vector cells; cells.reserve(32); - incident_cells(v, std::back_inserter(cells)); - std::vector boundary; boundary.reserve(cells.size()); - - for (typename std::vector::const_iterator it = cells.begin(); it != cells.end(); ++it) { - Cell_handle n = (*it)->neighbor((*it)->index(v)); - boundary.push_back(Facet(n, n->index(*it))); - } - cells.clear(); - - // Remove the vertex. -#ifdef CGAL_DELAUNAY_3_OLD_REMOVE - if (dimension() == 3 && !test_dim_down(v)) { - remove_3D_ear(v); - } else { -#endif - cDT tmp; - typename DT::Vertex_remover remover(tmp); - Tr::remove(v, remover); -#ifdef CGAL_DELAUNAY_3_OLD_REMOVE - } -#endif - - // Reinsert any conforming edges on the boundary of the retriangulated region. - for (typename std::vector::iterator it = boundary.begin(); it != boundary.end(); ++it) { - restore_conforming(mirror_facet(*it)); - } - - CGAL_triangulation_expensive_postcondition(is_valid()); - return true; -} - -template < class Gt, class Tds, class Itag > -bool Conforming_Delaunay_triangulation_3:: -remove_conforming(Cell_handle c, int li, int lj) { - // Mark the edge as non-conforming - set_conforming(c, li, lj, false); - return true; -} - -template < class Gt, class Tds, class Itag > -bool Conforming_Delaunay_triangulation_3:: -remove_conforming(Vertex_handle va, Vertex_handle vb) { - CGAL_triangulation_precondition(va != vb); - // All the parts of the segment ab that are constrained and - // not incident to a constrained facet are marked unconstrained. - // The return is true if the whole segment ab is covered by unconstrained edges; - // otherwise, returns false; - - // Check if first part of the segment is an edge in the triangulation. - Vertex_handle vi; - Cell_handle c; - int li, lj; - bool complete = true; - if (includes_edge(va, vb, vi, c, li, lj)) { - // Remove the first part. - if (!remove_conforming(c, li, lj)) - complete = false; - if (vi != vb) { - // Remove the next part. - if (!remove_conforming(vi, vb)) - complete = false; - } - return complete; - } - - return false; -} - -template < class Gt, class Tds, class Itag > -bool Conforming_Delaunay_triangulation_3:: -are_there_incident_conforming(Vertex_handle v) const { - if (dimension() > 0) { - std::vector edges; - edges.reserve(32); - finite_incident_edges(v, std::back_inserter(edges)); - for (typename std::vector::iterator it = edges.begin(); it != edges.end(); ++it) - if (is_conforming(*it)) - return true; - } - - return false; -} - -template < class Gt, class Tds, class Itag > -bool Conforming_Delaunay_triangulation_3:: -includes_edge(Vertex_handle va, Vertex_handle vb, Vertex_handle& vc, Cell_handle& c, int& li, int& lj) const { - CGAL_triangulation_precondition(!is_infinite(va) && !is_infinite(vb)); - - // Returns true if the line segment ab contains an edge e - // incident to a, false otherwise. - // If true, vc becomes the vertex of e distinct from a, - // c is a cell incident to e where e = (c, li, lj). - Vertex_handle v; - std::vector edges; edges.reserve(64); - finite_incident_edges(va, std::back_inserter(edges)); - for (typename std::vector::const_iterator it = edges.begin(); it != edges.end(); ++it) { - // Find the other vertex of the edge. - v = it->first->vertex(it->second); - if (v == va) v = it->first->vertex(it->third); - - if (is_infinite(v)) - continue; - - // Check if this edge is (part of) ab. - if (v == vb) { - vc = vb; - c = it->first; - li = it->second; - lj = it->third; - return true; - } - else if (sign((va->point() - v->point()) * (vb->point() - v->point())) == NEGATIVE && - collinear(va->point(), v->point(), vb->point())/* && - collinear_are_ordered_along_line(va->point(), v->point(), vb->point())*/) { - vc = v; - c = it->first; - li = it->second; - lj = it->third; - return true; - } - } - return false; -} - -template < class Gt, class Tds, class Itag > -bool Conforming_Delaunay_triangulation_3:: -is_valid(bool verbose, int level) const { - if (!tds().is_valid(verbose,level)) { - if (verbose) - std::cerr << "invalid data structure" << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - - if (infinite_vertex() == Vertex_handle()) { - if (verbose) - std::cerr << "no infinite vertex" << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - - int inf; - Cell_handle n; - Vertex_handle v1, v2; - switch (dimension()) { - case 3: { - for (All_cells_iterator it = this->all_cells_begin(); it != this->all_cells_end(); ++it) { - for (int i = 0; i < 4; ++i) { - if (!it->neighbor(i)->has_neighbor(it)) { - if (verbose) - std::cerr << "inconsistent neighbors " << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - } - if (it->has_vertex(infinite_vertex(), inf)) { - for (int i = 0; i < 4; ++i) { - if (i != inf) { - if (is_conforming(it, i, inf)) { - if (verbose) - std::cerr << "conforming infinite edge " << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - } - } - } - else { - is_valid_finite(it); - for (int i = 0; i < 4; ++i) { - n = it->neighbor(i); - if (!is_infinite(n->vertex(n->index(it)))) { - if (this->side_of_sphere(it, n->vertex(n->index(it))->point()) == ON_BOUNDED_SIDE) { - if (verbose) - std::cerr << "non-empty sphere " << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - } - } - } - } - for (Finite_edges_iterator it = finite_edges_begin(); it != finite_edges_end(); ++it) { - v1 = it->first->vertex(it->second); - v2 = it->first->vertex(it->third); - bool conforming = is_conforming(*it); - - Cell_circulator cit = incident_cells(it->first, it->second, it->third, it->first), start(cit); - do { - if (cit->is_marked(cit->index(v1), cit->index(v2))) { - if (verbose) - std::cerr << "marked edge " << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - - if (cit->is_conforming(cit->index(v1), cit->index(v2)) != conforming) { - if (verbose) - std::cerr << "inconsistent conforming edge " << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - cit++; - } while (cit != start); - } - break; - } - case 2: { - for (typename DT::All_facets_iterator it = this->all_facets_begin(); it != this->all_facets_end(); ++it) { - for (int i = 0; i < 3; ++i) { - if (!it->first->neighbor(i)->has_neighbor(it->first)) { - if (verbose) - std::cerr << "inconsistent neighbors " << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - } - if (it->first->has_vertex(infinite_vertex(), inf) && inf < 3) { - for (int i = 0; i < 3; ++i) { - if (i != inf) { - if (is_conforming(it->first, i, 3)) { - if (verbose) - std::cerr << "conforming infinite edge " << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - } - } - } - else { - is_valid_finite(it->first); - for (int i = 0; i < 3; ++i) { - n = it->first->neighbor(i); - if(!is_infinite(n->vertex(n->index(it->first)))) { - if (this->side_of_circle(it->first, 3, n->vertex(n->index(it->first))->point()) == ON_BOUNDED_SIDE) { - if (verbose) - std::cerr << "non-empty circle " << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - } - - Edge o = this->opposite_edge(it->first, i, 3); - if (o.first->is_marked(o.second, o.third)) { - if (verbose) - std::cerr << "marked edge " << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - - if (is_conforming(o) != - is_conforming(n, n->index(o.first->vertex(o.second)), n->index(o.first->vertex(o.third)))) { - if (verbose) - std::cerr << "inconsistent conforming edge " << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - } - } - } - break; - } - case 1: { - for (Finite_edges_iterator it = finite_edges_begin(); it != finite_edges_end(); ++it) - is_valid_finite(it->first); - break; - } - } - if (verbose) - std::cerr << "valid conforming Delaunay triangulation" << std::endl; - return true; -} - -template < class Gt, class Tds, class Itag > -bool Conforming_Delaunay_triangulation_3:: -is_valid(Cell_handle c, bool verbose, int level) const { - if (!DT::is_valid(c, verbose, level)) { - CGAL_triangulation_assertion(false); - return false; - } - - switch (dimension()) { - case 3: { - for (int i = 0; i < 4; ++i) { - for (int j = i+1; j < 4; ++j) { - Vertex_handle v1 = c->vertex(i), v2 = c->vertex(j); - if (is_infinite(v1) || is_infinite(v2)) { - if (is_conforming(c, i, j)) { - if (verbose) - std::cerr << "conforming infinite edge " << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - } - else { - bool conforming = is_conforming(c, i, j); - - Cell_circulator cit = incident_cells(c, i, j, c), start(cit); - do { - if (cit->is_marked(cit->index(v1), cit->index(v2))) { - if (verbose) - std::cerr << "marked edge " << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - - if (cit->is_conforming(cit->index(v1), cit->index(v2)) != conforming) { - if (verbose) - std::cerr << "inconsistent conforming edge " << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - cit++; - } while (cit != start); - } - } - } - break; - } - case 2: { - for (int v=0; v<3; ++v) { - int i = (v+1)%3, - j = (v+2)%3; - - Vertex_handle v1 = c->vertex(i), v2 = c->vertex(j); - if (is_infinite(v1) || is_infinite(v2)) { - if (is_conforming(c, i, j)) { - if (verbose) - std::cerr << "conforming infinite edge " << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - } - else { - if (c->is_marked(i, j)) { - if (verbose) - std::cerr << "marked edge " << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - - Cell_handle n = c->neighbor(v); - if (c->is_conforming(i, j) != n->is_conforming(n->index(v1), n->index(v2))) { - if (verbose) - std::cerr << "inconsistent conforming edge " << std::endl; - CGAL_triangulation_assertion(false); - return false; - } - } - } - break; - } - } - if (verbose) - std::cerr << "valid Delaunay cell" << std::endl; - return true; -} - -template < class Gt, class Tds, class Itag > -template -bool Conforming_Delaunay_triangulation_3:: -collect_encroaching(Vertex_handle va, Vertex_handle vb, OutputIterator encroaching, Vertex_handle& vc) { - CGAL_triangulation_precondition(va != vb); - CGAL_triangulation_precondition(!is_infinite(va) && !is_infinite(vb)); - CGAL_triangulation_precondition(dimension() >= 2); - - // Fill encroaching with the points that are encroaching upon the segment ab. - // vc is set to the first vertex on ab starting from a; - // if ab intersects an edge or facet the intersection is inserted and true returned. - - // This is done by traversing cells from a towards b until either b is found, - // or any other collinear vertex, or an intersecting constraint. - Locate_type lt; int li, lj; - Encroaching_collecter cit(va, vb->point(), this, encroaching); - while (cit.has_next()) { - // Walk to the next cell. - ++cit; - - cit.traversed(lt, li, lj); - - // Stop when either a new vertex is reached, or a barrier (a conforming edge) is encountered. - if (lt == DT::VERTEX) { - vc = cit->vertex(li); - if (va == vc) - continue; - return false; - } - - if (cit.barrier_hit()) { - vc = intersect(cit, lt, li, lj, va, vb); - return true; - } - } - - // The cell must contain the target. - CGAL_triangulation_assertion(cit->has_vertex(vb)); - vc = vb; - return false; -} - -template < class Gt, class Tds, class Itag > -typename Conforming_Delaunay_triangulation_3::Points_iterator -Conforming_Delaunay_triangulation_3:: -sort_encroaching(Points_iterator begin, Points_iterator end, const Point& a, const Point& b) const { - // Sorts the points encroaching upon ab to the front of the collection. - // Returns an iterator to the first point not encroaching upon ab. - Point temp; - Points_iterator middle = end; - for (Points_iterator it = begin; it != middle;) { - if (collinear(a, b, *it) || !is_encroaching(a, b, *it)) { - if (it != --middle) { - // Swap the current point to the back. - temp = *middle; - *middle = *it; - *it = temp; - } - } - else - ++it; - } - return middle; -} - -} //end of CGAL namespace - -#endif // CGAL_CONFORMING_DELAUNAY_TRIANGULATION_3_H diff --git a/Conforming_triangulation_3/include/CGAL/Conforming_triangulation_cell_base_3.h b/Conforming_triangulation_3/include/CGAL/Conforming_triangulation_cell_base_3.h deleted file mode 100644 index 936afbd0fe7..00000000000 --- a/Conforming_triangulation_3/include/CGAL/Conforming_triangulation_cell_base_3.h +++ /dev/null @@ -1,179 +0,0 @@ -//A tetrahedral cell in a conforming Delaunay triangulation. -//Copyright (C) 2012 Utrecht University -// -//This program is free software: you can redistribute it and/or modify -//it under the terms of the GNU General Public License as published by -//the Free Software Foundation, either version 3 of the License, or -//(at your option) any later version. -// -//This program is distributed in the hope that it will be useful, -//but WITHOUT ANY WARRANTY; without even the implied warranty of -//MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -//GNU General Public License for more details. -// -//You should have received a copy of the GNU General Public License -//along with this program. If not, see . -// -// Author(s): Thijs van Lankveld - -// Based on CGAL/Triangulation_cell_base_3.h -// - -// Note that a constrained 3D Delaunay triangulation is partially conforming and guided. -// 1-dimensional constraints are specified by segments. -// 2-demensional constraints are specified by closure complete polygons. -// The constrained 3D DT must contain a collection of edges and faces that exactly cover -// each constrained segment or polygon. For example, a constrained edge will not necessarily -// be maintained as vertices are inserted, but after updating the triangulation to the vertex -// insertion, the constrained edge must be the union of one or more constrained edges. - -#ifndef CGAL_CONFORMING_TRIANGULATION_CELL_BASE_3_H -#define CGAL_CONFORMING_TRIANGULATION_CELL_BASE_3_H - -#include -#include -#include - -namespace CGAL { - -template < typename Gt, typename Cb = Triangulation_cell_base_3 > -class Conforming_triangulation_cell_base_3: public Cb { -public: - typedef typename Cb::Vertex_handle Vertex_handle; - typedef typename Cb::Cell_handle Cell_handle; - - typedef Gt Geom_traits; - typedef typename Geom_traits::Point_3 Point; - - typedef Point* Point_container; - typedef Point* Point_iterator; - typedef const Point* Point_const_iterator; - - template < typename TDS2 > - struct Rebind_TDS { - typedef typename Cb::template Rebind_TDS::Other Cb2; - typedef Conforming_triangulation_cell_base_3 Other; - }; // Rebind_TDS - -protected: - // The edge states. - bool EF[6]; - bool mm[6]; - -public: - Conforming_triangulation_cell_base_3(): Cb() {clear_conforming();} - Conforming_triangulation_cell_base_3(const Vertex_handle& v0, const Vertex_handle& v1, - const Vertex_handle& v2, const Vertex_handle& v3) - : Cb(v0, v1, v2, v3) {clear_conforming();} - Conforming_triangulation_cell_base_3(const Vertex_handle& v0, const Vertex_handle& v1, - const Vertex_handle& v2, const Vertex_handle& v3, - const Cell_handle& n0, const Cell_handle& n1, - const Cell_handle& n2, const Cell_handle& n3) - : Cb(v0, v1, v2, v3, n0, n1, n2, n3) {clear_conforming();} - - bool is_conforming(int i, int j) const {return EF[getEdgeIndex(i,j)];} - bool is_marked(int i, int j) const {return mm[getEdgeIndex(i,j)];} - - void set_conforming(bool c0, bool c1, bool c2, bool c3, bool c4, bool c5) { - EF[0] = c0; - EF[1] = c1; - EF[2] = c2; - EF[3] = c3; - EF[4] = c4; - EF[5] = c5; - clear_marked(); - } - void set_conforming(int i, int j, bool c) { - int index = getEdgeIndex(i,j); - EF[index] = c; - mm[index] = false; - } - void mark(int i, int j) { - int index = getEdgeIndex(i,j); - EF[index] = true; - mm[index] = true; - } - - void clear_conforming() {set_conforming(false, false, false, false, false, false);} - void clear_marked() { - mm[0] = false; - mm[1] = false; - mm[2] = false; - mm[3] = false; - mm[4] = false; - mm[5] = false; - } - virtual void clear() {clear_conforming();} - - bool has_conforming() const {return EF[0] || EF[1] || EF[2] || EF[3] || EF[4] || EF[5];} - - virtual void reorient() { - bool tmp = EF[1]; - EF[1] = EF[3]; - EF[3] = tmp; - tmp = mm[1]; - mm[1] = mm[3]; - mm[3] = tmp; - } - - virtual std::istream& read_cell(std::istream& is); - virtual std::ostream& write_cell(std::ostream& os) const; - -private: - int getEdgeIndex(int i, int j) const { - CGAL_triangulation_precondition(i>=0 && i<=3); - CGAL_triangulation_precondition(j>=0 && j<=3); - CGAL_triangulation_precondition(i != j); - return (i==0 || j==0) ? i+j-1 : i+j; - } -}; // Conforming_triangulation_cell_base_3 - -template < class Gt, class Cb > -inline std::istream& operator>>(std::istream& is, Conforming_triangulation_cell_base_3& c) { - is >> static_cast(c); - return c.read_cell(is); -} - -template < class Gt, class Cb > -inline std::ostream& operator<<(std::ostream &os, const Conforming_triangulation_cell_base_3& c) { - os << static_cast(c); - return c.write_cell(os); -} - -template < class Gt, class Cb > -std::istream& Conforming_triangulation_cell_base_3::read_cell(std::istream& is) { - char s; - for (int i = 0; i < 4; ++i) { - for (int j = i+1; j < 4; ++j) { - if (is_ascii(is)) - is >> s; - else - read(is, s); - if (s == 'C') - set_conforming(i, j, true); - } - } - - return is; -} - -template < class Gt, class Cb > -std::ostream& Conforming_triangulation_cell_base_3::write_cell(std::ostream& os) const { - for (int i = 0; i < 4; ++i) { - for (int j = i+1; j < 4; ++j) { - if (is_conforming(i, j)) { - os << "C"; - } - else - os << "N"; - if (is_ascii(os)) - os << ' '; - } - } - - return os; -} - -} //end of CGAL namespace - -#endif // CGAL_CONFORMING_TRIANGULATION_CELL_BASE_3_H \ No newline at end of file diff --git a/Conforming_triangulation_3/include/CGAL/Conforming_triangulation_vertex_base_3.h b/Conforming_triangulation_3/include/CGAL/Conforming_triangulation_vertex_base_3.h deleted file mode 100644 index afd635b24c9..00000000000 --- a/Conforming_triangulation_3/include/CGAL/Conforming_triangulation_vertex_base_3.h +++ /dev/null @@ -1,69 +0,0 @@ -//A vertex in the conforming Delaunay triangulation structure. -//Copyright (C) 2012 Utrecht University -// -//This program is free software: you can redistribute it and/or modify -//it under the terms of the GNU General Public License as published by -//the Free Software Foundation, either version 3 of the License, or -//(at your option) any later version. -// -//This program is distributed in the hope that it will be useful, -//but WITHOUT ANY WARRANTY; without even the implied warranty of -//MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -//GNU General Public License for more details. -// -//You should have received a copy of the GNU General Public License -//along with this program. If not, see . -// -// Author(s): Thijs van Lankveld - -#ifndef CGAL_CONFORMING_TRIANGULATION_VERTEX_BASE_3_H -#define CGAL_CONFORMING_TRIANGULATION_VERTEX_BASE_3_H - -#include - -namespace CGAL { - -template < typename GT, typename Vb = Triangulation_vertex_base_3 > -class Conforming_triangulation_vertex_base_3: public Vb { - bool _s; - -public: - typedef typename Vb::Cell_handle Cell_handle; - typedef typename Vb::Point Point; - - template < typename TDS2 > struct Rebind_TDS { - typedef typename Vb::template Rebind_TDS::Other Vb2; - typedef Conforming_triangulation_vertex_base_3 Other; - }; - - Conforming_triangulation_vertex_base_3(): Vb(), _s(false) {} - Conforming_triangulation_vertex_base_3(const Point& p): Vb(p), _s(false) {} - Conforming_triangulation_vertex_base_3(const Point& p, bool s): Vb(p), _s(s) {} - Conforming_triangulation_vertex_base_3(const Point& p, Cell_handle c): Vb(p, c), _s(false) {} - Conforming_triangulation_vertex_base_3(const Point& p, bool s, Cell_handle c): Vb(p, c), _s(s) {} - Conforming_triangulation_vertex_base_3(Cell_handle c): Vb(c), _s(false) {} - - const bool& steiner() const {return _s;} - bool& steiner() {return _s;} -}; - -template < class GT, class Vb > -std::istream& -operator>>(std::istream &is, Conforming_triangulation_vertex_base_3 &v) { - is >> static_cast(v); - if (is_ascii(is)) is >> v.steiner(); - else read(is, v.steiner()); - return is; -} - -template < class GT, class Vb > -std::ostream& -operator<<(std::ostream &os, const Conforming_triangulation_vertex_base_3 &v) { - os << static_cast(v); - if (is_ascii(os)) os << ' '; - return os << v.steiner(); -} - -} //end of CGAL namespace - -#endif // CGAL_CONFORMING_TRIANGULATION_VERTEX_BASE_3_H diff --git a/Conforming_triangulation_3/include/CGAL/Delaunay_triangulation_utils_3.h b/Conforming_triangulation_3/include/CGAL/Delaunay_triangulation_utils_3.h deleted file mode 100644 index 59874fc1fc1..00000000000 --- a/Conforming_triangulation_3/include/CGAL/Delaunay_triangulation_utils_3.h +++ /dev/null @@ -1,406 +0,0 @@ -//Some additional utilities for the Delaunay triangulation structure. -//Copyright (C) 2012 Utrecht University -// -//This program is free software: you can redistribute it and/or modify -//it under the terms of the GNU General Public License as published by -//the Free Software Foundation, either version 3 of the License, or -//(at your option) any later version. -// -//This program is distributed in the hope that it will be useful, -//but WITHOUT ANY WARRANTY; without even the implied warranty of -//MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -//GNU General Public License for more details. -// -//You should have received a copy of the GNU General Public License -//along with this program. If not, see . -// -// Author(s): Thijs van Lankveld - -// Based on CGAL/Delaunay_triangulation_3.h -// -// For the love of something, why is the 3D triangulation package so badly written? - -#ifndef CGAL_DELAUNAY_TRIANGULATION_3_UTILS_H -#define CGAL_DELAUNAY_TRIANGULATION_3_UTILS_H - -#include - -#include -#include - -#include -#include - -namespace CGAL { - -template < class DT > class Natural_neighbors_3; - -template < class Gt, class Tds > class Triangulation_cell_traverser_3; - -template < class Gt, - class Tds_ = Triangulation_data_structure_3 < Triangulation_vertex_base_3, Triangulation_cell_base_3 > > -class Delaunay_triangulation_utils_3: public Delaunay_triangulation_3 { - //typedef Triangulation_data_structure Tds; - typedef Tds_ Tds; - typedef Delaunay_triangulation_utils_3 Self; - typedef Delaunay_triangulation_3 DT; - typedef Triangulation_3 Tr; - - friend class Natural_neighbors_3; - -public: - typedef Tds Triangulation_data_structure; - typedef Gt Geom_traits; - - typedef typename Gt::FT FT; - - typedef typename Gt::Point_3 Point; - typedef typename Gt::Vector_3 Vector; - typedef typename Gt::Segment_3 Segment; - typedef typename Gt::Line_3 Line; - typedef typename Gt::Triangle_3 Triangle; - typedef typename Gt::Tetrahedron_3 Tetrahedron; - typedef typename Gt::Plane_3 Plane; - - typedef typename DT::Cell_handle Cell_handle; - typedef typename DT::Vertex_handle Vertex_handle; - - typedef typename DT::Cell Cell; - typedef typename DT::Vertex Vertex; - typedef typename DT::Facet Facet; - typedef typename DT::Edge Edge; - - typedef std::pair Bi_vertex; - typedef Triple Tri_vertex; - - typedef std::pair Constraint_1; - typedef std::vector Constraint_2; - - typedef typename DT::Cell_circulator Cell_circulator; - typedef typename DT::Facet_circulator Facet_circulator; - typedef typename DT::Cell_iterator Cell_iterator; - typedef typename DT::Facet_iterator Facet_iterator; - typedef typename DT::Edge_iterator Edge_iterator; - typedef typename DT::Vertex_iterator Vertex_iterator; - - typedef typename DT::Finite_vertices_iterator Finite_vertices_iterator; - typedef typename DT::Finite_cells_iterator Finite_cells_iterator; - typedef typename DT::Finite_facets_iterator Finite_facets_iterator; - typedef typename DT::Finite_edges_iterator Finite_edges_iterator; - - typedef typename DT::All_cells_iterator All_cells_iterator; - - typedef Triangulation_segment_traverser_3 Cell_traverser; - - typedef typename DT::Locate_type Locate_type; - -#ifndef CGAL_CFG_USING_BASE_MEMBER_BUG_2 - using DT::cw; - using DT::ccw; - using DT::coplanar; - using DT::coplanar_orientation; - using DT::dimension; - using DT::finite_facets_begin; - using DT::finite_facets_end; - using DT::finite_vertices_begin; - using DT::finite_vertices_end; - using DT::finite_cells_begin; - using DT::finite_cells_end; - using DT::finite_edges_begin; - using DT::finite_edges_end; - using DT::geom_traits; - using DT::infinite_vertex; - using DT::is_valid; - using DT::mirror_vertex; - using DT::next_around_edge; - using DT::number_of_vertices; - using DT::orientation; - using DT::remove; - using DT::coplanar_side_of_bounded_circle; - using DT::side_of_oriented_sphere; - using DT::side_of_segment; - using DT::tds; - using DT::vertex_triple_index; - using DT::locate; - using DT::insert_in_conflict; - using DT::is_infinite; - using DT::insert; -#endif - -protected: - typedef std::list Tri_vertex_collection; - typedef std::list Bi_vertex_collection; - -protected: - // Test whether a newly inserted point conflicts with the existing cells. - class Conflict_tester_3 { - protected: - const Point& _p; - const Self* _tr; - - public: - Conflict_tester_3(const Point& pt, const Self* tr): _p(pt), _tr(tr) {} - bool operator()(const Cell_handle c) const {return _tr->side_of_sphere(c, _p, true) == ON_BOUNDED_SIDE;} - Oriented_side compare_weight(const Point&, const Point&) const {return ZERO;} - bool test_initial_cell(Cell_handle) const {return true;} - }; // class Conflict_tester_3 - - // Test whether a newly inserted point conflicts with the existing cells. - class Conflict_tester_2 { - protected: - const Point& _p; - const Self* _tr; - - public: - Conflict_tester_2(const Point& pt, const Self* tr): _p(pt), _tr(tr) {} - bool operator()(const Cell_handle c) const {return _tr->side_of_circle(c, 3, _p, true) == ON_BOUNDED_SIDE;} - Oriented_side compare_weight(const Point&, const Point&) const {return ZERO;} - bool test_initial_cell(Cell_handle) const {return true;} - }; // class Conflict_tester_2 - - class Hidden_point_visitor { - public: - Hidden_point_visitor() {} - - template < class InputIterator > - void process_cells_in_conflict(InputIterator, InputIterator) const {} - void reinsert_vertices(Vertex_handle) {} - Vertex_handle replace_vertex(Cell_handle c, int index, const Point&) {return c->vertex(index);} - void hide_point(Cell_handle, const Point&) {} - }; // class Hidden_point_visitor - - Hidden_point_visitor hidden_point_visitor; - -public: - // Constructors for bi-vertices. - inline Bi_vertex to_bi_vertex(Cell_handle c, int li, int lj) const {return Bi_vertex(c->vertex(li), c->vertex(lj));} - inline Bi_vertex to_bi_vertex(const Edge& e) const {return to_bi_vertex(e.first, e.second, e.third);} - inline Bi_vertex sort_bi_vertex(Vertex_handle v1, Vertex_handle v2) const {if (v1 < v2) return Bi_vertex(v1, v2); return Bi_vertex(v2, v1);} - inline Bi_vertex sort_bi_vertex(const Bi_vertex& bv) const {return sort_bi_vertex(bv.first, bv.second);} - inline Bi_vertex sort_bi_vertex(Cell_handle c, int li, int lj) const {return sort_bi_vertex(c->vertex(li), c->vertex(lj));} - inline Bi_vertex sort_bi_vertex(const Edge& e) const {return sort_bi_vertex(e.first, e.second, e.third);} - - // Constructors for tri-vertices. - inline Tri_vertex to_tri_vertex(Cell_handle c, int li) const {return Tri_vertex(c->vertex((li+1)&3), c->vertex((li+2)&3), c->vertex((li+3)&3));} - inline Tri_vertex to_tri_vertex(const Facet& f) const {return to_tri_vertex(f.first, f.second);} - inline Tri_vertex sort_tri_vertex(Vertex_handle v1, Vertex_handle v2, Vertex_handle v3) const { - CGAL_triangulation_precondition(v1 != v2 && v1 != v3 && v2 != v3); - if (v1 < v2) { - if (v2 < v3) - return Tri_vertex(v3, v2, v1); - if (v3 < v1) - return Tri_vertex(v2, v1, v3); - return Tri_vertex(v2, v3, v1); - } - else { // v2 < v1 - if (v3 < v2) - return Tri_vertex(v1, v2, v3); - if (v1 < v3) - return Tri_vertex(v3, v1, v2); - return Tri_vertex(v1, v3, v2); - } - } - inline Tri_vertex sort_tri_vertex(const Tri_vertex& tv) const {return sort_tri_vertex(tv.first, tv.second, tv.third);} - inline Tri_vertex sort_tri_vertex(Cell_handle c, int li) const {return sort_tri_vertex(to_tri_vertex(c, li));} - inline Tri_vertex sort_tri_vertex(const Facet& f) const {return sort_tri_vertex(f.first, f.second);} - inline Tri_vertex oriented_tri_vertex(Cell_handle c, int li) const { - CGAL_triangulation_precondition(dimension() == 2 || dimension() == 3); - CGAL_triangulation_precondition((dimension() == 2 && li == 3) || (dimension() == 3 && li >= 0 && li <= 3) ); - CGAL_triangulation_precondition(!is_infinite(c, li)); - if ((li&1) == 0) - return Tri_vertex(c->vertex((li+2)&3), - c->vertex((li+1)&3), - c->vertex((li+3)&3)); - return Tri_vertex(c->vertex((li+1)&3), - c->vertex((li+2)&3), - c->vertex((li+3)&3)); - } - inline Tri_vertex oriented_tri_vertex(const Facet& f) const {return oriented_tri_vertex(f.first, f.second);} - -protected: - // tri_vertex methods similar to Triangulation_3. - Tri_vertex make_tri_vertex(const Facet& f) const { - return Tri_vertex(f.first->vertex(vertex_triple_index(f.second, 0)), - f.first->vertex(vertex_triple_index(f.second, 1)), - f.first->vertex(vertex_triple_index(f.second, 2))); - } - void make_canonical(Tri_vertex& t) const { - int i = (&*(t.first) < &*(t.second))? 0 : 1; - if (i == 0) - i = (&*(t.first) < &*(t.third))? 0 : 2; - else - i = (&*(t.second) < &*(t.third))? 1 : 2; - Vertex_handle tmp; - switch (i) { - case 0: return; - case 1: - tmp = t.first; - t.first = t.second; - t.second = t.third; - t.third = tmp; - return; - default: - tmp = t.first; - t.first = t.third; - t.third = t.second; - t.second = tmp; - } - } - - // PREDICATES - Orientation coplanar_orientation(const Point& p0, const Point& p1, const Point& p2, const Point& p) const { - return Gt().coplanar_orientation_3_object()(p0, p1, p2, p); - } - -public: - // CONSTRUCTORS - Delaunay_triangulation_utils_3(const Gt& gt = Gt()): DT(gt) {} - - // Copy constructor duplicates vertices and cells. - Delaunay_triangulation_utils_3(const Delaunay_triangulation_utils_3& tr): DT(tr) { - CGAL_triangulation_postcondition(is_valid()); - } - - // Create a 3D Delaunay triangulation from a number of points. - template < typename InputIterator > - Delaunay_triangulation_utils_3(InputIterator first, InputIterator last, const Gt& gt = Gt()): DT(gt) { - insert(first, last); - CGAL_triangulation_postcondition(is_valid()); - } - -public: - // INSERT POINT - - Vertex_handle insert(const Point& p, Cell_handle start = Cell_handle()); - virtual Vertex_handle insert(const Point& p, Locate_type lt, Cell_handle c, int li, int lj); - - template < class InputIterator > - int insert(InputIterator first, InputIterator last) { - int n = number_of_vertices(); - - std::vector points(first, last); - std::random_shuffle(points.begin(), points.end()); - spatial_sort(points.begin(), points.end(), geom_traits()); - - Cell_handle hint; - for (typename std::vector::const_iterator it = points.begin(), end = points.end(); it != end; ++it) - hint = insert(*it, hint)->cell(); - - return number_of_vertices() - n; - } - -public: - // MOVE - - Vertex_handle move_point(Vertex_handle v, const Point& p); - -public: - // Construct the edge of c opposite to ij. - Edge opposite_edge(Cell_handle c, int li, int lj) const { - CGAL_triangulation_precondition(li >= 0 && li < 4); - CGAL_triangulation_precondition(lj >= 0 && lj < 4); - CGAL_triangulation_precondition(li != lj); - - switch (6-li-lj) { // i + j + missing indices = 6. - case 1: return Edge(c, 0, 1); - case 2: return Edge(c, 0, 2); - case 3: return (li == 0 || lj == 0) - ? Edge(c, 1, 2) - : Edge(c, 0, 3); - case 4: return Edge(c, 1, 3); - case 5: return Edge(c, 2, 3); - } - - CGAL_triangulation_assertion(false); - return Edge(); - } - Edge opposite_edge(const Edge& e) const {return opposite_edge(e.first, e.second, e.third);} - - // Give the same facet as seen from the other side. - Facet mirror_facet(Cell_handle c, int li) const {return Facet(c->neighbor(li), c->neighbor(li)->index(c));} - Facet mirror_facet(const Facet& f) const {return mirror_facet(f.first, f.second);} - - // Construct a plane of a facet. - Plane plane(Cell_handle c, int li) const { - /* This should be removed to reduce orientation errors with the inexact kernel*/ - CGAL_triangulation_precondition(dimension() >= 2); - CGAL_triangulation_precondition(dimension() == 3 || li == 3); - CGAL_triangulation_precondition(li >= 0 && li <= 3); - CGAL_triangulation_precondition(!is_infinite(c, li)); - if ((li&1) == 0) - return Plane(c->vertex((li+2)&3)->point(), - c->vertex((li+1)&3)->point(), - c->vertex((li+3)&3)->point()); - return Plane(c->vertex((li+1)&3)->point(), - c->vertex((li+2)&3)->point(), - c->vertex((li+3)&3)->point()); - } - Plane plane(const Facet& f) const {return plane(f.first, f.second);} - - // Construct a line of an edge. - Line line(Cell_handle c, int li, int lj) const { - CGAL_triangulation_precondition(dimension() >= 1); - CGAL_triangulation_precondition(li >= 0 && li <= 3); - CGAL_triangulation_precondition(lj >= 0 && lj <= 3); - CGAL_triangulation_precondition(li != lj); - CGAL_triangulation_precondition(li + lj < dimension() * 2); - CGAL_triangulation_precondition(!is_infinite(c, li, lj)); - return Line(c->vertex(li)->point(), - c->vertex(lj)->point()); - } - Line line(const Edge& e) const {return line(e.first, e.second, e.third);} -}; // Delaunay_triangulation_utils_3 - -template < class Gt, class Tds > -typename Delaunay_triangulation_utils_3::Vertex_handle -Delaunay_triangulation_utils_3:: -insert(const Point& p, Cell_handle start) { - Locate_type lt; - int li, lj; - Cell_handle c = locate(p, lt, li, lj, start); - return insert(p, lt, c, li, lj); -} - -template < class Gt, class Tds > -typename Delaunay_triangulation_utils_3::Vertex_handle -Delaunay_triangulation_utils_3:: -insert(const Point& p, Locate_type lt, Cell_handle c, int li, int lj) { - switch (dimension()) { - case 3: { - Conflict_tester_3 tester(p, this); - Vertex_handle v = insert_in_conflict(p, lt, c, li, lj, tester, hidden_point_visitor); - return v; - }// dim 3 - case 2: { - Conflict_tester_2 tester(p, this); - return insert_in_conflict(p, lt, c, li, lj, tester, hidden_point_visitor); - }//dim 2 - default : - // dimension <= 1 - // Do not use the generic insert. - return Tr::insert(p, c); - } -} - -template < class Gt, class Tds > -typename Delaunay_triangulation_utils_3::Vertex_handle -Delaunay_triangulation_utils_3::move_point(Vertex_handle v, const Point& p) { - CGAL_triangulation_precondition(!is_infinite(v)); - CGAL_triangulation_expensive_precondition(is_vertex(v)); - - // Remember an incident vertex to restart the point location after the removal. - Cell_handle c = v->cell(); - Vertex_handle old_neighbor = c->vertex(c->index(v) == 0 ? 1 : 0); - CGAL_triangulation_assertion(old_neighbor != v); - - if (!remove(v)) - return v; - - if (dimension() <= 0) - return insert(p); - return insert(p, old_neighbor->cell()); -} - -} //end of CGAL namespace - -#endif // CGAL_DELAUNAY_TRIANGULATION_3_UTILS_H diff --git a/Conforming_triangulation_3/include/CGAL/Encroaching_collecter_3.h b/Conforming_triangulation_3/include/CGAL/Encroaching_collecter_3.h deleted file mode 100644 index f371916bb11..00000000000 --- a/Conforming_triangulation_3/include/CGAL/Encroaching_collecter_3.h +++ /dev/null @@ -1,135 +0,0 @@ -//A traverser that collects points encroaching upon its path. -//Copyright (C) 2012 Utrecht University -// -//This program is free software: you can redistribute it and/or modify -//it under the terms of the GNU General Public License as published by -//the Free Software Foundation, either version 3 of the License, or -//(at your option) any later version. -// -//This program is distributed in the hope that it will be useful, -//but WITHOUT ANY WARRANTY; without even the implied warranty of -//MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -//GNU General Public License for more details. -// -//You should have received a copy of the GNU General Public License -//along with this program. If not, see . -// -// Author(s): Thijs van Lankveld - -// A class that traverses the cells of a triangulation by following a segment from source to target, -// while collecting points that encroach upon the segment. - -#ifndef CGAL_ENCROACHING_COLLECTER_3_H -#define CGAL_ENCROACHING_COLLECTER_3_H - -#include "Triangulation_segment_traverser_3.h" - -namespace CGAL { - -#ifndef CGAL_CONSTRAINED_TRIANGULATION_2_H -/// \todo factorize with CDT2 -struct No_intersection_tag{}; -struct Exact_intersections_tag{}; // To be used with an exact number type. -struct Exact_predicates_tag{}; // To be used with filtered exact number type. -#endif - -template < class Gt, class Tds, class Itag > -class Conforming_Delaunay_triangulation_3; - -template < class Gt, class Tds, class Itag = No_intersection_tag, class Out = std::back_insert_iterator > > > -class Encroaching_collecter_3: public Triangulation_segment_traverser_3 { - typedef Encroaching_collecter_3 Self; - typedef Triangulation_segment_traverser_3 Base; - - typedef typename Gt::Plane_3 Plane; - -public: - typedef Conforming_Delaunay_triangulation_3 CDT; - - typedef typename Tds::Vertex Vertex; - typedef typename Tds::Cell Cell; - typedef typename Tds::Edge Edge; - typedef typename Tds::Facet Facet; - typedef typename Tds::Vertex_handle Vertex_handle; - typedef typename Tds::Cell_handle Cell_handle; - typedef typename Tds::Cell_iterator Cell_iterator; - typedef typename Tds::Cell_circulator Cell_circulator; - typedef typename Tds::Facet_iterator Facet_iterator; - typedef typename Tds::Facet_circulator Facet_circulator; - - typedef typename Gt::Point_3 Point; - typedef typename CDT::Locate_type Locate_type; - typedef typename CDT::Bi_vertex Bi_vertex; - using Base::_lt; - using Base::_tr; - using Base::_pos; - using Base::_li; - using Base::_lj; - using Base::_source; - using Base::_target; -protected: - Out _out; - - Bi_vertex se; // These are used when the source lies on an edge. - -public: - Encroaching_collecter_3(): Base() {} - Encroaching_collecter_3(Vertex_handle v, const Point& t, const CDT* tr, Out output) - : Base(v, t, tr), _out(output) {if (_lt == CDT::EDGE) se = _tr->sort_bi_vertex(_pos, _li, _lj);} - Encroaching_collecter_3(const Point& s, const Point& t, const CDT* tr, Out output, Cell_handle hint = Cell_handle()) - : Base(s, t, tr, hint), _out(output) {if (_lt == CDT::EDGE) se = _tr->sort_bi_vertex(_pos, _li, _lj);} - - virtual bool barrier_hit() const {return _lt == CDT::EDGE && _pos->is_conforming(_li, _lj) && _tr->sort_bi_vertex(_pos, _li, _lj) != se;} - -protected: - // Traverse to the next cell along the line. - virtual void increment(); -}; - -template < class Gt, class Tds, class Itag, class Out > -void Encroaching_collecter_3::increment() { - // Walk to the next cell. - Base::increment(); - - // Check for encroaching points. - // Which points are checked is based on the type of simplex traversed. - const CDT* _cdt = dynamic_cast(_tr); - switch (_lt) { - case CDT::FACET: - for (int i = 0; i < 4; ++i) - if (i != _li && !_cdt->is_infinite(_pos->vertex(i)) && _cdt->is_encroaching(_source, _target, _pos->vertex(i)->point())) - *_out++ = _pos->vertex(i)->point(); - break; - case CDT::EDGE: { - Vertex_handle vi = _pos->vertex(_li), vj = _pos->vertex(_lj), vk; - Cell_handle c; - - // Check the vertices in the star around the edge. - Facet_circulator fit = _cdt->incident_facets(_pos, _li, _lj), start(fit); - do { - c = fit->first; - vk = c->vertex(6 - c->index(vi) - c->index(vj) - fit->second); - if (!_cdt->is_infinite(vk) && !collinear(_source, _target, vk->point()) && _cdt->is_encroaching(_source, _target, vk->point())) - *_out++ = vk->point(); - ++fit; - } while (fit != start); - - // Check the vertices of the edge. - if (_cdt->is_encroaching(_source, _target, vi->point())) - *_out++ = vi->point(); - if (_cdt->is_encroaching(_source, _target, vj->point())) - *_out++ = vj->point(); - break; - } - case CDT::VERTEX: - CGAL_triangulation_assertion(collinear(_source, _target, _pos->vertex(_li)->point())); - break; - default: - CGAL_triangulation_assertion(false); - break; - } -} - -} //end of CGAL namespace - -#endif // CGAL_ENCROACHING_COLLECTER_3_H \ No newline at end of file diff --git a/Conforming_triangulation_3/include/CGAL/Triangulation_segment_traverser_3.h b/Conforming_triangulation_3/include/CGAL/Triangulation_segment_traverser_3.h deleted file mode 100644 index fe2c13a2915..00000000000 --- a/Conforming_triangulation_3/include/CGAL/Triangulation_segment_traverser_3.h +++ /dev/null @@ -1,695 +0,0 @@ -//A class that follows a straight line through a Delaunay triangulation structure. -//Copyright (C) 2012 Utrecht University -// -//This program is free software: you can redistribute it and/or modify -//it under the terms of the GNU General Public License as published by -//the Free Software Foundation, either version 3 of the License, or -//(at your option) any later version. -// -//This program is distributed in the hope that it will be useful, -//but WITHOUT ANY WARRANTY; without even the implied warranty of -//MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -//GNU General Public License for more details. -// -//You should have received a copy of the GNU General Public License -//along with this program. If not, see . -// -// Author(s): Thijs van Lankveld - -// A class that traverses the cells of a triangulation by following a segment from source to target. - -#ifndef CGAL_TRIANGULATION_SEGMENT_TRAVERSER_3_H -#define CGAL_TRIANGULATION_SEGMENT_TRAVERSER_3_H - -#include "Delaunay_triangulation_utils_3.h" - -#include - -namespace CGAL { - -template < class Gt, class Tds> -class Delaunay_triangulation_utils_3; - -template < class Gt, class Tds > -class Triangulation_segment_traverser_3: public Forward_circulator_base< typename Tds::Cell, std::ptrdiff_t, std::size_t > { - typedef Triangulation_segment_traverser_3 Self; - -public: - typedef Delaunay_triangulation_utils_3 Tr; - - typedef typename Tds::Vertex Vertex; - typedef typename Tds::Cell Cell; - typedef typename Tds::Edge Edge; - typedef typename Tds::Facet Facet; - typedef typename Tds::Vertex_handle Vertex_handle; - typedef typename Tds::Cell_handle Cell_handle; - typedef typename Tds::Cell_iterator Cell_iterator; - typedef typename Tds::Cell_circulator Cell_circulator; - typedef typename Tds::Facet_iterator Facet_iterator; - typedef typename Tds::Facet_circulator Facet_circulator; - - typedef typename Gt::Point_3 Point; - typedef typename Tr::Locate_type Locate_type; - -protected: - // The triangulation being traversed. - const Tr* _tr; - - // The source and target points of the traversal. - Point _source, _target; - - // The cell currently traversed and the previous one. - Cell_handle _pos, _prev; - - // Otherwise, they indicate the simplex through which this cell was entered, - // or the location of the source if it is in this cell. - int _li, _lj; - Locate_type _lt; - - // This bit signifies when a cell containing the target is found. - bool _done; - - // Where possible, facets are checked in a random order. - mutable Random rng; - -public: - Triangulation_segment_traverser_3(): _tr(NULL), _pos(), _prev(), _li(-1), _lj(-1), _lt(Tr::CELL), _done(true) {} - Triangulation_segment_traverser_3(Vertex_handle s, const Point& t, const Tr* tr); - Triangulation_segment_traverser_3(const Point& s, const Point& t, const Tr* tr, Cell_handle hint = Cell_handle()); - - Self& operator++(); - Self operator++(int); - Cell* operator->() {return &*_pos;} - Cell& operator*() {return *_pos;} - Cell_handle handle() {return _pos;} - operator const Cell_handle() const {return _pos;} - bool operator==(const Self& ct) const; - bool operator!=(const Self& ct) const; - - bool operator==(const Cell_handle& ch) const {return ch == _pos;} - bool operator!=(const Cell_handle& ch) const {return ch != _pos;} - - bool operator==(Nullptr_t CGAL_triangulation_assertion_code(n)) const; - bool operator!=(Nullptr_t n) const; - - // Check if a cell containing the target has been reached. - bool has_next() const {return !_done;} - - // Access the current and previous cells. - Cell_handle previous() const {return _prev;} - - // Traverse to a cell containing the target. - Cell_handle traverse(); - - // Get the type of simplex traversed last. - void traversed(Locate_type& lt, int& li, int& lj) const {lt = _lt; li = _li; lj = _lj;} - -protected: - // Traverse to the next cell along the line. - virtual void increment(); - -private: - inline int edgeIndex(int i, int j) const { - CGAL_triangulation_precondition(i>=0 && i<=3); - CGAL_triangulation_precondition(j>=0 && j<=3); - CGAL_triangulation_precondition(i != j); - return (i==0 || j==0) ? i+j-1 : i+j; - } - - inline Orientation coplanar_orientation(const Point& p, const Point& q, const Point& r, const Point& s) { - return Gt().coplanar_orientation_3_object()(p, q, r, s); - } -}; // class Triangulation_segment_traverser_3 - -template < class Gt, class Tds > -inline bool operator==(typename Tds::Cell_handle ch, Triangulation_segment_traverser_3 tci) {return tci == ch;} - -template < class Gt, class Tds > -inline bool operator!=(typename Tds::Cell_handle ch, Triangulation_segment_traverser_3 tci) {return tci != ch;} - -template < class Gt, class Tds > -Triangulation_segment_traverser_3::Triangulation_segment_traverser_3(Vertex_handle s, const Point& t, const Tr* tr) -: _tr(tr), _pos(), _prev(), _lj(-1), _lt(Tr::VERTEX), _done(false) { - CGAL_triangulation_precondition(!_tr->is_infinite(s)); - CGAL_triangulation_precondition(s->point() != t); - CGAL_triangulation_precondition(_tr->dimension() >= 2); - CGAL_triangulation_assertion(_tr->dimension() == 3 || _tr->plane(*_tr->finite_facets_begin()).has_on(_target)); - - _source = s->point(); - _target = t; - - _pos = s->cell(); - int inf; - if (_pos->has_vertex(_tr->infinite_vertex(), inf)) - _pos = _pos->neighbor(inf); - _li = _pos->index(s); - - CGAL_triangulation_postcondition(_pos != Cell_handle()); -} - -template < class Gt, class Tds > -Triangulation_segment_traverser_3::Triangulation_segment_traverser_3(const Point& s, const Point& t, const Tr* tr, Cell_handle hint) -: _tr(tr), _pos(), _prev(), _li(-1), _lj(-1), _done(false) { - CGAL_triangulation_precondition(s != t); - CGAL_triangulation_precondition(_tr->dimension() >= 2); - CGAL_triangulation_assertion(_tr->dimension() == 3 || _tr->plane(*_tr->finite_facets_begin()).has_on(_target)); - - _source = s; - _target = t; - - _pos = _tr->locate(s, _lt, _li, _lj, hint); - - CGAL_triangulation_postcondition(_pos != Cell_handle()); -} - -template < class Gt, class Tds > -inline Triangulation_segment_traverser_3& -Triangulation_segment_traverser_3::operator++() { - CGAL_triangulation_precondition(_pos != Cell_handle()); - increment(); - return *this; -} - -template < class Gt, class Tds > -inline Triangulation_segment_traverser_3 -Triangulation_segment_traverser_3::operator++(int) { - Self tmp(*this); - ++(*this); - return tmp; -} - -template < class Gt, class Tds > -inline bool Triangulation_segment_traverser_3::operator==(const Self& ct) const { - CGAL_triangulation_precondition(_pos != Cell_handle() && ct._pos != Cell_handle()); - return (_pos == ct._pos && - _prev == ct._prev && - _tr == ct._tr && - _li == ct._li && - _lj == ct._lj && - _lt == ct._lt && - _source == ct._source && - _target == ct._target && - _done == ct._done); -} - -template < class Gt, class Tds > -inline bool Triangulation_segment_traverser_3::operator!=(const Self& ct) const { - return !(*this == ct); -} - -template < class Gt, class Tds > -inline bool Triangulation_segment_traverser_3::operator==(Nullptr_t CGAL_triangulation_assertion_code(n)) const { - CGAL_triangulation_assertion(n == NULL); - return _pos == Cell_handle(); -} - -template < class Gt, class Tds > -inline bool Triangulation_segment_traverser_3::operator!=(Nullptr_t n) const { - CGAL_triangulation_assertion(n == NULL); - return !(*this == n); -} - -template < class Gt, class Tds > -inline typename Triangulation_segment_traverser_3::Cell_handle -Triangulation_segment_traverser_3::traverse() { - while (!_done) - increment(); - return _pos; -} - -template < class Gt, class Tds > -void Triangulation_segment_traverser_3::increment() { - CGAL_triangulation_precondition(!_done); - - // Walks to the next cell over a facet intersected by the line from source to target. - // This method is based on Triangulation_3::locate(). - switch (_tr->dimension()) { - case 3: { - // Infinite cells should be handled differently. - int inf; - if (_pos->has_vertex(_tr->infinite_vertex(), inf)) { - // If this cell was reached by traversal from a finite one, it must be the final cell. - Cell_handle fin = _pos->neighbor(inf); - if (fin == _prev) { - _done = true; - // Here, I do not change the _lt etc. because they were already set. - return; - } - - const Point* pts[4]; - for (int i = 0; i != 4; ++i) - if (i != inf) - pts[i] = &(_pos->vertex(i)->point()); - pts[inf] = &_target; - Orientation o[4]; - - // For the remembering stochastic walk, we start trying with a random index: - int li = rng.template get_bits<2>(); - - // Check if the target lies outside the convex hull. - if (orientation(*pts[0], *pts[1], *pts[2], *pts[3]) == POSITIVE) { - // The target lies in an infinite cell. - _done = true; - return; - - // Note that we do not traverse to other infinite cells. - } - - pts[inf] = &_source; - CGAL_triangulation_assertion(orientation(*pts[0], *pts[1], *pts[2], *pts[3]) == POSITIVE); - - // Check if the line enters an adjacent infinite cell. - // This occurs if the target lies on the other side of - // a plane through one of the finite edges and the source point. - for (int j = 0; j != 4; ++j, li = (li+1)&3) { - if (li == inf) { - o[li] = COPLANAR; - continue; - } - - Cell_handle next = _pos->neighbor(li); - if (next == _prev) { - o[li] = POSITIVE; - continue; - } - - const Point* backup = pts[li]; - pts[li] = &_target; - o[li] = orientation(*pts[0], *pts[1], *pts[2], *pts[3]); - - if (o[li] != NEGATIVE) { - pts[li] = backup; - continue; - } - - // The target lies behind the plane through the source and two finite vertices. - // Traverse to the incident infinite cell. - _prev = _pos; - _pos = next; - CGAL_triangulation_assertion(_tr->is_infinite(_pos)); - - _lt = Tr::FACET; - _li = _pos->index(_prev); - return; - } - - // The line enters the convex hull here (or lies on the finite facet). - _prev = _pos; - _pos = fin; - - // Check through which simplex the line traverses. - switch (o[0]+o[1]+o[2]+o[3]) { - case 3: - _lt = Tr::FACET; - _li = _pos->index(_prev); - return; - case 2: - _lt = Tr::EDGE; - for (int i = 0; i < 4; ++i) { - if (o[i] == COPLANAR && i != inf) { - Edge opp = _tr->opposite_edge(_prev, inf, i); - _li = _pos->index(_prev->vertex(opp.second)); - _lj = _pos->index(_prev->vertex(opp.third)); - return; - } - } - CGAL_triangulation_assertion(false); - return; - case 1: - _lt = Tr::VERTEX; - for (int i = 0; i < 4; ++i) { - if (o[i] == POSITIVE) { - _li = _pos->index(_prev->vertex(i)); - return; - } - } - CGAL_triangulation_assertion(false); - return; - default: - CGAL_triangulation_assertion(false); - return; - } - } - - const Point* pts[4] = {&(_pos->vertex(0)->point()), - &(_pos->vertex(1)->point()), - &(_pos->vertex(2)->point()), - &(_pos->vertex(3)->point())}; - - // We check in which direction the target lies - // by comparing its position relative to the planes through the - // source and the edges of the cell. - Orientation o[6]; - bool calc[6] = {false, false, false, false, false, false}; - - if (_lt == Tr::VERTEX) { - // The three planes through the vertex are set to coplanar. - for (int j = 0; j < 4; ++j) { - if (_li != j) { - int ij = edgeIndex(_li, j); - o[ij] = COPLANAR; - calc[ij] = true; - } - } - } - else if (_lt == Tr::EDGE) { - // The plane through the edge is set to coplanar. - int ij = edgeIndex(_li, _lj); - o[ij] = COPLANAR; - calc[ij] = true; - } - - // For the remembering stochastic walk, we start trying with a random facet: - int li = rng.template get_bits<2>(); - - CGAL_triangulation_assertion_code(bool incell = true;) - for (int k = 0; k < 4; ++k, li = (li+1)&3) { - Cell_handle next = _pos->neighbor(li); - if (next == _prev) - continue; - - // Check if the target is outside the cell. - const Point* backup = pts[li]; - pts[li] = &_target; - if (orientation(*pts[0], *pts[1], *pts[2], *pts[3]) != NEGATIVE) { - pts[li] = backup; - continue; - } - - CGAL_triangulation_assertion_code(incell = false;) - - // Check if the target is inside the pyramid. - int lj = rng.template get_bits<2>(); - int Or = 0; - for (int l = 0; l < 4; ++l, lj = (lj+1)&3) { - if (li == lj) - continue; - - // We check the orientation of the target compared to the plane - // Through the source and the edge opposite of ij. - int oij = 5 - edgeIndex(li, lj); - if (!calc[oij]) { - const Point* backup2 = pts[lj]; - pts[lj] = &_source; - o[oij] = orientation(*pts[0], *pts[1], *pts[2], *pts[3]); - pts[lj] = backup2; - calc[oij] = true; - } - - Or -= o[oij]; - if (o[oij] == POSITIVE) { - // The target is not inside the pyramid. - // Invert the planes. - // This can be safely done because either - // they were not calculated yet, - // or they will no longer be used. - for (int j = 0; j < 4; ++j) { - if (li == j) continue; - int oij = 5 - edgeIndex(li, j); - o[oij] = -o[oij]; - } - Or = 0; - break; - } - } - - if (Or == 0) { - // Either the target is not inside the pyramid, - // or the pyramid is degenerate. - pts[li] = backup; - continue; - } - - // The target is inside the pyramid. - _prev = _pos; - _pos = next; - - switch (Or) { - case 3: - _lt = Tr::FACET; - _li = _pos->index(_prev); - return; - case 2: - _lt = Tr::EDGE; - for (int j = 0; j < 4; ++j) { - if (li != j && o[5-edgeIndex(li, j)] == COPLANAR) { - Edge opp = _tr->opposite_edge(_prev, li, j); - _li = _pos->index(_prev->vertex(opp.second)); - _lj = _pos->index(_prev->vertex(opp.third)); - return; - } - } - CGAL_triangulation_assertion(false); - return; - case 1: - _lt = Tr::VERTEX; - for (int j = 0; j < 4; ++j) { - if (li != j && o[5-edgeIndex(li, j)] == NEGATIVE) { - _li = _pos->index(_prev->vertex(j)); - return; - } - } - CGAL_triangulation_assertion(false); - return; - default: - CGAL_triangulation_assertion(false); - return; - } - } - - // The target lies inside this cell. - CGAL_triangulation_assertion(incell); - _done = true; - return; - } - case 2: { - int inf; - if (_pos->has_vertex(_tr->infinite_vertex(), inf) && inf < 3) { - // If this cell was reached by traversal from a finite one, it must be the final cell. - Cell_handle fin = _pos->neighbor(inf); - if (fin == _prev) { - _done = true; - return; - } - - // Check the neighboring cells. - if (coplanar_orientation(_source, _pos->vertex(_tr->ccw(inf))->point(), _pos->vertex(_tr->cw(inf))->point(), _target) == NEGATIVE) { - Cell_handle tmp = _pos->neighbor(_tr->cw(inf)); - _prev = _pos; - _pos = tmp; - return; - } - if (coplanar_orientation(_source, _pos->vertex(_tr->cw(inf))->point(), _pos->vertex(_tr->ccw(inf))->point(), _target) == NEGATIVE) { - Cell_handle tmp = _pos->neighbor(_tr->ccw(inf)); - _prev = _pos; - _pos = tmp; - return; - } - if (coplanar_orientation(_pos->vertex(_tr->ccw(inf))->point(), _pos->vertex(_tr->cw(inf))->point(), _source, _target) != POSITIVE) { - _prev = _pos; - _pos = fin; - return; - } - - // The target lies in this infinite cell. - _done = true; - return; - } - - const Point* pts[3] = {&(_pos->vertex(0)->point()), - &(_pos->vertex(1)->point()), - &(_pos->vertex(2)->point())}; - - switch (_lt) { - case Tr::VERTEX: { - // First we try the incident edges. - Orientation ocw = coplanar_orientation(*pts[_li], *pts[(_li+2)%3], *pts[(_li+1)%3], _target); - if (_pos->neighbor((_li+1)%3) != _prev && - ocw == NEGATIVE) { - Cell_handle tmp = _pos->neighbor((_li+1)%3); - _prev = _pos; - _pos = tmp; - _li = _pos->index(_prev->vertex(_li)); - return; - } - Orientation occw = coplanar_orientation(*pts[_li], *pts[(_li+1)%3], *pts[(_li+2)%3], _target); - if (_pos->neighbor((_li+2)%3) != _prev && - occw == NEGATIVE) { - Cell_handle tmp = _pos->neighbor((_li+2)%3); - _prev = _pos; - _pos = tmp; - _li = _pos->index(_prev->vertex(_li)); - return; - } - - // Then we try the opposite edge. - if (coplanar_orientation(*pts[(_li+1)%3], *pts[(_li+2)%3], *pts[_li], _target) == NEGATIVE) { - Cell_handle tmp = _pos->neighbor(_li); - _prev = _pos; - _pos = tmp; - - switch (ocw+occw) { - case 2: { - _lt = Tr::EDGE; - _li = _pos->index(_prev->vertex((_li+1)%3)); - _lj = _pos->index(_prev->vertex((_li+2)%3)); - return; - } - case 1: { - _lt = Tr::VERTEX; - if (ocw == COLLINEAR) _li = _pos->index(_prev->vertex((_li+2)%3)); - else _li = _pos->index(_prev->vertex((_li+1)%3)); - return; - } - default: - CGAL_triangulation_assertion(false); - return; - } - } - - // The target lies in this cell. - _done = true; - return; - } - case Tr::EDGE: { - int lk = 3-_li-_lj; - - if (_pos->neighbor(lk) != _prev) { - // Check the edge itself - switch (coplanar_orientation(*pts[_li], *pts[_lj], *pts[lk], _target)) { - case COLLINEAR: { - // The target lies in this cell. - _done = true; - return; - } - case NEGATIVE: { - // The target lies opposite of the edge. - Cell_handle tmp = _pos->neighbor(lk); - _prev = _pos; - _pos = tmp; - _li = _pos->index(_prev->vertex(_li)); - _lj = _pos->index(_prev->vertex(_lj)); - return; - } - default: - break; - } - } - - Orientation o = coplanar_orientation(_source, *pts[lk], *pts[_li], _target); - switch (o) { - case POSITIVE: { - // The ray passes through the edge ik. - if (coplanar_orientation(*pts[lk], *pts[_li], _source, _target) == NEGATIVE) { - Cell_handle tmp = _pos->neighbor(_lj); - _prev = _pos; - _pos = tmp; - - if (collinear(_source, *pts[_li], _target)) { - _lt = Tr::VERTEX; - _li = _pos->index(_prev->vertex(_li)); - } - else { - _lt = Tr::EDGE; - _li = _pos->index(_prev->vertex(_li)); - _lj = _pos->index(_prev->vertex(lk)); - } - return; - } - break; - } - default: { - // The ray passes through the edge jk. - if (coplanar_orientation(*pts[lk], *pts[_lj], _source, _target) == NEGATIVE) { - Cell_handle tmp = _pos->neighbor(_li); - _prev = _pos; - _pos = tmp; - - if (collinear(_source, *pts[_lj], _target)) { - _lt = Tr::VERTEX; - _li = _pos->index(_prev->vertex(_lj)); - } - else if (o == COLLINEAR) { - _lt = Tr::VERTEX; - _li = _pos->index(_prev->vertex(lk)); - } - else { - _lt = Tr::EDGE; - _li = _pos->index(_prev->vertex(lk)); - _lj = _pos->index(_prev->vertex(_lj)); - } - return; - } - break; - } - } - - // The target lies in this cell. - _done = true; - return; - } - case Tr::FACET: { - // We test its edges in a random order until we find a neighbor to go further - int li = rng.get_int(0, 3); - - Orientation o[3]; - bool calc[3] = {false, false, false}; - - for (int j = 0; j != 3; ++j, li = (li+1)%3) { - Cell_handle next = _pos->neighbor(li); - if (next == _prev) - continue; - - // The target should lie on the other side of the edge. - if (coplanar_orientation(*pts[(li+1)%3], *pts[(li+2)%3], *pts[li], _target) != NEGATIVE) - continue; - - // The target should lie inside the wedge. - if (!calc[(li+1)%3]) { - o[(li+1)%3] = coplanar_orientation(_source, *pts[(li+1)%3], *pts[(li+2)%3], _target); - calc[(li+1)%3] = true; - } - if (o[(li+1)%3] == NEGATIVE) continue; - - if (!calc[(li+2)%3]) { - o[(li+2)%3] = coplanar_orientation(_source, *pts[(li+2)%3], *pts[(li+1)%3], _target); - calc[(li+2)%3] = true; - } - if (o[(li+2)%3] == POSITIVE) continue; - - _prev = _pos; - _pos = next; - - switch (o[(li+1)%3]+o[(li+2)%3]) { - case 2: { - _lt = Tr::EDGE; - _li = _pos->index(_prev->vertex((li+1)%3)); - _lj = _pos->index(_prev->vertex((li+2)%3)); - return; - } - case 1: { - _lt = Tr::VERTEX; - if (o[(li+1)%3] == COLLINEAR) _li = _pos->index(_prev->vertex((li+1)%3)); - else _li = _pos->index(_prev->vertex((li+2)%3)); - return; - } - default: - CGAL_triangulation_assertion(false); - return; - } - } - - // The target lies in this cell. - _done = true; - return; - } - default: - CGAL_triangulation_assertion(false); - } - } - } -} - -} //end of CGAL namespace - -#endif // CGAL_TRIANGULATION_SEGMENT_TRAVERSER_3_H \ No newline at end of file diff --git a/Conforming_triangulation_3/package_info/Conforming_triangulation_3/copyright b/Conforming_triangulation_3/package_info/Conforming_triangulation_3/copyright deleted file mode 100644 index dea96f89cdd..00000000000 --- a/Conforming_triangulation_3/package_info/Conforming_triangulation_3/copyright +++ /dev/null @@ -1,2 +0,0 @@ -Utrecht University (The Netherlands) & INRIA Sophia-Antipolis (France) - diff --git a/Conforming_triangulation_3/package_info/Conforming_triangulation_3/description.txt b/Conforming_triangulation_3/package_info/Conforming_triangulation_3/description.txt deleted file mode 100644 index 64a95593a98..00000000000 --- a/Conforming_triangulation_3/package_info/Conforming_triangulation_3/description.txt +++ /dev/null @@ -1,3 +0,0 @@ -Package Conforming_triangulation_3 : -conforming triangulation data structure -this can be wrapped around any CGAL::Triangulation_3, usually CGAL::Delaunay_triangulation_3 diff --git a/Conforming_triangulation_3/package_info/Conforming_triangulation_3/license.txt b/Conforming_triangulation_3/package_info/Conforming_triangulation_3/license.txt deleted file mode 100644 index 8bb8efcb72b..00000000000 --- a/Conforming_triangulation_3/package_info/Conforming_triangulation_3/license.txt +++ /dev/null @@ -1 +0,0 @@ -GPL (v3 or later) diff --git a/Conforming_triangulation_3/package_info/Conforming_triangulation_3/maintainer b/Conforming_triangulation_3/package_info/Conforming_triangulation_3/maintainer deleted file mode 100644 index 8dc617e8bb8..00000000000 --- a/Conforming_triangulation_3/package_info/Conforming_triangulation_3/maintainer +++ /dev/null @@ -1 +0,0 @@ -Thijs van Lankveld diff --git a/Triangulation_3/doc/Triangulation_3/examples.txt b/Triangulation_3/doc/Triangulation_3/examples.txt index 44f09d428a1..d6818e3bfdd 100644 --- a/Triangulation_3/doc/Triangulation_3/examples.txt +++ b/Triangulation_3/doc/Triangulation_3/examples.txt @@ -7,6 +7,7 @@ \example Triangulation_3/info_insert_with_transform_iterator.cpp \example Triangulation_3/simplex.cpp \example Triangulation_3/fast_location_3.cpp +\example Triangulation_3/segment_traverser_3.cpp \example Triangulation_3/find_conflicts_3.cpp \example Triangulation_3/regular_3.cpp \example Triangulation_3/copy_triangulation_3.cpp diff --git a/Triangulation_3/examples/Triangulation_3/README b/Triangulation_3/examples/Triangulation_3/README index d67b701d3d5..968c5a53b65 100644 --- a/Triangulation_3/examples/Triangulation_3/README +++ b/Triangulation_3/examples/Triangulation_3/README @@ -22,6 +22,11 @@ in a triangulation, when he needs to add handles in it. This example shows the use of the Fast_location policy to speed up point location queries in a Delaunay triangulation. +------- segment_traverser_3 -------------------------------------------- + +This example shows how to use a Triangulation_segment_traverser_3 +to traverse the triangulation in a straight line. + ------- tds -------------------------------------------------- Shows how to construct a 3D triangulation data structure by diff --git a/Triangulation_3/examples/Triangulation_3/segment_traverser_3.cpp b/Triangulation_3/examples/Triangulation_3/segment_traverser_3.cpp new file mode 100644 index 00000000000..691faa15751 --- /dev/null +++ b/Triangulation_3/examples/Triangulation_3/segment_traverser_3.cpp @@ -0,0 +1,45 @@ +#include +#include +#include + +#include +#include + +// Define the kernel. +typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel; +typedef Kernel::Point_3 Point_3; + +// Define the structure. +typedef CGAL::Delaunay_triangulation_3< Kernel > DT; + +typedef DT::Cell_handle Cell_handle; + +typedef CGAL::Triangulation_segment_traverser_3< DT > Traverser; + +void main() { + std::vector< Point_3 > points; + points.reserve( 6 ); + points.push_back( Point_3(-.9063, .4226, 1.0 ) ); + points.push_back( Point_3( .8192, .5736, 1.0 ) ); + points.push_back( Point_3( .0872,-.9992, 1.0 ) ); + points.push_back( Point_3(-.8192, .5736,-1.0 ) ); + points.push_back( Point_3( .9063, .4226,-1.0 ) ); + points.push_back( Point_3( .0872,-.9962,-1.0 ) ); + + // Construct the Delaunay triangulation. + DT dt( points.begin(), points.end() ); + assert( dt.is_valid() ); + + // Construct a traverser. + Traverser st( dt, Point_3(-3,0,0), Point_3(3,0,0) ); + + // Count the number of finite cells traversed. + unsigned int cells = 0; + while( st.has_next() ) { + if( !dt.is_infinite(st) ) + ++cells; + ++st; + } + + std::cout << "The number of finite cells between " << st.source() << " and " << st.target() << " is " << cells << std::endl; +} \ No newline at end of file diff --git a/Triangulation_3/include/CGAL/Triangulation_3.h b/Triangulation_3/include/CGAL/Triangulation_3.h index b56c13eeb47..1cbee55c57e 100644 --- a/Triangulation_3/include/CGAL/Triangulation_3.h +++ b/Triangulation_3/include/CGAL/Triangulation_3.h @@ -539,6 +539,7 @@ protected: GT _gt; Vertex_handle infinite; //infinite vertex +public: Comparison_result compare_xyz(const Point &p, const Point &q) const { @@ -557,6 +558,9 @@ protected: { return geom_traits().orientation_3_object()(p, q, r, s); } + + // Compute the orientation of a point compared to the oriented plane supporting a half-facet. + Orientation orientation( const Facet& f, const Point& p ) const; bool coplanar(const Point &p, const Point &q, @@ -565,6 +569,21 @@ protected: return orientation(p, q, r, s) == COPLANAR; } + // Check whether the points of four vertices `a`, `b`, `c`, and `d` are coplanar. + bool coplanar( Vertex_handle a, Vertex_handle b, Vertex_handle c, Vertex_handle d ) const { + return coplanar( a->point(), b->point(), c->point(), d->point() ); + } + + // Check whether the points of facet `f` and point `p` are coplanar. + bool coplanar( const Facet& f, const Point& p ) const { + return orientation(f, p) == COPLANAR; + } + + // Check whether the points of facet `f` and vertex `v` are coplanar. + bool coplanar( const Facet& f, Vertex_handle v ) const{ + return orientation(f, v->point()) == COPLANAR; + } + Orientation coplanar_orientation(const Point &p, const Point &q, const Point &r) const { @@ -577,6 +596,7 @@ protected: return coplanar_orientation(p, q, r) == COLLINEAR; } +protected: Segment construct_segment(const Point &p, const Point &q) const { @@ -1034,6 +1054,22 @@ public: Facet mirror_facet(Facet f) const { return _tds.mirror_facet(f);} + + // Gives the edge incident to the same cell that is not incident to any of the input vertices. + Edge opposite_edge( Cell_handle c, int li, int lj ) const; + + // Gives the edge incident to the same cell that is not incident to any of the vertices of the input edge. + Edge opposite_edge( const Edge& e ) const { return opposite_edge( e.first, e.second, e.third ); } + + // Gives the half-facet incident to the opposite cell. + Facet mirror_facet( Cell_handle c, int li ) const { return Facet( c->neighbor(li), c->neighbor(li)->index(c) ); } + + // Gives the opposite half-facet. + Facet mirror_facet( const Facet& f ) const { return mirror_facet( f.first, f.second ); } + + + + // MODIFIERS bool flip(const Facet &f) // returns false if the facet is not flippable @@ -2244,9 +2280,23 @@ operator<< (std::ostream& os, const Triangulation_3 &tr) return os ; } -template < class GT, class Tds, class Lds > -typename Triangulation_3::size_type -Triangulation_3:: + + + + +template +inline Orientation Triangulation_3:: +orientation( const Facet& f, const Point& p ) const { + return Gt().orientation_3_object()( f.first->vertex( vertex_triple_index(f.second, 0) )->point(), + f.first->vertex( vertex_triple_index(f.second, 1) )->point(), + f.first->vertex( vertex_triple_index(f.second, 2) )->point(), + p ); +} + + +template < class GT, class Tds > +typename Triangulation_3::size_type +Triangulation_3:: number_of_finite_cells() const { if ( dimension() < 3 ) return 0; @@ -3339,7 +3389,30 @@ side_of_edge(const Point & p, } } -template < class GT, class Tds, class Lds > + +template +inline +typename Triangulation_3::Edge +Triangulation_3:: +opposite_edge( Cell_handle c, int li, int lj ) const { + CGAL_triangulation_precondition( li >= 0 && li < 4 ); + CGAL_triangulation_precondition( lj >= 0 && lj < 4 ); + CGAL_triangulation_precondition( li != lj ); + + switch( 6-li-lj ) { // i + j + missing indices = 6. + case 1: return Edge( c, 0, 1 ); + case 2: return Edge( c, 0, 2 ); + case 3: return ( li == 0 || lj == 0 ) ? Edge( c, 1, 2 ) : Edge( c, 0, 3 ); + case 4: return Edge( c, 1, 3 ); + case 5: return Edge( c, 2, 3 ); + } + + CGAL_triangulation_assertion( false ); + return Edge(); +} + + +template < class GT, class Tds > bool Triangulation_3:: flip( Cell_handle c, int i ) diff --git a/Triangulation_3/include/CGAL/Triangulation_segment_traverser_3.h b/Triangulation_3/include/CGAL/Triangulation_segment_traverser_3.h new file mode 100644 index 00000000000..d2fb23730b1 --- /dev/null +++ b/Triangulation_3/include/CGAL/Triangulation_segment_traverser_3.h @@ -0,0 +1,356 @@ +//A class that follows a straight line through a Delaunay triangulation structure. +//Copyright (C) 2012 Utrecht University +// +//This program is free software: you can redistribute it and/or modify +//it under the terms of the GNU General Public License as published by +//the Free Software Foundation, either version 3 of the License, or +//(at your option) any later version. +// +//This program is distributed in the hope that it will be useful, +//but WITHOUT ANY WARRANTY; without even the implied warranty of +//MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +//GNU General Public License for more details. +// +//You should have received a copy of the GNU General Public License +//along with this program. If not, see . +// +// Author(s): Thijs van Lankveld + +#ifndef CGAL_TRIANGULATION_SEGMENT_TRAVERSER_3_H +#define CGAL_TRIANGULATION_SEGMENT_TRAVERSER_3_H + +#include +#include + +namespace CGAL { + +#ifndef DOXYGEN_RUNNING +/** + * \ingroup PkgConformingTriangulation3Classes + * + * The `Triangulation_segment_traverser_3` iterates over the cells + * of a `Triangulation_3` by following a straight line segment \f$ st \f$. + * + * This class is closely related to `Triangulation_3::locate(...)`. + * However, unlike this `locate(...)` method, all the cells traversed + * by the `Triangulation_segment_traverser_3` intersect the interior of the line + * segment \f$ st \f$. + * + * Traversal starts from a cell containing \f$ s \f$ and it ends in a cell containing + * \f$ t \f$. + * If \f$ st \f$ is coplanar with a facet or collinear with an edge, at most one of the + * incident cells is traversed. + * If \f$ st \f$ intersects an edge or vertex, at most two incident cells are traversed: + * the cells intersecting \f$ st \f$ strictly in their interior. + * + * If \f$ s \f$ lies on the convex hull, traversal starts in an incident cell inside + * the convex hull. Similarly, if \f$ t \f$ lies on the convex hull, traversal ends in + * an adjacent cell inside the convex hull. + * + * Both \f$ s \f$ and \f$ t \f$ may lie outside the convex hull of the triangulation, + * but they must lie within the affine hull of the triangulation. In either case, the + * finite facet of any infinite cells traversed must intersect \f$ st \f$. + * + * The traverser may be applied to any triangulation of dimension > 0. + * However, for triangulations of dimension 1, the functionality is somewhat trivial. + * + * The traverser becomes invalid whenever the triangulation is changed. + * + * \tparam Tr_ is the triangulation type to traverse. + * + * \sa `Triangulation_3` + * \sa `Forward_circulator_base` + */ +#endif // DOXYGEN_RUNNING +template < class Tr_ = void > +class Triangulation_segment_traverser_3 +: public Forward_circulator_base< typename Tr_::Triangulation_data_structure::Cell, + std::ptrdiff_t, std::size_t > { + typedef Tr_ Tr; + typedef typename Tr::Triangulation_data_structure Tds; + typedef typename Tr::Geom_traits Gt; + +public: +/// \name Types +/// \{ + typedef Tr Triangulation; ///< The triangulation type. + typedef Triangulation_segment_traverser_3 TST; ///< A triangulation segment traverser. + typedef typename Tr::Point Point; ///< A point embedded in the vertices of the triangulation. + typedef typename Tr::Segment Segment; ///< A segment connecting two points. + typedef typename Tr::Cell Cell; ///< A cell in the triangulation. + typedef typename Tr::Edge Edge; ///< An edge in the triangulation. + typedef typename Tr::Vertex_handle Vertex_handle; ///< A handle for a vertex in the triangulation. + typedef typename Tr::Cell_handle Cell_handle; ///< A handle for a cell in the triangulation. + typedef typename Tr::Locate_type Locate_type; ///< The simplex containing a geometric object. +/// \} + + template < class Tr2 > + struct Rebind_Tr { + typedef Triangulation_segment_traverser_3 Other; + }; + +protected: +/// \internal \name Protected Attributes +/// \{ + /// \internal The triangulation to traverse. + const Tr& _tr; +/// \} + + // The source and target points of the traversal. + Point _source; + Point _target; + Vertex_handle _s_vertex; + Vertex_handle _t_vertex; + + // The cell currently traversed and the previous one. + Cell_handle _pos, _prev; + + // Otherwise, they indicate the simplex through which this cell was entered, + // or the location of the source if it is in this cell. + int _li, _lj; + Locate_type _lt; + + // This bit signifies when a cell containing the target is found. + bool _done; + + // Where possible, facets are checked in a random order. + mutable Random rng; + +private: + Tds _tds2; + Vertex_handle _s_vert, _t_vert; + +public: +/// \name Constructors +/// \{ + /// Constructs a traverser. + /** \param tr the triangulation to traverse. This triangulation must have dimension > 0. + * \param s the source vertex. This vertex must be initialized and cannot be the infinite vertex. + * \param t the target vertex. This vertex must be initialized and cannot be the infinite vertex. + * It cannot equal `s`. + */ + Triangulation_segment_traverser_3( const Tr& tr, Vertex_handle s, Vertex_handle t ); + + /// Constructs a traverser. + /** \param tr the triangulation to traverse. This triangulation must have dimension > 0. + * \param s the source vertex. This vertex must be initialized and cannot be the infinite vertex. + * \param t the target point. This point must be initialized and it cannot be be at the same location as `s`. + * If `tr` has dimension < 3, `t` must lie inside the affine hull of `tr`. + */ + Triangulation_segment_traverser_3( const Tr& tr, Vertex_handle s, const Point& t ); + + /// Constructs a traverser. + /** \param tr the triangulation to traverse. This triangulation must have dimension > 0. + * \param s the source point. This point must be initialized. If `tr` has dimension < 3, `s` must lie inside + * the affine hull of `tr`. + * \param t the target point. This point must be initialized and it cannot be be at the same location as `s`. + * If `tr` has dimension < 3, `t` must lie inside the affine hull of `tr`. + * \param hint the starting point to search for `s`. + */ + Triangulation_segment_traverser_3( const Tr& tr, const Point& s, const Point& t, Cell_handle hint = Cell_handle() ); + + /// Constructs a traverser. + /** \param tr the triangulation to traverse. This triangulation must have dimension > 0. + * \param S the segment to be traversed. If `tr` has dimension < 3, `S` must lie inside + * the affine hull of `tr`. `S` must not be degenerate, i.e. its source and target must not be equal. + * \param hint the starting point to search for `S`. + */ + Triangulation_segment_traverser_3( const Tr& tr, const Segment& S, Cell_handle hint = Cell_handle() ); +/// \} + +public: +/// \name Accessors +/// \{ + /// Gives the source point of the traversed segment. + /** \return the source point. + */ + const Point& source() const { return _source; } + + /// Gives the target point of the traversed segment. + /** \return the target point. + */ + const Point& target() const { return _target; } + + /// Gives the cell currently traversed. + /** By invariance, this cell is intersected by the segment + * between `source()` and `target()`. + * \return the cell currently traversed. + * \sa `handle()`. + */ + Cell_handle cell() const { return _pos; } + + /// Gives a handle to the cell currently traversed. + /** By invariance, this cell is intersected by the segment + * between `source()` and `target()`. + * \return a handle to the cell currently traversed. + * \sa `cell()`. + */ + Cell_handle handle() { return _pos; } + + /// Gives the cell traversed previously. + /** This cell is uninitialized until the traverser leaves the initial + * cell containing the `source()`. + * By invariance, once initialized, this cell must be intersected by the segment + * between `source()` and `target()`. + * \return the cell traversed previously. + * \sa `handle()`. + */ + Cell_handle previous() const { return _prev; } + + /// Dereference operator. + /** \return a pointer to the cell currently traversed. + */ + Cell* operator->() { return &*_pos; } + + /// Indirection operator. + /** \return the cell currently traversed. + */ + Cell& operator*() { return *_pos; } + + /// Conversion operator. + /** \return a handle to the cell currently traversed. + */ + operator const Cell_handle() const { return _pos; } + + /// Checks if the traverser has reached the final cell, which contains the `target()`. + /** If the `target()` lies on a facet, edge, or vertex, the final cell is the cell containing + * the interior of the segment between `source()` and `target()`. + * \return true iff the current cell contains the `target()`. + */ + bool has_next() const { return !_done; } + + /// Gives the type of simplex traversed last. + /** This simplex indicates where traversal has entered the cell. + * For the first cell, containing the `source()` \f$ s \f$, this indicates the location of \f$ s \f$ in this cell. + */ + void entry( Locate_type& lt, int& li, int& lj ) const { lt = _lt; li = _li; lj = _lj; } +/// \} + +public: +/// \name Mutators +/// \{ + /// Increment postfix operator. + /** The current cell must not contain the `target()`. + * After incrementing the traverser, the current cell intersects the segment + * between `source()` and `target()` closer to the `target()` than the previous cell. + * \sa `operator++(int)`. + */ + TST& operator++(); + + /// Increment prefix operator. + /** The current cell must not contain the `target()`. + * After incrementing the traverser, the current cell intersects the segment + * between `source()` and `target()` closer to the `target()` than the previous cell. + * than the previous cell. + * \sa `operator++()`. + */ + TST operator++( int ); + + /// Traverse to the final cell, which contains the `target()`. + /** \return the final cell. + */ + Cell_handle traverse(); +/// \} + +public: +/// \name Comparison +/// \{ + /// Compares this traverser with `ct`. + /** \param ct the other traverser. + * \return true iff the other traverser traverses the same triangulation along the same line segment + * and has the same current cell. + * \sa `operator!=( const TST& t )`. + */ + bool operator==( const TST& ct ) const; + + /// Compares this traverser with `ct`. + /** \param ct the other traverser. + * \return `false` iff the other traverser traverses the same triangulation along the same line segment + * and has the same current cell. + * \sa `operator==( const TST& t ) const`. + */ + bool operator!=( const TST& ct ) const; + + /// Compares the cell currently traversed with `ch`. + /** \param ch a handle to the other cell. + * \return true iff the cell currently traversed is the same as the one pointed to by `ch`. + * \sa `operator!=( const Cell_handle& ch ) const`. + * \sa `operator==( typename TriangulationTraits_3::Cell_handle ch, Triangulation_segment_traverser_3 t )`. + */ + bool operator==( const Cell_handle& ch ) const { return ch == _pos; } + + /// Compares the cell currently traversed with `ch`. + /** \param ch a handle to the other cell. + * \return `false` iff the cell currently traversed is the same as the one pointed to by `ch`. + * \sa `operator==( const Cell_handle& ch )`. + * \sa `operator!=( typename TriangulationTraits_3::Cell_handle ch, Triangulation_segment_traverser_3 t )`. + */ + bool operator!=( const Cell_handle& ch ) const { return ch != _pos; } +/// \} + + bool operator==( Nullptr_t CGAL_triangulation_assertion_code(n) ) const; + bool operator!=( Nullptr_t n ) const; + +protected: +/// \internal \name Protected Member Functions +/// \{ + /// \internal Traverse to the next cell. + /** \internal \sa `traverse()`. + */ + virtual void increment(); +/// \} + +private: + // increment(), if the triangulation is 3D. + void increment_3(); + void increment_3_inf( int inf ); + + // increment(), if the triangulation is 2D. + void increment_2(); + void increment_2_inf( int inf ); + +private: + inline int edgeIndex( int i, int j ) const { + CGAL_triangulation_precondition( i>=0 && i<=3 ); + CGAL_triangulation_precondition( j>=0 && j<=3 ); + CGAL_triangulation_precondition( i != j ); + return ( i==0 || j==0 ) ? i+j-1 : i+j; + } +}; // class Triangulation_segment_traverser_3 + +/// Compares a handle to a cell to a traverser. +/** \param ch the handle to a cell. + * \param t the traverser. + * \return true iff the cell currently traversed by `t` is the same as the one pointed to by `ch`. + * \sa `operator!=( typename TriangulationTraits_3::Cell_handle ch, Triangulation_segment_traverser_3 t )`. + * \sa `Triangulation_segment_traverser_3::operator==( const Cell_handle& ch )`. + */ +template +inline bool operator==( typename Tr::Cell_handle ch, Triangulation_segment_traverser_3 tci ) { return tci == ch; } + +/// Compares a handle to a cell to a traverser. +/** \param ch the handle to a cell. + * \param t the traverser. + * \return `false` iff the cell currently traversed by `t` is the same as the one pointed to by `ch`. + * \sa `operator==( typename TriangulationTraits_3::Cell_handle ch, Triangulation_segment_traverser_3 t )`. + * \sa `Triangulation_segment_traverser_3::operator!=( const Cell_handle& ch )`. + */ +template +inline bool operator!=( typename Tr::Cell_handle ch, Triangulation_segment_traverser_3 tci ) { return tci != ch; } + +// Specialization for void +template <> +class Triangulation_segment_traverser_3 { + typedef internal::Dummy_gt Gt; + typedef internal::Dummy_tds_3 Tds; +public: + typedef Triangulation_3 Triangulation; + template + struct Rebind_Tr { typedef Triangulation_segment_traverser_3 Other; }; +}; // class Triangulation_segment_traverser_3 + +} // namespace CGAL + +#include + +#endif // CGAL_TRIANGULATION_SEGMENT_TRAVERSER_3_H diff --git a/Triangulation_3/include/CGAL/Triangulation_segment_traverser_3_impl.h b/Triangulation_3/include/CGAL/Triangulation_segment_traverser_3_impl.h new file mode 100644 index 00000000000..3134dbef015 --- /dev/null +++ b/Triangulation_3/include/CGAL/Triangulation_segment_traverser_3_impl.h @@ -0,0 +1,665 @@ +//A class that follows a straight line through a Delaunay triangulation structure. +//Copyright (C) 2012 Utrecht University +// +//This program is free software: you can redistribute it and/or modify +//it under the terms of the GNU General Public License as published by +//the Free Software Foundation, either version 3 of the License, or +//(at your option) any later version. +// +//This program is distributed in the hope that it will be useful, +//but WITHOUT ANY WARRANTY; without even the implied warranty of +//MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +//GNU General Public License for more details. +// +//You should have received a copy of the GNU General Public License +//along with this program. If not, see . +// +// Author(s): Thijs van Lankveld + +#ifndef CGAL_TRIANGULATION_SEGMENT_TRAVERSER_3_IMPL_H +#define CGAL_TRIANGULATION_SEGMENT_TRAVERSER_3_IMPL_H + +namespace CGAL { + +template +Triangulation_segment_traverser_3:: +Triangulation_segment_traverser_3( const Tr& tr, Vertex_handle s, const Point& t ) +: _tr(tr), _pos(), _prev(), _lj(-1), _lt(Tr::VERTEX), _done(false) { + CGAL_triangulation_precondition( !_tr.is_infinite(s) ); + CGAL_triangulation_precondition( s->point() != t ); + CGAL_triangulation_precondition( _tr.dimension() >= 2 ); + CGAL_triangulation_precondition( _tr.dimension() == 3 || + _tr.orientation( *_tr.finite_facets_begin(), t ) == COPLANAR ); + + _source = s->point(); + _target = t; + _s_vertex = s; + _t_vertex = Vertex_handle(); + _s_vert = s; + _t_vert = _tds2.create_vertex( Vertex(_target) ); + + _pos = s->cell(); + int inf; + if( _pos->has_vertex( _tr.infinite_vertex(), inf ) ) + _pos = _pos->neighbor(inf); + _li = _pos->index(s); + + CGAL_triangulation_postcondition( _pos != Cell_handle() ); +} + +template < class Tr> +Triangulation_segment_traverser_3:: +Triangulation_segment_traverser_3( const Tr& tr, Vertex_handle s, Vertex_handle t ) +: _tr(tr), _pos(), _prev(), _lj(-1), _lt(Tr::VERTEX), _done(false) { + CGAL_triangulation_precondition( !_tr.is_infinite(s) ); + CGAL_triangulation_precondition( !_tr.is_infinite(t) ); + CGAL_triangulation_precondition( s->point() != t->point() ); + CGAL_triangulation_precondition( _tr.dimension() >= 2 ); + + _source = s->point(); + _target = t->point(); + _s_vertex = s; + _t_vertex = t; + _s_vert = s; + _t_vert = t; + + _pos = s->cell(); + int inf; + if( _pos->has_vertex( _tr.infinite_vertex(), inf ) ) + _pos = _pos->neighbor(inf); + _li = _pos->index(s); +} + + +template +Triangulation_segment_traverser_3:: +Triangulation_segment_traverser_3( const Tr& tr, const Point& s, const Point& t, Cell_handle hint ) +: _tr(tr), _pos(), _prev(), _li(-1), _lj(-1), _done(false) { + CGAL_triangulation_precondition( s != t ); + CGAL_triangulation_precondition( _tr.dimension() >= 2 ); + CGAL_triangulation_precondition( _tr.dimension() == 3 || + _tr.coplanar( *_tr.finite_facets_begin(), _target ) ); + + _source = s; + _target = t; + _s_vertex = Vertex_handle(); + _t_vertex = Vertex_handle(); + _s_vert = _tds2.create_vertex( Tr::Vertex(_source) ); + _t_vert = _tds2.create_vertex( Tr::Vertex(_target) ); + + _pos = _tr.locate( s, _lt, _li, _lj, hint ); + + CGAL_triangulation_postcondition( _pos != Cell_handle() ); +} + +template +Triangulation_segment_traverser_3:: +Triangulation_segment_traverser_3( const Tr& tr, const Segment& s, Cell_handle hint ) +: Triangulation_segment_traverser_3( tr, s.source(), s.target(), hint ) {} + +template +inline Triangulation_segment_traverser_3& +Triangulation_segment_traverser_3::operator++() { + CGAL_triangulation_precondition( _pos != Cell_handle() ); + increment(); + return *this; +} + +template +inline Triangulation_segment_traverser_3 +Triangulation_segment_traverser_3::operator++( int ) { + TST tmp( *this ); + ++( *this ); + return tmp; +} + +template +inline typename Triangulation_segment_traverser_3::Cell_handle +Triangulation_segment_traverser_3::traverse() { + while( !_done ) + increment(); + return _pos; +} + +template +inline bool Triangulation_segment_traverser_3:: +operator==( const TST& ct ) const { + CGAL_triangulation_precondition( _pos != Cell_handle() ); + CGAL_triangulation_precondition( ct._pos != Cell_handle() ); + return ( _pos == ct._pos && + _prev == ct._prev && + _tr == ct._tr && + _li == ct._li && + _lj == ct._lj && + _lt == ct._lt && + _source == ct._source && + _target == ct._target && + _s_vertex == ct._s_vertex && + _t_vertex == ct._t_vertex && + _done == ct._done ); +} + +template +inline bool Triangulation_segment_traverser_3:: +operator!=( const TST& ct ) const { + return !( *this == ct ); +} + +template +inline bool Triangulation_segment_traverser_3:: +operator==( Nullptr_t CGAL_triangulation_assertion_code(n) ) const { + CGAL_triangulation_assertion( n == NULL ); + return _pos == Cell_handle(); +} + +template +inline bool Triangulation_segment_traverser_3:: +operator!=( Nullptr_t n ) const { + return !( *this == n ); +} + +template +void Triangulation_segment_traverser_3:: +increment() { + CGAL_triangulation_precondition( !_done ); + // Walks to the next cell over a facet intersected by the line from source to target. + // This method is based on Triangulation_3::locate(). + int inf; + switch( _tr.dimension() ) { + case 3: { + // Infinite cells should be handled differently. + if( _pos->has_vertex( _tr.infinite_vertex(), inf ) ) + increment_3_inf( inf ); + else + increment_3(); + break; + } + case 2: { + if( _pos->has_vertex( _tr.infinite_vertex(), inf ) ) + increment_2_inf( inf ); + else + increment_2(); + break; + } + } +} + +template +void Triangulation_segment_traverser_3:: +increment_3() { + Vertex_handle vert[4] = { _pos->vertex(0), + _pos->vertex(1), + _pos->vertex(2), + _pos->vertex(3) }; + + // We check in which direction the target lies + // by comparing its position relative to the planes through the + // source and the edges of the cell. + Orientation o[6]; + // We keep track of which orientations are calculated. + bool calc[6] = { false, false, false, false, false, false }; + + if( _lt == Tr::VERTEX ) { + // The three planes through the vertex are set to coplanar. + for( int j = 0; j < 4; ++j ) { + if( _li != j ) { + int ij = edgeIndex( _li, j ); + o[ij] = COPLANAR; + calc[ij] = true; + } + } + } + else if( _lt == Tr::EDGE ) { + // The plane through the edge is set to coplanar. + int ij = edgeIndex( _li, _lj ); + o[ij] = COPLANAR; + calc[ij] = true; + } + + // For the remembering stochastic walk, we start trying with a random facet. + int li = rng.template get_bits<2>(); + + CGAL_triangulation_assertion_code( bool incell = true; ) + for( int k = 0; k < 4; ++k, li = _tr.incr_ind(li) ) { + Cell_handle next = _pos->neighbor(li); + if( next == _prev ) + continue; + + if( _t_vertex == _pos->vertex(0) || _t_vertex == _pos->vertex(1) || + _t_vertex == _pos->vertex(2) || _t_vertex == _pos->vertex(3) ) { + // The target is inside the cell. + _done = true; + return; + } + + // Check if the target is outside the cell. + Vertex_handle backup = vert[li]; + vert[li] = _t_vert; + + if( _tr.orientation( vert[0], vert[1], vert[2], vert[3] ) != NEGATIVE ) { + vert[li] = backup; + continue; + } + CGAL_triangulation_assertion_code( incell = false; ) + + // Check if the target is inside the pyramid. + int lj = rng.template get_bits<2>(); + int Or = 0; + for( int l = 0; l < 4; ++l, lj = _tr.incr_ind(lj) ) { + if( li == lj ) + continue; + + // We check the orientation of the target compared to the plane + // Through the source and the edge opposite of ij. + int oij = 5 - edgeIndex( li, lj ); + if( !calc[oij] ) { + Vertex_handle backup2 = vert[lj]; + vert[lj] = _s_vert; + o[oij] = _tr.orientation( vert[0], vert[1], vert[2], vert[3] ); + vert[lj] = backup2; + calc[oij] = true; + } + + Or -= o[oij]; + if( o[oij] == POSITIVE ) { + // The target is not inside the pyramid. + // Invert the planes. + // This can be safely done because either + // they were not calculated yet, + // or they will no longer be used. + for( int j = 0; j < 4; ++j ) { + if( li == j ) continue; + int oij = 5 - edgeIndex( li, j ); + o[oij] = -o[oij]; + } + Or = 0; + break; + } + } + + if( Or == 0 ) { + // Either the target is not inside the pyramid, + // or the pyramid is degenerate. + vert[li] = backup; + continue; + } + + // The target is inside the pyramid. + _prev = _pos; + _pos = next; + + switch( Or ) { + case 3: + _lt = Tr::FACET; + _li = _pos->index(_prev); + return; + case 2: + _lt = Tr::EDGE; + for( int j = 0; j < 4; ++j ) { + if( li != j && o[ 5 - edgeIndex(li, j) ] == COPLANAR) { + Edge opp = _tr.opposite_edge( _prev, li, j ); + _li = _pos->index( _prev->vertex(opp.second) ); + _lj = _pos->index( _prev->vertex(opp.third) ); + return; + } + } + CGAL_triangulation_assertion( false ); + return; + case 1: + _lt = Tr::VERTEX; + for( int j = 0; j < 4; ++j ) { + if( li != j && o[ 5 - edgeIndex(li, j) ] == NEGATIVE ) { + _li = _pos->index( _prev->vertex(j) ); + return; + } + } + CGAL_triangulation_assertion( false ); + return; + default: + CGAL_triangulation_assertion( false ); + return; + } + } + + // The target lies inside this cell. + CGAL_triangulation_assertion( incell ); + _done = true; + return; +} + +template +void Triangulation_segment_traverser_3:: +increment_3_inf( int inf ) { + CGAL_triangulation_precondition( _tr.is_infinite( _pos->vertex(inf) ) ); + + // If this cell was reached by traversal from a finite one, it must be the final cell. + Cell_handle fin = _pos->neighbor(inf); + if( fin == _prev ) { + _done = true; + // Here, I do not change the _lt etc. because they were already set. + return; + } + + Vertex_handle vert[4]; + for( int i = 0; i != 4; ++i ) + if( i != inf ) + vert[i] = _pos->vertex(i); + vert[inf] = _t_vert; + Orientation o[4]; + + // For the remembering stochastic walk, we start trying with a random index: + int li = rng.template get_bits<2>(); + + // Check if the target lies outside the convex hull. + if( _tr.orientation( vert[0], vert[1], vert[2], vert[3] ) == POSITIVE ) { + // The target lies in an infinite cell. + // Note that we do not traverse to other infinite cells. + _done = true; + return; + } + + vert[inf] = _s_vert; + CGAL_triangulation_assertion( _tr.orientation( vert[0], vert[1], vert[2], vert[3] ) == POSITIVE ); + + // Check if the line enters an adjacent infinite cell. + // This occurs if the target lies on the other side of + // a plane through one of the finite edges and the source point. + for( int j = 0; j != 4; ++j, li = _tr.incr_ind(li) ) { + if( li == inf ) { + o[li] = COPLANAR; + continue; + } + + Cell_handle next = _pos->neighbor(li); + if( next == _prev ) { + o[li] = POSITIVE; + continue; + } + + Vertex_handle backup = vert[li]; + vert[li] = _t_vert; + o[li] = _tr.orientation( vert[0], vert[1], vert[2], vert[3] ); + + if( o[li] != NEGATIVE ) { + vert[li] = backup; + continue; + } + + // The target lies behind the plane through the source and two finite vertices. + // Traverse to the incident infinite cell. + _prev = _pos; + _pos = next; + CGAL_triangulation_assertion( _tr.is_infinite(_pos) ); + + _lt = Tr::FACET; + _li = _pos->index(_prev); + return; + } + + // The line enters the convex hull here (or lies on the finite facet). + _prev = _pos; + _pos = fin; + + // Check through which simplex the line traverses. + switch( o[0]+o[1]+o[2]+o[3] ) { + case 3: + _lt = Tr::FACET; + _li = _pos->index(_prev); + return; + case 2: + _lt = Tr::EDGE; + for( int i = 0; i < 4; ++i ) { + if( o[i] == COPLANAR && i != inf ) { + Edge opp = _tr.opposite_edge( _prev, inf, i ); + _li = _pos->index( _prev->vertex(opp.second) ); + _lj = _pos->index( _prev->vertex(opp.third) ); + return; + } + } + CGAL_triangulation_assertion( false ); + return; + case 1: + _lt = Tr::VERTEX; + for( int i = 0; i < 4; ++i ) { + if( o[i] == POSITIVE ) { + _li = _pos->index( _prev->vertex(i) ); + return; + } + } + CGAL_triangulation_assertion( false ); + return; + default: + CGAL_triangulation_assertion( false ); + return; + } +} + +template +void Triangulation_segment_traverser_3:: +increment_2() { + Vertex_handle vert[3] = { _pos->vertex(0), + _pos->vertex(1), + _pos->vertex(2) }; + + switch( _lt ) { + case Tr::VERTEX: { + // First we try the incident edges. + Orientation ocw = _tr.coplanar_orientation( vert[_li], vert[_tr.cw(_li)], vert[_tr.ccw(_li)], _t_vert ); + if( _pos->neighbor( _tr.ccw(_li) ) != _prev && ocw == NEGATIVE) { + Cell_handle tmp = _pos->neighbor( _tr.ccw(_li) ); + _prev = _pos; + _pos = tmp; + _li = _pos->index( _prev->vertex(_li) ); + return; + } + Orientation occw = _tr.coplanar_orientation( vert[_li], vert[_tr.ccw(_li)], vert[_tr.cw(_li)], _t_vert ); + if( _pos->neighbor( _tr.cw(_li) ) != _prev && occw == NEGATIVE) { + Cell_handle tmp = _pos->neighbor( _tr.cw(_li) ); + _prev = _pos; + _pos = tmp; + _li = _pos->index( _prev->vertex(_li) ); + return; + } + + // Then we try the opposite edge. + if( _tr.coplanar_orientation( vert[_tr.ccw(_li)], vert[_tr.cw(_li)], vert[_li], _t_vert ) == NEGATIVE) { + Cell_handle tmp = _pos->neighbor(_li); + _prev = _pos; + _pos = tmp; + + switch( ocw+occw ) { + case 2: + _lt = Tr::EDGE; + _li = _pos->index( _prev->vertex( _tr.ccw(_li) ) ); + _lj = _pos->index( _prev->vertex( _tr.cw(_li) ) ); + return; + case 1: + _lt = Tr::VERTEX; + if( ocw == COLLINEAR ) _li = _pos->index( _prev->vertex( _tr.cw(_li) ) ); + else _li = _pos->index( _prev->vertex( _tr.ccw(_li) ) ); + return; + default: + CGAL_triangulation_assertion(false); + return; + } + } + + // The target lies in this cell. + _done = true; + return; + } + case Tr::EDGE: { + int lk = 3-_li-_lj; + + if( _pos->neighbor(lk) != _prev ) { + // Check the edge itself + switch( _tr.coplanar_orientation( vert[_li], vert[_lj], vert[lk], _t_vert ) ) { + case COLLINEAR: + // The target lies in this cell. + _done = true; + return; + case NEGATIVE: { + // The target lies opposite of the edge. + Cell_handle tmp = _pos->neighbor(lk); + _prev = _pos; + _pos = tmp; + _li = _pos->index( _prev->vertex(_li) ); + _lj = _pos->index( _prev->vertex(_lj) ); + return; + } + default: + break; + } + } + + Orientation o = _tr.coplanar_orientation( _s_vert, vert[lk], vert[_li], _t_vert ); + switch( o ) { + case POSITIVE: { + // The ray passes through the edge ik. + if( _tr.coplanar_orientation( vert[lk], vert[_li], _s_vert, _t_vert ) == NEGATIVE ) { + Cell_handle tmp = _pos->neighbor(_lj); + _prev = _pos; + _pos = tmp; + + if( _tr.collinear( _s_vert, vert[_li], _t_vert ) ) { + _lt = Tr::VERTEX; + _li = _pos->index( _prev->vertex(_li) ); + } + else { + _lt = Tr::EDGE; + _li = _pos->index( _prev->vertex(_li) ); + _lj = _pos->index( _prev->vertex(lk) ); + } + return; + } + break; + } + default: { + // The ray passes through the edge jk. + if( _tr.coplanar_orientation( vert[lk], vert[_lj], _s_vert, _t_vert ) == NEGATIVE ) { + Cell_handle tmp = _pos->neighbor(_li); + _prev = _pos; + _pos = tmp; + + if( _tr.collinear( _s_vert, vert[_lj], _t_vert ) ) { + _lt = Tr::VERTEX; + _li = _pos->index( _prev->vertex(_lj) ); + } + else if( o == COLLINEAR ) { + _lt = Tr::VERTEX; + _li = _pos->index( _prev->vertex(lk) ); + } + else { + _lt = Tr::EDGE; + _li = _pos->index( _prev->vertex(lk) ); + _lj = _pos->index( _prev->vertex(_lj) ); + } + return; + } + break; + } + } + + // The target lies in this cell. + _done = true; + return; + } + case Tr::FACET: { + // We test its edges in a random order until we find a neighbor to go further + int li = rng.get_int(0, 3); + + Orientation o[3]; + bool calc[3] = { false, false, false }; + + for( int j = 0; j != 3; ++j, li = _tr.ccw(li) ) { + Cell_handle next = _pos->neighbor(li); + if( next == _prev ) + continue; + + // The target should lie on the other side of the edge. + if( _tr.coplanar_orientation( vert[_tr.ccw(li)], vert[_tr.cw(li)], vert[li], _t_vert ) != NEGATIVE ) + continue; + + // The target should lie inside the wedge. + if( !calc[_tr.ccw(li)] ) { + o[_tr.ccw(li)] = _tr.coplanar_orientation( _s_vert, vert[_tr.ccw(li)], vert[_tr.cw(li)], _t_vert ); + calc[_tr.ccw(li)] = true; + } + if( o[_tr.ccw(li)] == NEGATIVE ) continue; + + if( !calc[_tr.cw(li)] ) { + o[_tr.cw(li)] = _tr.coplanar_orientation( _s_vert, vert[_tr.cw(li)], vert[_tr.ccw(li)], _t_vert ); + calc[_tr.cw(li)] = true; + } + if( o[_tr.cw(li)] == POSITIVE ) continue; + + _prev = _pos; + _pos = next; + + switch( o[_tr.ccw(li)] + o[_tr.cw(li)] ) { + case 2: + _lt = Tr::EDGE; + _li = _pos->index( _prev->vertex( _tr.ccw(li) ) ); + _lj = _pos->index( _prev->vertex( _tr.cw(li) ) ); + return; + case 1: + _lt = Tr::VERTEX; + if( o[_tr.ccw(li)] == COLLINEAR ) _li = _pos->index( _prev->vertex( _tr.ccw(li) ) ); + else _li = _pos->index( _prev->vertex( _tr.cw(li) ) ); + return; + default: + CGAL_triangulation_assertion( false ); + return; + } + } + + // The target lies in this cell. + _done = true; + return; + } + default: + CGAL_triangulation_assertion( false ); + } +} + +template +void Triangulation_segment_traverser_3:: +increment_2_inf( int inf ) { + CGAL_triangulation_precondition( _tr.is_infinite( _pos->vertex(3) ) ); + CGAL_triangulation_precondition( _tr.is_infinite( _pos->vertex(inf) ) ); + + // If this cell was reached by traversal from a finite one, it must be the final cell. + Cell_handle fin = _pos->neighbor(inf); + if (fin == _prev) { + _done = true; + return; + } + + // Check the neighboring cells. + Gt::Coplanar_orientation_3 coplanar_orientation = Gt().coplanar_orientation_3_object(); + if( _tr.coplanar_orientation( _s_vert, _pos->vertex( _tr.ccw(inf)), _pos->vertex(_tr.cw(inf)), _t_vert ) == NEGATIVE ) { + Cell_handle tmp = _pos->neighbor(_tr.cw(inf)); + _prev = _pos; + _pos = tmp; + return; + } + if( _tr.coplanar_orientation( _s_vert, _pos->vertex( _tr.cw(inf)), _pos->vertex(_tr.ccw(inf)), _t_vert ) == NEGATIVE ) { + Cell_handle tmp = _pos->neighbor(_tr.ccw(inf)); + _prev = _pos; + _pos = tmp; + return; + } + if( _tr.coplanar_orientation( _pos->vertex( _tr.ccw(inf)), _pos->vertex(_tr.cw(inf)), _s_vert, _t_vert ) != POSITIVE ) { + _prev = _pos; + _pos = fin; + return; + } + + // The target lies in this infinite cell. + _done = true; + return; +} + +} //end of CGAL namespace + +#endif // CGAL_TRIANGULATION_SEGMENT_TRAVERSER_3_IMPL_H diff --git a/Triangulation_3/include/CGAL/Triangulation_utils_3.h b/Triangulation_3/include/CGAL/Triangulation_utils_3.h index 8695c445461..0ac0c642ffe 100644 --- a/Triangulation_3/include/CGAL/Triangulation_utils_3.h +++ b/Triangulation_3/include/CGAL/Triangulation_utils_3.h @@ -29,55 +29,81 @@ namespace CGAL { // We use the following template class in order to avoid having a static data // member of a non-template class which would require src/Triangulation_3.C . template < class T = void > -struct Triangulation_utils_base_3 -{ - static const char tab_next_around_edge[4][4]; - static const int tab_vertex_triple_index[4][3]; +struct Triangulation_utils_base_3 { + static const char tab_next_around_edge[4][4]; + static const int tab_vertex_triple_index[4][3]; + + static const int index_increment_map[4]; + static const int index_jump_map[4]; + static const int index_decrement_map[4]; }; template < class T > const char Triangulation_utils_base_3::tab_next_around_edge[4][4] = { - {5, 2, 3, 1}, - {3, 5, 0, 2}, - {1, 3, 5, 0}, - {2, 0, 1, 5} }; + { 5, 2, 3, 1 }, + { 3, 5, 0, 2 }, + { 1, 3, 5, 0 }, + { 2, 0, 1, 5 } +}; template < class T > const int Triangulation_utils_base_3::tab_vertex_triple_index[4][3] = { - {1, 3, 2}, - {0, 2, 3}, - {0, 3, 1}, - {0, 1, 2} + { 1, 3, 2 }, + { 0, 2, 3 }, + { 0, 3, 1 }, + { 0, 1, 2 } }; +template < class T > +const int Triangulation_utils_base_3::index_increment_map[4] = { 1, 2, 3, 0 }; + +template < class T > +const int Triangulation_utils_base_3::index_jump_map[4] = { 2, 3, 0, 1 }; + +template < class T > +const int Triangulation_utils_base_3::index_decrement_map[4] = { 3, 0, 1, 2 }; + // We derive from Triangulation_cw_ccw_2 because we still use cw() and ccw() // in the 2D part of the code. Ideally, this should go away when we re-use // T2D entirely. struct Triangulation_utils_3 - : public Triangulation_cw_ccw_2, - public Triangulation_utils_base_3<> -{ - static int next_around_edge(const int i, const int j) - { - // index of the next cell when turning around the - // oriented edge vertex(i) vertex(j) in 3d - CGAL_triangulation_precondition( ( i >= 0 && i < 4 ) && - ( j >= 0 && j < 4 ) && - ( i != j ) ); - return tab_next_around_edge[i][j]; - } +: public Triangulation_cw_ccw_2, + public Triangulation_utils_base_3<> { + static int next_around_edge( const int i, const int j ) { + // index of the next cell when turning around the + // oriented edge vertex(i) vertex(j) in 3d + CGAL_triangulation_precondition( ( i >= 0 && i < 4 ) && + ( j >= 0 && j < 4 ) && + ( i != j ) ); + return tab_next_around_edge[i][j]; + } + static int vertex_triple_index( const int i, const int j ) { + // indexes of the jth vertex of the facet of a cell + // opposite to vertx i + CGAL_triangulation_precondition( ( i >= 0 && i < 4 ) && + ( j >= 0 && j < 3 ) ); + return tab_vertex_triple_index[i][j]; + } - static int vertex_triple_index(const int i, const int j) - { - // indexes of the jth vertex of the facet of a cell - // opposite to vertx i - CGAL_triangulation_precondition( ( i >= 0 && i < 4 ) && - ( j >= 0 && j < 3 ) ); - return tab_vertex_triple_index[i][j]; - } + // Get the index of the next vertex or facet. + static int increment_index( int li ) { + CGAL_triangulation_precondition( li >= 0 && li < 4 ); + return index_increment_map[ li ]; + } + // Get the index of the vertex or facet two places further. + static int jump_index( int li ) { + CGAL_triangulation_precondition( li >= 0 && li < 4 ); + return index_jump_map[ li ]; + } + + // Get the index of the previous vertex or facet. + static int decrement_index( int li ) { + CGAL_triangulation_precondition( li >= 0 && li < 4 ); + return index_decrement_map[ li ]; + } }; } //namespace CGAL diff --git a/Triangulation_3/include/CGAL/internal/Dummy_gt.h b/Triangulation_3/include/CGAL/internal/Dummy_gt.h new file mode 100644 index 00000000000..77c7596fac1 --- /dev/null +++ b/Triangulation_3/include/CGAL/internal/Dummy_gt.h @@ -0,0 +1,37 @@ +// A dummy geometry traits class. +// Copyright (c) 2003 INRIA Sophia-Antipolis (France). +// All rights reserved. +// +//This program is free software: you can redistribute it and/or modify +//it under the terms of the GNU General Public License as published by +//the Free Software Foundation, either version 3 of the License, or +//(at your option) any later version. +// +//This program is distributed in the hope that it will be useful, +//but WITHOUT ANY WARRANTY; without even the implied warranty of +//MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +//GNU General Public License for more details. +// +//You should have received a copy of the GNU General Public License +//along with this program. If not, see . +// +// Author(s): Thijs van Lankveld + +#ifndef CGAL_INTERNAL_TRIANGULATION_DUMMY_GT_H +#define CGAL_INTERNAL_TRIANGULATION_DUMMY_GT_H + +namespace CGAL { namespace internal { + +// Dummy triangulation which provides all types that a triangulation can use. +struct Dummy_gt { + struct Point_3 {}; + struct Segment_3 {}; + struct Triangle_3 {}; + struct Tetrahedron_3 {}; +}; + +} // namespace internal + +} // namespace CGAL + +#endif // CGAL_INTERNAL_TRIANGULATION_DUMMY_GT_H